Add support for reading pixel_quality ancillary variables, FCI reader no longer logs warnings

satpy/etc/readers/fci_l1c_fdhsi.yaml

@@ -283,6 +283,101 @@ datasets:
         units: W m-2 um-1 sr-1
     file_type: fci_l1c_fdhsi
 
+  vis_04_pixel_quality:
+    name: vis_04_pixel_quality
+    sensor: fci
+    resolution: 1000
+    file_type: fci_l1c_fdhsi
+
+  vis_05_pixel_quality:
+    name: vis_05_pixel_quality
+    sensor: fci
+    resolution: 1000
+    file_type: fci_l1c_fdhsi
+
+  vis_06_pixel_quality:
+    name: vis_06_pixel_quality
+    sensor: fci
+    resolution: 1000
+    file_type: fci_l1c_fdhsi
+
+  vis_08_pixel_quality:
+    name: vis_08_pixel_quality
+    sensor: fci
+    resolution: 1000
+    file_type: fci_l1c_fdhsi
+
+  vis_09_pixel_quality:
+    name: vis_09_pixel_quality
+    sensor: fci
+    resolution: 1000
+    file_type: fci_l1c_fdhsi
+
+  nir_13_pixel_quality:
+    name: nir_13_pixel_quality
+    sensor: fci
+    resolution: 1000
+    file_type: fci_l1c_fdhsi
+
+  nir_16_pixel_quality:
+    name: nir_16_pixel_quality
+    sensor: fci
+    resolution: 1000
+    file_type: fci_l1c_fdhsi
+
+  nir_22_pixel_quality:
+    name: nir_22_pixel_quality
+    sensor: fci
+    resolution: 1000
+    file_type: fci_l1c_fdhsi
+
+  ir_38_pixel_quality:
+    name: ir_38_pixel_quality
+    sensor: fci
+    resolution: 2000
+    file_type: fci_l1c_fdhsi
+
+  wv_63_pixel_quality:
+    name: wv_63_pixel_quality
+    sensor: fci
+    resolution: 2000
+    file_type: fci_l1c_fdhsi
+
+  wv_73_pixel_quality:
+    name: wv_73_pixel_quality
+    sensor: fci
+    resolution: 2000
+    file_type: fci_l1c_fdhsi
+
+  ir_87_pixel_quality:
+    name: ir_87_pixel_quality
+    sensor: fci
+    resolution: 2000
+    file_type: fci_l1c_fdhsi
+
+  ir_97_pixel_quality:
+    name: ir_97_pixel_quality
+    sensor: fci
+    resolution: 2000
+    file_type: fci_l1c_fdhsi
+
+  ir_105_pixel_quality:
+    name: ir_105_pixel_quality
+    sensor: fci
+    resolution: 2000
+    file_type: fci_l1c_fdhsi
+
+  ir_123_pixel_quality:
+    name: ir_123_pixel_quality
+    sensor: fci
+    resolution: 2000
+    file_type: fci_l1c_fdhsi
+
+  ir_133_pixel_quality:
+    name: ir_133_pixel_quality
+    sensor: fci
+    resolution: 2000
+    file_type: fci_l1c_fdhsi
 
 # Source: FCI L1 Dataset User Guide [FCIL1DUG]
 # ftp://ftp.eumetsat.int/pub/OPS/out/test-data/FCI_L1C_Format_Familiarisation/FCI_L1_Dataset_User_Guide_[FCIL1DUG].pdf

satpy/readers/fci_l1c_fdhsi.py

@@ -56,6 +56,18 @@
 The brightness temperature calculation is based on the formulas indicated in
 `PUG`_ and `RADTOBR`_.
 
+The reading routine supports channel data in counts, radiances, and (depending
+on channel) brightness temperatures or reflectances.  For each channel, it also
+supports the pixel quality, obtained by prepending the channel name such as
+``"vis_04_pixel_quality"``.
+
+.. warning::
+    The API for the direct reading of pixel quality is temporary and likely to
+    change.  Currently, for each channel, the pixel quality is available by
+    ``<chan>_pixel_quality``.  In the future, they will likely all be called
+    ``pixel_quality`` and disambiguated by a to-be-decided property in the
+    `DatasetID`.
+
 .. _RADTOBR: https://www.eumetsat.int/website/wcm/idc/idcplg?IdcService=GET_FILE&dDocName=PDF_EFFECT_RAD_TO_BRIGHTNESS&RevisionSelectionMethod=LatestReleased&Rendition=Web
 .. _PUG: http://www.eumetsat.int/website/wcm/idc/idcplg?IdcService=GET_FILE&dDocName=PDF_DMT_719113&RevisionSelectionMethod=LatestReleased&Rendition=Web
 .. _EUMETSAT: https://www.eumetsat.int/website/home/Satellites/FutureSatellites/MeteosatThirdGeneration/MTGDesign/index.html#fci  # noqa: E501
@@ -126,6 +138,21 @@ def end_time(self):
     def get_dataset(self, key, info=None):
         """Load a dataset."""
         logger.debug('Reading {} from {}'.format(key.name, self.filename))
+        if "pixel_quality" in key.name:
+            return self._get_dataset_quality(key, info=info)
+        elif any(lb in key.name for lb in {"vis_", "ir_", "nir_", "wv_"}):
+            return self._get_dataset_measurand(key, info=info)
+        else:
+            raise ValueError("Unknown dataset key, not a channel or quality: "
+                             f"{key.name:s}")
+
+    def _get_dataset_measurand(self, key, info=None):
+        """Load dataset corresponding to channel measurement.
+
+        Load a dataset when the key refers to a measurand, whether uncalibrated
+        (counts) or calibrated in terms of brightness temperature, radiance, or
+        reflectance.
+        """
         # Get the dataset
         # Get metadata for given dataset
         measured = self.get_channel_measured_group_path(key.name)
@@ -152,6 +179,24 @@ def get_dataset(self, key, info=None):
         info.pop("units")
         attrs.pop("units")
 
+        # For each channel, the effective_radiance contains in the
+        # "ancillary_variables" attribute the value "pixel_quality".  In
+        # FileYAMLReader._load_ancillary_variables, satpy will try to load
+        # "pixel_quality" but is lacking the context from what group to load
+        # it.  Until we can have multiple pixel_quality variables defined (for
+        # example, with https://github.com/pytroll/satpy/pull/1088), rewrite
+        # the ancillary variable to include the channel.  See also
+        # https://github.com/pytroll/satpy/issues/1171.
+        if "pixel_quality" in attrs["ancillary_variables"]:
+            attrs["ancillary_variables"] = attrs["ancillary_variables"].replace(
+                    "pixel_quality", key.name + "_pixel_quality")
+        else:
+            raise ValueError(
+                "Unexpected value for attribute ancillary_variables, "
+                "which the FCI file handler intends to rewrite (see "
+                "https://github.com/pytroll/satpy/issues/1171 for why). "
+                f"Expected 'pixel_quality', got {attrs['ancillary_variables']:s}")
+
         res.attrs.update(key.to_dict())
         res.attrs.update(info)
         res.attrs.update(attrs)
@@ -171,6 +216,26 @@ def get_dataset(self, key, info=None):
 
         return res
 
+    def _get_dataset_quality(self, key, info=None):
+        """Load quality for channel.
+
+        Load a quality field for an FCI channel.  This is a bit involved in
+        case of FCI because each channel group (data/<channel>/measured) has
+        its own data variable 'pixel_quality', referred to in ancillary
+        variables (see also Satpy issue 1171), so some special treatment in
+        necessary.
+        """
+        # FIXME: replace by .removesuffix after we drop support for Python < 3.9
+        if key.name.endswith("_pixel_quality"):
+            chan_lab = key.name[:-len("_pixel_quality")]
+        else:
+            raise ValueError("Quality label must end with pixel_quality, got "
+                             f"{key.name:s}")
+        grp_path = self.get_channel_measured_group_path(chan_lab)
+        dv_path = grp_path + "/pixel_quality"
+        data = self[dv_path]
+        return data
+
     def get_channel_measured_group_path(self, channel):
         """Get the channel's measured group path."""
         measured_group_path = 'data/{}/measured'.format(channel)

satpy/tests/reader_tests/test_fci_l1c_fdhsi.py

@@ -23,6 +23,7 @@
 import dask.array as da
 import numpy.testing
 import pytest
+import logging
 from unittest import mock
 from satpy.tests.reader_tests.test_netcdf_utils import FakeNetCDF4FileHandler
 
@@ -58,6 +59,7 @@ def _get_test_content_for_channel(self, pat, ch):
         chroot = "data/{:s}"
         meas = chroot + "/measured"
         rad = meas + "/effective_radiance"
+        qual = meas + "/pixel_quality"
         pos = meas + "/{:s}_position_{:s}"
         shp = rad + "/shape"
         x = meas + "/x"
@@ -66,6 +68,13 @@ def _get_test_content_for_channel(self, pat, ch):
         ch_str = pat.format(ch)
         ch_path = rad.format(ch_str)
 
+        common_attrs = {
+                "scale_factor": 5,
+                "add_offset": 10,
+                "long_name": "Effective Radiance",
+                "units": "mW.m-2.sr-1.(cm-1)-1",
+                "ancillary_variables": "pixel_quality"
+                }
         if ch == 38:
             fire_line = da.ones((1, ncols), dtype="uint16", chunks=1024) * 5000
             data_without_fires = da.ones((nrows-1, ncols), dtype="uint16", chunks=1024)
@@ -74,12 +83,9 @@ def _get_test_content_for_channel(self, pat, ch):
                 dims=("y", "x"),
                 attrs={
                     "valid_range": [0, 8191],
-                    "scale_factor": 5,
-                    "add_offset": 10,
                     "warm_scale_factor": 2,
                     "warm_add_offset": -300,
-                    "long_name": "Effective Radiance",
-                    "units": "mW.m-2.sr-1.(cm-1)-1",
+                    **common_attrs
                 }
             )
         else:
@@ -88,14 +94,11 @@ def _get_test_content_for_channel(self, pat, ch):
                 dims=("y", "x"),
                 attrs={
                     "valid_range": [0, 4095],
-                    "scale_factor": 5,
-                    "add_offset": 10,
                     "warm_scale_factor": 1,
                     "warm_add_offset": 0,
-                    "long_name": "Effective Radiance",
-                    "units": "mW.m-2.sr-1.(cm-1)-1",
-                }
-            )
+                    **common_attrs
+                    }
+                )
 
         data[ch_path] = d
         data[x.format(ch_str)] = xrda(
@@ -114,6 +117,9 @@ def _get_test_content_for_channel(self, pat, ch):
                     "add_offset": 0.155619515845,
                     }
                 )
+        data[qual.format(ch_str)] = xrda(
+                da.arange(nrows*ncols, dtype="uint8").reshape(nrows, ncols) % 128,
+                dims=("y", "x"))
 
         data[pos.format(ch_str, "start", "row")] = xrda(0)
         data[pos.format(ch_str, "start", "column")] = xrda(0)
@@ -210,6 +216,13 @@ def _get_test_calib_for_channel_ir(self, chroot, meas):
         data[meas + "/radiance_to_bt_conversion_constant_c2"] = v
         return data
 
+    def _get_test_calib_for_channel_vis(self, chroot, meas):
+        data = super()._get_test_calib_for_channel_vis(chroot, meas)
+        from netCDF4 import default_fillvals
+        v = xr.DataArray(default_fillvals["f4"])
+        data[meas + "/channel_effective_solar_irradiance"] = v
+        return data
+
 
 @pytest.fixture
 def reader_configs():
@@ -368,7 +381,7 @@ def test_load_reflectance(self, reader_configs):
             assert res[ch].attrs["units"] == "%"
             numpy.testing.assert_array_equal(res[ch], 100 * 15 * 1 * np.pi / 50)
 
-    def test_load_bt(self, reader_configs):
+    def test_load_bt(self, reader_configs, caplog):
         """Test loading with bt."""
         from satpy import DatasetID
         from satpy.readers import load_reader
@@ -381,10 +394,11 @@ def test_load_bt(self, reader_configs):
         reader = load_reader(reader_configs)
         loadables = reader.select_files_from_pathnames(filenames)
         reader.create_filehandlers(loadables)
-        res = reader.load(
-                [DatasetID(name=name, calibration="brightness_temperature") for
-                    name in self._chans["terran"]])
-        assert 8 == len(res)
+        with caplog.at_level(logging.WARNING):
+            res = reader.load(
+                    [DatasetID(name=name, calibration="brightness_temperature") for
+                        name in self._chans["terran"]])
+            assert caplog.text == ""
         for ch in self._chans["terran"]:
             assert res[ch].shape == (200, 11136)
             assert res[ch].dtype == np.float64
@@ -429,14 +443,40 @@ def test_platform_name(self, reader_configs):
         res = reader.load(["ir_123"])
         assert res["ir_123"].attrs["platform_name"] == "MTG-I1"
 
+    def test_excs(self, reader_configs, caplog):
+        """Test that exceptions are raised where expected."""
+        from satpy import DatasetID
+        from satpy.readers import load_reader
+
+        filenames = [
+            "W_XX-EUMETSAT-Darmstadt,IMG+SAT,MTI1+FCI-1C-RRAD-FDHSI-FD--"
+            "CHK-BODY--L2P-NC4E_C_EUMT_20170410114434_GTT_DEV_"
+            "20170410113925_20170410113934_N__C_0070_0067.nc",
+        ]
+
+        reader = load_reader(reader_configs)
+        loadables = reader.select_files_from_pathnames(filenames)
+        fhs = reader.create_filehandlers(loadables)
+
+        with pytest.raises(ValueError):
+            fhs["fci_l1c_fdhsi"][0].get_dataset(DatasetID(name="invalid"), {})
+        with pytest.raises(ValueError):
+            fhs["fci_l1c_fdhsi"][0]._get_dataset_quality(DatasetID(name="invalid"),
+                                                         {})
+        with caplog.at_level(logging.ERROR):
+            fhs["fci_l1c_fdhsi"][0].get_dataset(
+                    DatasetID(name="ir_123", calibration="unknown"),
+                    {"units": "unknown"})
+            assert "unknown calibration key" in caplog.text
+
 
 class TestFCIL1CFDHSIReaderBadData(TestFCIL1CFDHSIReader):
     """Test the FCI L1C FDHSI Reader for bad data input."""
 
     _alt_handler = FakeNetCDF4FileHandler3
 
     def test_handling_bad_data_ir(self, reader_configs, caplog):
-        """Test handling of bad data."""
+        """Test handling of bad IR data."""
         from satpy import DatasetID
         from satpy.readers import load_reader
 
@@ -454,3 +494,23 @@ def test_handling_bad_data_ir(self, reader_configs, caplog):
                     name="ir_123",
                     calibration="brightness_temperature")])
             assert "cannot produce brightness temperature" in caplog.text
+
+    def test_handling_bad_data_vis(self, reader_configs, caplog):
+        """Test handling of bad VIS data."""
+        from satpy import DatasetID
+        from satpy.readers import load_reader
+
+        filenames = [
+            "W_XX-EUMETSAT-Darmstadt,IMG+SAT,MTI1+FCI-1C-RRAD-FDHSI-FD--"
+            "CHK-BODY--L2P-NC4E_C_EUMT_20170410114434_GTT_DEV_"
+            "20170410113925_20170410113934_N__C_0070_0067.nc",
+        ]
+
+        reader = load_reader(reader_configs)
+        loadables = reader.select_files_from_pathnames(filenames)
+        reader.create_filehandlers(loadables)
+        with caplog.at_level("ERROR"):
+            reader.load([DatasetID(
+                    name="vis_04",
+                    calibration="reflectance")])
+            assert "cannot produce reflectance" in caplog.text