From 12f7792b48c1b30786dac99eb7d37532b3815b29 Mon Sep 17 00:00:00 2001 From: dwfncar Date: Wed, 16 Oct 2024 13:48:30 -0600 Subject: [PATCH 1/3] Getting started on a VIIRS reader. --- monetio/sat/_viirs_l2_mm.py | 114 ++++++++++++++++++++++++++++++++++++ 1 file changed, 114 insertions(+) create mode 100644 monetio/sat/_viirs_l2_mm.py diff --git a/monetio/sat/_viirs_l2_mm.py b/monetio/sat/_viirs_l2_mm.py new file mode 100644 index 00000000..b09c5c09 --- /dev/null +++ b/monetio/sat/_viirs_l2_mm.py @@ -0,0 +1,114 @@ +import logging +import sys +from collections import OrderedDict +from datetime import datetime, timezone +from glob import glob + +import numpy as np +import xarray as xr + + +def read_dataset(fname, variable_dict): + """ + Parameters + ---------- + fname : str + Input file path. + + Returns + ------- + xarray.Dataset + """ + from monetio.sat.hdfio import hdf_close, hdf_list, hdf_open, hdf_read + + epoch_1993 = int(datetime(1993, 1, 1, tzinfo=timezone.utc).timestamp()) + + print("reading " + fname) + + ds = xr.Dataset() + + f = hdf_open(fname) + hdf_list(f) + # Geolocation and Time Parameters + latitude = hdf_read(f, "Latitude") + longitude = hdf_read(f, "Longitude") + # convert seconds since 1993 to since 1970 + start_time = hdf_read(f, "Scan_Start_Time") + epoch_1993 + for varname in variable_dict: + print(varname) + values = hdf_read(f, varname) + print("min, max: ", values.min(), " ", values.max()) + if "scale" in variable_dict[varname]: + values = variable_dict[varname]["scale"] * values + if "minimum" in variable_dict[varname]: + minimum = variable_dict[varname]["minimum"] + values[values < minimum] = np.nan + if "maximum" in variable_dict[varname]: + maximum = variable_dict[varname]["maximum"] + values[values > maximum] = np.nan + ds[varname] = xr.DataArray(values) + if "quality_flag" in variable_dict[varname]: + ds.attrs["quality_flag"] = varname + ds.attrs["quality_thresh"] = variable_dict[varname]["quality_flag"] + hdf_close(f) + + ds = ds.assign_coords( + { + "lon": (["dim_0", "dim_1"], longitude), + "lat": (["dim_0", "dim_1"], latitude), + "time": (["dim_0", "dim_1"], start_time), + } + ) + ds = ds.rename_dims({"dim_0": "Cell_Along_Swath", "dim_1": "Cell_Across_Swath"}) + + return ds + + +def apply_quality_flag(ds): + """ + Parameters + ---------- + ds : xarray.Dataset + """ + if "quality_flag" in ds.attrs: + quality_flag = ds[ds.attrs["quality_flag"]] + quality_thresh = ds.attrs["quality_thresh"] + for varname in ds: + if varname != ds.attrs["quality_flag"]: + logging.debug(varname) + values = ds[varname].values + values[quality_flag >= quality_thresh] = np.nan + ds[varname].values = values + + +def read_mfdataset(fnames, variable_dict, debug=False): + """ + Parameters + ---------- + fnames : str + Regular expression for input file paths. + + Returns + ------- + xarray.Dataset + """ + if debug: + logging_level = logging.DEBUG + logging.basicConfig(stream=sys.stdout, level=logging_level) + + if isinstance(fnames, str): + files = sorted(glob(fnames)) + else: + files = fnames + + granules = OrderedDict() + + for file in files: + granule = read_dataset(file, variable_dict) + apply_quality_flag(granule) + granule_str = file.split("/")[-1] + granule_info = granule_str.split(".") + datetime_str = granule_info[1][1:] + granule_info[2] + granules[datetime_str] = granule + + return granules From 318c41359e56dd6e776583bfc5d02144bad37545 Mon Sep 17 00:00:00 2001 From: dwfncar Date: Wed, 16 Oct 2024 14:44:26 -0600 Subject: [PATCH 2/3] Added VIIRS reader. --- monetio/sat/_viirs_l2_mm.py | 33 +++++++-------------------------- 1 file changed, 7 insertions(+), 26 deletions(-) diff --git a/monetio/sat/_viirs_l2_mm.py b/monetio/sat/_viirs_l2_mm.py index b09c5c09..85eb673d 100644 --- a/monetio/sat/_viirs_l2_mm.py +++ b/monetio/sat/_viirs_l2_mm.py @@ -19,25 +19,16 @@ def read_dataset(fname, variable_dict): ------- xarray.Dataset """ - from monetio.sat.hdfio import hdf_close, hdf_list, hdf_open, hdf_read - - epoch_1993 = int(datetime(1993, 1, 1, tzinfo=timezone.utc).timestamp()) - print("reading " + fname) - ds = xr.Dataset() + ds_subset = xr.Dataset() + + ds = xr.open_dataset(fname) + print(ds) - f = hdf_open(fname) - hdf_list(f) - # Geolocation and Time Parameters - latitude = hdf_read(f, "Latitude") - longitude = hdf_read(f, "Longitude") - # convert seconds since 1993 to since 1970 - start_time = hdf_read(f, "Scan_Start_Time") + epoch_1993 for varname in variable_dict: print(varname) - values = hdf_read(f, varname) - print("min, max: ", values.min(), " ", values.max()) + values = ds[varname].values if "scale" in variable_dict[varname]: values = variable_dict[varname]["scale"] * values if "minimum" in variable_dict[varname]: @@ -50,18 +41,8 @@ def read_dataset(fname, variable_dict): if "quality_flag" in variable_dict[varname]: ds.attrs["quality_flag"] = varname ds.attrs["quality_thresh"] = variable_dict[varname]["quality_flag"] - hdf_close(f) - - ds = ds.assign_coords( - { - "lon": (["dim_0", "dim_1"], longitude), - "lat": (["dim_0", "dim_1"], latitude), - "time": (["dim_0", "dim_1"], start_time), - } - ) - ds = ds.rename_dims({"dim_0": "Cell_Along_Swath", "dim_1": "Cell_Across_Swath"}) - - return ds + + return ds_subset def apply_quality_flag(ds): From 2c6f117edc5e3299929f3cc1d4f2a648339c3620 Mon Sep 17 00:00:00 2001 From: dwfncar Date: Wed, 16 Oct 2024 14:55:00 -0600 Subject: [PATCH 3/3] Fixed subsetting. --- monetio/sat/_viirs_l2_mm.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/monetio/sat/_viirs_l2_mm.py b/monetio/sat/_viirs_l2_mm.py index 85eb673d..bc3518a6 100644 --- a/monetio/sat/_viirs_l2_mm.py +++ b/monetio/sat/_viirs_l2_mm.py @@ -24,7 +24,7 @@ def read_dataset(fname, variable_dict): ds_subset = xr.Dataset() ds = xr.open_dataset(fname) - print(ds) + # print(ds) for varname in variable_dict: print(varname) @@ -37,10 +37,10 @@ def read_dataset(fname, variable_dict): if "maximum" in variable_dict[varname]: maximum = variable_dict[varname]["maximum"] values[values > maximum] = np.nan - ds[varname] = xr.DataArray(values) + ds_subset[varname] = xr.DataArray(values) if "quality_flag" in variable_dict[varname]: - ds.attrs["quality_flag"] = varname - ds.attrs["quality_thresh"] = variable_dict[varname]["quality_flag"] + ds_subset.attrs["quality_flag"] = varname + ds_subset.attrs["quality_thresh"] = variable_dict[varname]["quality_flag"] return ds_subset