Implement goodtimes culling of entire pointing when gain config does …

imap_processing/hi/hi_goodtimes.py

@@ -71,6 +71,8 @@ def hi_goodtimes(
 
     0. mark_bad_voltage - Remove times with invalid ESA or detector HV
        configuration (e.g. gain test intervals)
+       mark_unmatched_gain_config - Remove the full pointing if its detector
+       gain state matches no gain_config_id in the cal-prod config
     1. mark_incomplete_spin_sets - Remove incomplete 8-spin histogram periods
     2. mark_drf_times - Remove times during spacecraft drift restabilization
     3. mark_bad_tdc_cal - Remove times with failed TDC calibration
@@ -281,6 +283,11 @@ def _apply_goodtimes_filters(
     logger.info("Applying filter: mark_bad_voltage")
     mark_bad_voltage(goodtimes_ds, current_l1b_de)
 
+    # 0b. Mark the full pointing bad if its gain state matches no cal-prod
+    # gain_config_id
+    logger.info("Applying filter: mark_unmatched_gain_config")
+    mark_unmatched_gain_config(goodtimes_ds, current_l1b_de, cal_product_config)
+
     # 1. Mark incomplete spin sets
     logger.info("Applying filter: mark_incomplete_spin_sets")
     mark_incomplete_spin_sets(goodtimes_ds, current_l1b_de)
@@ -1028,6 +1035,69 @@ def mark_bad_voltage(
     )
 
 
+def mark_unmatched_gain_config(
+    goodtimes_ds: xr.Dataset,
+    l1b_de: xr.Dataset,
+    cal_product_config: pd.DataFrame,
+    cull_code: int = CullCode.BAD_HV_VALUE,
+) -> None:
+    """
+    Mark the full pointing bad if its gain state matches no gain_config_id.
+
+    The pointing's detector gain state is recorded by
+    hi_l1b.de_gain_test_filter() as L1B DE global attributes (one per
+    CalibrationProductConfig.GAIN_MATCH_FIELDS). L1C matches these against
+    the cal-prod config to select the geometric factor and the calibration
+    products used for counts (see hi_l1c.generate_pset_dataset()). When no
+    single gain_config_id matches, L1C produces zero counts and fill
+    geometric factors, but exposure time would still be non-zero. Culling
+    every time and spin bin here makes the L1C exposure zero as well, so the
+    pointing is excluded from L2 maps rather than biasing intensities low.
+
+    Parameters
+    ----------
+    goodtimes_ds : xarray.Dataset
+        Goodtimes dataset to update with cull flags.
+    l1b_de : xarray.Dataset
+        L1B Direct Event data for the current pointing, with the gain match
+        global attributes set by hi_l1b.de_gain_test_filter().
+    cal_product_config : pandas.DataFrame
+        Calibration product configuration DataFrame. Use
+        CalibrationProductConfig.from_csv() to load. This should be the same
+        cal-prod file version that L1C uses so that both agree on which
+        pointings are excluded.
+    cull_code : int, optional
+        Cull code to use for marking bad times (default: CullCode.BAD_HV_VALUE).
+
+    Notes
+    -----
+    match_gain_config_id() returns None (and the pointing is culled) when
+    there is no match, more than one match, or any of the gain match values
+    is NaN. A missing global attribute is treated as NaN.
+    """
+    logger.info("Running mark_unmatched_gain_config culling")
+
+    hv_deltas = {
+        field: float(l1b_de.attrs.get(field, np.nan))
+        for field in CalibrationProductConfig.GAIN_MATCH_FIELDS
+    }
+    gain_config_id = cal_product_config.cal_prod_config.match_gain_config_id(hv_deltas)
+    if gain_config_id is not None:
+        logger.info(f"Pointing gain state matches gain_config_id={gain_config_id}")
+        return
+
+    if any(np.isnan(value) for value in hv_deltas.values()):
+        reason = "one or more gain match values are NaN or missing"
+    else:
+        reason = "zero or multiple gain_config_ids matched"
+    logger.warning(
+        f"Pointing gain state {hv_deltas} does not match exactly one "
+        f"gain_config_id in the cal-prod config ({reason}). "
+        f"Marking all times as bad."
+    )
+    goodtimes_ds["cull_flags"].values[:, :] |= np.uint8(cull_code)
+
+
 def mark_incomplete_spin_sets(
     goodtimes_ds: xr.Dataset,
     l1b_de: xr.Dataset,

imap_processing/tests/hi/test_hi_goodtimes.py

@@ -32,9 +32,20 @@
     mark_statistical_filter_0,
     mark_statistical_filter_1,
     mark_statistical_filter_2,
+    mark_unmatched_gain_config,
 )
+from imap_processing.hi.utils import CalibrationProductConfig
 from imap_processing.quality_flags import ImapHiL1bDeFlags
 
+# HV deltas within tolerance of the test cal-prod config's gain_config_id=0.
+# See imap_processing/tests/hi/data/l1/imap_hi_90sensor-cal-prod_20240101_v001.csv
+NOMINAL_HV_DELTAS = {
+    "mcp_delta_v": 875.0,
+    "cem_a_delta_v": 2150.0,
+    "cem_b_delta_v": 2150.0,
+    "tof_v": -8000.0,
+}
+
 
 @pytest.fixture
 def mock_l1b_de():
@@ -4279,6 +4290,7 @@ def test_loads_cal_config(self, tmp_path):
                 "imap_processing.hi.utils.CalibrationProductConfig.from_csv"
             ) as mock_cal_load,
             patch("imap_processing.hi.hi_goodtimes.mark_bad_voltage"),
+            patch("imap_processing.hi.hi_goodtimes.mark_unmatched_gain_config"),
             patch("imap_processing.hi.hi_goodtimes.mark_incomplete_spin_sets"),
             patch("imap_processing.hi.hi_goodtimes.mark_drf_times"),
             patch("imap_processing.hi.hi_goodtimes.mark_overflow_packets"),
@@ -4301,7 +4313,7 @@ def test_loads_cal_config(self, tmp_path):
             mock_cal_load.assert_called_once_with(cal_path)
 
     def test_calls_all_filters(self, tmp_path):
-        """Test that all 8 filters are called."""
+        """Test that all 9 filters are called."""
         mock_goodtimes = MagicMock()
         mock_goodtimes.goodtimes.get_cull_statistics.return_value = {
             "good_bins": 100,
@@ -4317,6 +4329,9 @@ def test_calls_all_filters(self, tmp_path):
                 return_value=mock_cal,
             ),
             patch("imap_processing.hi.hi_goodtimes.mark_bad_voltage") as mock_f0,
+            patch(
+                "imap_processing.hi.hi_goodtimes.mark_unmatched_gain_config"
+            ) as mock_f0b,
             patch(
                 "imap_processing.hi.hi_goodtimes.mark_incomplete_spin_sets"
             ) as mock_f1,
@@ -4343,6 +4358,7 @@ def test_calls_all_filters(self, tmp_path):
             )
 
             mock_f0.assert_called_once()
+            mock_f0b.assert_called_once()
             mock_f1.assert_called_once()
             mock_f2.assert_called_once()
             mock_f3.assert_called_once()
@@ -4368,6 +4384,7 @@ def test_raises_statistical_filter_0_errors(self, tmp_path):
                 return_value=mock_cal,
             ),
             patch("imap_processing.hi.hi_goodtimes.mark_bad_voltage"),
+            patch("imap_processing.hi.hi_goodtimes.mark_unmatched_gain_config"),
             patch("imap_processing.hi.hi_goodtimes.mark_incomplete_spin_sets"),
             patch("imap_processing.hi.hi_goodtimes.mark_drf_times"),
             patch("imap_processing.hi.hi_goodtimes.mark_bad_tdc_cal"),
@@ -4404,6 +4421,7 @@ def test_raises_statistical_filter_1_errors(self, tmp_path):
                 return_value=mock_cal,
             ),
             patch("imap_processing.hi.hi_goodtimes.mark_bad_voltage"),
+            patch("imap_processing.hi.hi_goodtimes.mark_unmatched_gain_config"),
             patch("imap_processing.hi.hi_goodtimes.mark_incomplete_spin_sets"),
             patch("imap_processing.hi.hi_goodtimes.mark_drf_times"),
             patch("imap_processing.hi.hi_goodtimes.mark_bad_tdc_cal"),
@@ -4726,3 +4744,95 @@ def test_mark_bad_voltage_custom_cull_code(
         custom_code = 200
         mark_bad_voltage(goodtimes_for_esa, l1b_de_with_zero, cull_code=custom_code)
         assert np.all(goodtimes_for_esa["cull_flags"].values[1, :] == custom_code)
+
+
+class TestMarkUnmatchedGainConfig:
+    """Tests for mark_unmatched_gain_config culling function."""
+
+    @pytest.fixture
+    def cal_config(self, hi_test_cal_prod_config_path):
+        """Cal-prod config containing only gain_config_id=0."""
+        return CalibrationProductConfig.from_csv(hi_test_cal_prod_config_path)
+
+    @pytest.fixture
+    def goodtimes(self):
+        """Create an all-good goodtimes dataset."""
+        met_values = np.array([1000.0, 1050.0, 1100.0])
+        return xr.Dataset(
+            {
+                "cull_flags": xr.DataArray(
+                    np.zeros((3, 90), dtype=np.uint8),
+                    dims=["met", "spin_bin"],
+                ),
+                "esa_step": xr.DataArray(np.array([1, 2, 3]), dims=["met"]),
+            },
+            coords={"met": met_values, "spin_bin": np.arange(90)},
+        )
+
+    @staticmethod
+    def _l1b_de(hv_deltas):
+        return xr.Dataset(attrs=dict(hv_deltas))
+
+    def test_matching_gain_config_unchanged(self, goodtimes, cal_config):
+        """Test that a pointing matching a gain_config_id is not culled."""
+        mark_unmatched_gain_config(
+            goodtimes, self._l1b_de(NOMINAL_HV_DELTAS), cal_config
+        )
+        assert np.all(goodtimes["cull_flags"].values == CullCode.GOOD)
+
+    def test_no_match_culls_all(self, goodtimes, cal_config):
+        """Test that a pointing matching no gain_config_id is fully culled."""
+        # New gain configuration: CEM A/B back voltages raised by ~200 V
+        hv_deltas = dict(NOMINAL_HV_DELTAS, cem_a_delta_v=2350.0, cem_b_delta_v=2350.0)
+        mark_unmatched_gain_config(goodtimes, self._l1b_de(hv_deltas), cal_config)
+        assert np.all(goodtimes["cull_flags"].values == CullCode.BAD_HV_VALUE)
+        assert goodtimes.goodtimes.get_cull_statistics()["good_bins"] == 0
+        assert len(goodtimes.goodtimes.get_good_intervals()) == 0
+
+    def test_multiple_matches_culls_all(self, goodtimes, cal_config):
+        """Test that a pointing matching multiple gain_config_ids is culled."""
+        duplicate = cal_config.rename(index={0: 1}, level="gain_config_id")
+        ambiguous_config = pd.concat([cal_config, duplicate])
+        mark_unmatched_gain_config(
+            goodtimes, self._l1b_de(NOMINAL_HV_DELTAS), ambiguous_config
+        )
+        assert np.all(goodtimes["cull_flags"].values == CullCode.BAD_HV_VALUE)
+
+    @pytest.mark.parametrize(
+        "field", ["mcp_delta_v", "cem_a_delta_v", "cem_b_delta_v", "tof_v"]
+    )
+    def test_nan_delta_culls_all(self, goodtimes, cal_config, field):
+        """Test that a NaN gain match value culls the full pointing."""
+        hv_deltas = dict(NOMINAL_HV_DELTAS, **{field: np.nan})
+        mark_unmatched_gain_config(goodtimes, self._l1b_de(hv_deltas), cal_config)
+        assert np.all(goodtimes["cull_flags"].values == CullCode.BAD_HV_VALUE)
+
+    def test_missing_attribute_culls_all(self, goodtimes, cal_config):
+        """Test that a missing gain match attribute culls the full pointing."""
+        hv_deltas = dict(NOMINAL_HV_DELTAS)
+        del hv_deltas["tof_v"]
+        mark_unmatched_gain_config(goodtimes, self._l1b_de(hv_deltas), cal_config)
+        assert np.all(goodtimes["cull_flags"].values == CullCode.BAD_HV_VALUE)
+
+    def test_combines_with_existing_culls(self, goodtimes, cal_config):
+        """Test that existing cull codes are preserved via bitwise OR."""
+        goodtimes["cull_flags"].values[0, :] = CullCode.DRF
+        mark_unmatched_gain_config(
+            goodtimes,
+            self._l1b_de({k: np.nan for k in NOMINAL_HV_DELTAS}),
+            cal_config,
+        )
+        assert np.all(
+            goodtimes["cull_flags"].values[0, :] == CullCode.DRF | CullCode.BAD_HV_VALUE
+        )
+        assert np.all(goodtimes["cull_flags"].values[1:, :] == CullCode.BAD_HV_VALUE)
+
+    def test_custom_cull_code(self, goodtimes, cal_config):
+        """Test using a custom cull code."""
+        mark_unmatched_gain_config(
+            goodtimes,
+            self._l1b_de({k: np.nan for k in NOMINAL_HV_DELTAS}),
+            cal_config,
+            cull_code=CullCode.DRF,
+        )
+        assert np.all(goodtimes["cull_flags"].values == CullCode.DRF)