Merge pull request #2802 from bkremmli/read_error

Adapt reader mviri_l1b_fiduceo_nc

AUTHORS.md

@@ -47,6 +47,7 @@ The following people have made contributions to this project:
 - [Sauli Joro (sjoro)](https://github.com/sjoro)
 - [Pouria Khalaj](https://github.com/pkhalaj)
 - [Janne Kotro (jkotro)](https://github.com/jkotro)
+- [Beke Kremmling (bkremmli)](https://github.com/bkremmli) - Deutscher Wetterdienst
 - [Ralph Kuehn (ralphk11)](https://github.com/ralphk11)
 - [Panu Lahtinen (pnuu)](https://github.com/pnuu)
 - [Clement Laplace (ClementLaplace)](https://github.com/ClementLaplace)

satpy/etc/readers/mviri_l1b_fiduceo_nc.yaml

@@ -20,14 +20,18 @@ file_types:
     nc_easy:
         file_reader: !!python/name:satpy.readers.mviri_l1b_fiduceo_nc.FiduceoMviriEasyFcdrFileHandler
         file_patterns: [
-          'FIDUCEO_FCDR_{level}_{sensor}_{platform}-{projection_longitude:f}_{start_time:%Y%m%d%H%M}_{end_time:%Y%m%d%H%M}_EASY_{processor_version}_{format_version}.nc'
+          'FIDUCEO_FCDR_{level}_{sensor}_{platform}-{projection_longitude:f}_{start_time:%Y%m%d%H%M}_{end_time:%Y%m%d%H%M}_EASY_{processor_version}_{format_version}.nc',
           # Example: FIDUCEO_FCDR_L15_MVIRI_MET7-57.0_201701201000_201701201030_EASY_v2.6_fv3.1.nc
+          '{sensor}_FCDR-EASY_{level}_{platform}-E{projection_longitude:s}_{start_time:%Y%m%d%H%M}_{end_time:%Y%m%d%H%M}_{release}.nc'
+          # Example: MVIRI_FCDR-EASY_L15_MET7-E0000_200607060600_200607060630_0200.nc
         ]
     nc_full:
         file_reader: !!python/name:satpy.readers.mviri_l1b_fiduceo_nc.FiduceoMviriFullFcdrFileHandler
         file_patterns: [
-          'FIDUCEO_FCDR_{level}_{sensor}_{platform}-{projection_longitude:f}_{start_time:%Y%m%d%H%M}_{end_time:%Y%m%d%H%M}_FULL_{processor_version}_{format_version}.nc'
+          'FIDUCEO_FCDR_{level}_{sensor}_{platform}-{projection_longitude:f}_{start_time:%Y%m%d%H%M}_{end_time:%Y%m%d%H%M}_FULL_{processor_version}_{format_version}.nc',
           # Example: FIDUCEO_FCDR_L15_MVIRI_MET7-57.0_201701201000_201701201030_FULL_v2.6_fv3.1.nc
+          '{sensor}_FCDR-FULL_{level}_{platform}-E{projection_longitude:s}_{start_time:%Y%m%d%H%M}_{end_time:%Y%m%d%H%M}_{release}.nc'
+          # Example: MVIRI_FCDR-FULL_L15_MET7-E0000_200607060600_200607060630_0200.nc
         ]
 
 datasets:

satpy/readers/mviri_l1b_fiduceo_nc.py

@@ -162,9 +162,9 @@
 
 from satpy.readers._geos_area import get_area_definition, get_area_extent, sampling_to_lfac_cfac
 from satpy.readers.file_handlers import BaseFileHandler
-from satpy.utils import get_legacy_chunk_size
+from satpy.utils import get_chunk_size_limit
 
-CHUNK_SIZE = get_legacy_chunk_size()
+CHUNK_SIZE = get_chunk_size_limit()
 EQUATOR_RADIUS = 6378140.0
 POLE_RADIUS = 6356755.0
 ALTITUDE = 42164000.0 - EQUATOR_RADIUS
@@ -186,6 +186,8 @@
 ]
 HIGH_RESOL = 2250
 
+warnings.filterwarnings("ignore", message="^.*We do not yet support duplicate dimension names, but "
+                                          "we do allow initial construction of the object.*$")
 
 class IRWVCalibrator:
     """Calibrate IR & WV channels."""
@@ -459,6 +461,10 @@ def __init__(self, nc):
         """Wrap the given dataset."""
         self.nc = nc
 
+        self._decode_cf()
+        self._fix_duplicate_dimensions(self.nc)
+
+
     @property
     def attrs(self):
         """Exposes dataset attributes."""
@@ -509,6 +515,31 @@ def _cleanup_attrs(self, ds):
         # satpy warnings.
         ds.attrs.pop("ancillary_variables", None)
 
+    def _decode_cf(self):
+        """Decode data."""
+        # time decoding with decode_cf results in error - decode separately!
+        time_dims, time = self._decode_time()
+        self.nc = self.nc.drop_vars(time.name)
+        self.nc = xr.decode_cf(self.nc)
+        self.nc[time.name] = (time_dims, time.values)
+
+    def _decode_time(self):
+        """Decode time using fill value and offset."""
+        time = self.get_time()
+        time_dims = self.nc[time.name].dims
+        time = xr.where(time == time.attrs["_FillValue"], np.datetime64("NaT"),
+                        (time + time.attrs["add_offset"]).astype("datetime64[s]").astype("datetime64[ns]"))
+
+        return (time_dims, time)
+
+    def _fix_duplicate_dimensions(self, nc):
+        """Rename dimensions as duplicate dimensions names are not supported by xarray."""
+        nc.variables["covariance_spectral_response_function_vis"].dims = ("srf_size_1", "srf_size_2")
+        self.nc = nc.drop_dims("srf_size")
+        nc.variables["channel_correlation_matrix_independent"].dims = ("channel_1", "channel_2")
+        nc.variables["channel_correlation_matrix_structured"].dims = ("channel_1", "channel_2")
+        self.nc = nc.drop_dims("channel")
+
     def get_time(self):
         """Get time coordinate.
 
@@ -558,13 +589,17 @@ def __init__(self, filename, filename_info, filetype_info,  # noqa: D417
             chunks={"x": CHUNK_SIZE,
                     "y": CHUNK_SIZE,
                     "x_ir_wv": CHUNK_SIZE,
-                    "y_ir_wv": CHUNK_SIZE}
+                    "y_ir_wv": CHUNK_SIZE},
+            # see dataset wrapper for why decoding is disabled
+            decode_cf=False,
+            decode_times=False,
+            mask_and_scale=False,
         )
+
         self.nc = DatasetWrapper(nc_raw)
 
-        # Projection longitude is not provided in the file, read it from the
-        # filename.
-        self.projection_longitude = float(filename_info["projection_longitude"])
+        self.projection_longitude = self._get_projection_longitude(filename_info)
+
         self.calib_coefs = self._get_calib_coefs()
 
         self._get_angles = functools.lru_cache(maxsize=8)(
@@ -574,6 +609,13 @@ def __init__(self, filename, filename_info, filetype_info,  # noqa: D417
             self._get_acq_time_uncached
         )
 
+    def _get_projection_longitude(self, filename_info):
+        """Read projection longitude from filename as it is not provided in the file."""
+        if "." in str(filename_info["projection_longitude"]):
+            return float(filename_info["projection_longitude"])
+
+        return float(filename_info["projection_longitude"]) / 100
+
     def get_dataset(self, dataset_id, dataset_info):
         """Get the dataset."""
         name = dataset_id["name"]