VGAC decode the time variable

satpy/etc/readers/viirs_vgac_l1c_nc.yaml

@@ -247,8 +247,13 @@ datasets:
         units: degrees_east
         nc_key: lon
 
-  time:
-    name: time
+  scanline_timestamps:
+    name: scanline_timestamps
     resolution: 5000
     file_type:  vgac_nc
     nc_key: time
+
+  proj_time0:
+    name: proj_time0
+    file_type: vgac_nc
+    nc_key: proj_time0

satpy/readers/viirs_vgac_l1c_nc.py

@@ -63,7 +63,7 @@
 
     def fix_radiances_not_in_percent(self, data):
         """Scale radiances to percent. This was not done in first version of data."""
-        return 100*data
+        return 100 * data
 
     def set_time_attrs(self, data):
         """Set time from attributes."""
@@ -73,6 +73,33 @@
             self._end_time = data.attrs["end_time"]
             self._start_time = data.attrs["start_time"]
 
+    def dt64_to_datetime(self, dt64):
+        """Conversion of numpy.datetime64 to datetime objects."""
+        if isinstance(dt64, np.datetime64):
+            return dt64.astype(datetime)
+        return dt64
+
+    def extract_time_data(self, data, nc):
+        """Decode time data."""
+        reference_time = np.datetime64(datetime.strptime(nc["proj_time0"].attrs["units"],
+                                                         "days since %d/%m/%YT%H:%M:%S"))
+        delta_part_of_day, delta_full_days = np.modf(nc["proj_time0"].values[0])
+        delta_full_days = np.timedelta64(delta_full_days.astype(np.int64), "D").astype("timedelta64[us]")
+        delta_part_of_day = delta_part_of_day * np.timedelta64(1, "D").astype("timedelta64[us]")
+        delta_hours = data.values * np.timedelta64(1, "h").astype("timedelta64[us]")
+        time_data = xr.DataArray(reference_time + delta_full_days + delta_part_of_day + delta_hours,
+                                 coords=data.coords, attrs={"long_name": "Scanline time"})
+        return time_data
+
+    def decode_time_variable(self, data, file_key, nc):
+        """Decide if time data should be decoded."""
+        if file_key != "time":
+            return data
+        if data.attrs["units"] == "hours since proj_time0":
+            return self.extract_time_data(data, nc)
+        else:
+            raise AttributeError('Unit of time variable in VGAC nc file is not "hours since proj_time0"')
+
     def get_dataset(self, key, yaml_info):
         """Get dataset."""
         logger.debug("Getting data for: %s", yaml_info["name"])
@@ -82,6 +109,7 @@
         file_key = yaml_info.get("nc_key", name)
         data = nc[file_key]
         data = self.calibrate(data, yaml_info, file_key, nc)
+        data = self.decode_time_variable(data, file_key, nc)
         data.attrs.update(nc.attrs)  # For now add global attributes to all datasets
         data.attrs.update(yaml_info)
         self.set_time_attrs(data)

satpy/tests/reader_tests/test_viirs_vgac_l1c_nc.py

@@ -22,29 +22,33 @@
 """
 
 
-import datetime
+from datetime import datetime
 
 import numpy as np
 import pytest
+import xarray as xr
 from netCDF4 import Dataset
 
 
 @pytest.fixture()
 def nc_filename(tmp_path):
     """Create an nc test data file and return its filename."""
-    now = datetime.datetime.utcnow()
+    now = datetime.utcnow()
     filename = f"VGAC_VJ10XMOD_A{now:%Y%j_%H%M}_n004946_K005.nc"
     filename_str = str(tmp_path / filename)
     # Create test data
     with Dataset(filename_str, "w") as nc:
         nscn = 7
         npix = 800
         n_lut = 12000
+        start_time_srting = "2023-03-28T09:08:07"
+        end_time_string = "2023-03-28T10:11:12"
         nc.createDimension("npix", npix)
         nc.createDimension("nscn", nscn)
         nc.createDimension("n_lut", n_lut)
-        nc.StartTime = "2023-03-28T09:08:07"
-        nc.EndTime = "2023-03-28T10:11:12"
+        nc.createDimension("one", 1)
+        nc.StartTime = start_time_srting
+        nc.EndTime = end_time_string
         for ind in range(1, 11, 1):
             ch_name = "M{:02d}".format(ind)
             r_a = nc.createVariable(ch_name, np.int16, dimensions=("nscn", "npix"))
@@ -61,6 +65,23 @@ def nc_filename(tmp_path):
                 setattr(tb_b, attr, attrs[attr])
             tb_lut = nc.createVariable(ch_name + "_LUT", np.float32, dimensions=("n_lut"))
             tb_lut[:] = np.array(range(0, n_lut)) * 0.5
+            tb_lut.units = "Kelvin"
+        reference_time = np.datetime64("2010-01-01T00:00:00")
+        start_time = np.datetime64("2023-03-28T09:08:07") + np.timedelta64(123000, "us")
+        delta_days = start_time - reference_time
+        delta_full_days = delta_days.astype("timedelta64[D]")
+        hidden_reference_time = reference_time + delta_full_days
+        delta_part_of_days = start_time - hidden_reference_time
+        proj_time0 = nc.createVariable("proj_time0", np.float64, ("one",))
+        proj_time0[:] = (delta_full_days.astype(np.int64) +
+                         0.000001 * delta_part_of_days.astype("timedelta64[us]").astype(np.int64) / (60 * 60 * 24))
+        proj_time0.units = "days since 01/01/2010T00:00:00"
+        time_v = nc.createVariable("time", np.float64, ("nscn",))
+        delta_h = np.datetime64(end_time_string) - start_time
+        delta_hours = 0.000001 * delta_h.astype("timedelta64[us]").astype(np.int64) / (60 * 60)
+        time_v[:] = np.linspace(0, delta_hours, num=nscn).astype(np.float64)
+        time_v.units = "hours since proj_time0"
+
     return filename_str
 
 
@@ -75,10 +96,43 @@ def test_read_vgac(self, nc_filename):
         scn_ = Scene(
             reader="viirs_vgac_l1c_nc",
             filenames=[nc_filename])
-        scn_.load(["M05", "M15"])
+        scn_.load(["M05", "M15", "scanline_timestamps"])
+        diff_s = (scn_["scanline_timestamps"][0].values.astype("datetime64[us]") -
+                  np.datetime64("2023-03-28T09:08:07.123000").astype("datetime64[us]"))
+        diff_e = (np.datetime64("2023-03-28T10:11:12.000000").astype("datetime64[us]") -
+                  scn_["scanline_timestamps"][-1].values.astype("datetime64[us]"))
+        assert (diff_s < np.timedelta64(5, "us"))
+        assert (diff_s > np.timedelta64(-5, "us"))
+        assert (diff_e < np.timedelta64(5, "us"))
+        assert (diff_e > np.timedelta64(-5, "us"))
         assert (scn_["M05"][0, 0] == 100)
         assert (scn_["M15"][0, 0] == 400)
-        assert scn_.start_time == datetime.datetime(year=2023, month=3, day=28,
-                                                    hour=9, minute=8, second=7)
-        assert scn_.end_time == datetime.datetime(year=2023, month=3, day=28,
-                                                  hour=10, minute=11, second=12)
+        assert scn_.start_time == datetime(year=2023, month=3, day=28,
+                                           hour=9, minute=8, second=7)
+        assert scn_.end_time == datetime(year=2023, month=3, day=28,
+                                         hour=10, minute=11, second=12)
+
+    def test_dt64_to_datetime(self):
+        """Test datetime conversion branch."""
+        from satpy.readers.viirs_vgac_l1c_nc import VGACFileHandler
+        fh = VGACFileHandler(filename="",
+                             filename_info={"start_time": "2023-03-28T09:08:07"},
+                             filetype_info="")
+        in_dt = datetime(year=2023, month=3, day=28,
+                         hour=9, minute=8, second=7)
+        out_dt = fh.dt64_to_datetime(in_dt)
+        assert out_dt == in_dt
+
+    def test_decode_time_variable(self):
+        """Test decode time variable branch."""
+        from satpy.readers.viirs_vgac_l1c_nc import VGACFileHandler
+        fh = VGACFileHandler(filename="",
+                             filename_info={"start_time": "2023-03-28T09:08:07"},
+                             filetype_info="")
+        data = xr.DataArray(
+            [[1, 2],
+             [3, 4]],
+            dims=("y", "x"),
+            attrs={"units": "something not expected"})
+        with pytest.raises(AttributeError):
+            fh.decode_time_variable(data, "time", None)