Merge pull request #432 from lsst/tickets/DM-52435

DM-52435: Diffim bad mask plane handling

Author: isullivan

python/lsst/ip/diffim/detectAndMeasure.py

@@ -286,11 +286,13 @@ class DetectAndMeasureConfig(pipeBase.PipelineTaskConfig,
         dtype=str,
         doc="Sources with any of these flags set are removed before writing the output catalog.",
         default=("base_PixelFlags_flag_offimage",
+                 "base_PixelFlags_flag_edge",
                  "base_PixelFlags_flag_interpolatedCenterAll",
                  "base_PixelFlags_flag_badCenterAll",
                  "base_PixelFlags_flag_edgeCenterAll",
                  "base_PixelFlags_flag_nodataCenterAll",
                  "base_PixelFlags_flag_saturatedCenterAll",
+                 "base_PixelFlags_flag_saturated_templateCenterAll",
                  ),
     )
     clearMaskPlanes = lsst.pex.config.ListField(
@@ -371,9 +373,7 @@ def setDefaults(self):
         self.detection.thresholdValue = 5.0
         self.detection.reEstimateBackground = False
         self.detection.thresholdType = "pixel_stdev"
-        self.detection.excludeMaskPlanes = ["EDGE",
-                                            "BAD",
-                                            ]
+        self.detection.excludeMaskPlanes = []
 
         # Copy configs for binned streak detection from the base detection task
         self.streakDetection.thresholdType = self.detection.thresholdType
@@ -416,9 +416,9 @@ def setDefaults(self):
 
         # Keep track of which footprints contain streaks
         self.measurement.plugins["base_PixelFlags"].masksFpAnywhere = [
-            "STREAK", "INJECTED", "INJECTED_TEMPLATE", "HIGH_VARIANCE"]
+            "STREAK", "INJECTED", "INJECTED_TEMPLATE", "HIGH_VARIANCE", "SATURATED_TEMPLATE"]
         self.measurement.plugins["base_PixelFlags"].masksFpCenter = [
-            "STREAK", "INJECTED", "INJECTED_TEMPLATE", "HIGH_VARIANCE"]
+            "STREAK", "INJECTED", "INJECTED_TEMPLATE", "HIGH_VARIANCE", "SATURATED_TEMPLATE"]
         self.skySources.avoidMask = ["DETECTED", "DETECTED_NEGATIVE", "BAD", "NO_DATA", "EDGE"]
 
     def validate(self):
@@ -963,6 +963,9 @@ def addSkySources(self, diaSources, mask, seed,
         """
         if subtask is None:
             subtask = self.skySources
+        if subtask.config.nSources <= 0:
+            self.metadata[f"n_{subtask.getName()}"] = 0
+            return
         skySourceFootprints = subtask.run(mask=mask, seed=seed, catalog=diaSources)
         self.metadata[f"n_{subtask.getName()}"] = len(skySourceFootprints)
 

python/lsst/ip/diffim/utils.py

@@ -410,7 +410,7 @@ def footprint_mean(sources, sky=0):
     else:
         sky_mean = np.nan
         sky_std = np.nan
-        sky_difference = 0
+        sky_difference = np.nanmedian(np.abs(difference.image.array))
     science_footprints, difference_footprints, ratio = footprint_mean(selectStars, sky_difference)
     return lsst.pipe.base.Struct(differenceFootprintRatioMean=ratio.mean(),
                                  differenceFootprintRatioStdev=ratio.std(),

tests/test_detectAndMeasure.py

@@ -111,7 +111,8 @@ def _check_values(self, values, minValue=None, maxValue=None):
             self.assertTrue(np.all(values <= maxValue))
 
     def _setup_detection(self, doSkySources=True, nSkySources=5,
-                         doSubtractBackground=False, run_sattle=False, **kwargs):
+                         doSubtractBackground=False, run_sattle=False,
+                         badSourceFlags=None, **kwargs):
         """Setup and configure the detection and measurement PipelineTask.
 
         Parameters
@@ -143,12 +144,15 @@ def _setup_detection(self, doSkySources=True, nSkySources=5,
             detector=12,
             universe=lsst.daf.butler.DimensionUniverse(),
         )
+        if badSourceFlags is None:
+            badSourceFlags = ["base_PixelFlags_flag_offimage", ]
         config.idGenerator.packer.name = "observation"
         config.idGenerator.packer["observation"].n_observations = 10000
         config.idGenerator.packer["observation"].n_detectors = 99
         config.idGenerator.n_releases = 8
         config.idGenerator.release_id = 2
         config.doSubtractBackground = doSubtractBackground
+        config.badSourceFlags = badSourceFlags
         self.idGenerator = config.idGenerator.apply(dataId)
 
         return self.detectionTask(config=config)
@@ -542,8 +546,12 @@ def test_exclude_mask_detections(self):
         matchedTemplate, _ = makeTestImage(noiseLevel=noiseLevel/4, noiseSeed=7, **kwargs)
 
         # Configure the detection Task
-        detectionTask = self._setup_detection()
-        excludeMaskPlanes = detectionTask.config.detection.excludeMaskPlanes
+        badSourceFlags = ["base_PixelFlags_flag_offimage",
+                          "base_PixelFlags_flag_edgeCenterAll",
+                          "base_PixelFlags_flag_badCenterAll",
+                          "base_PixelFlags_flag_saturatedCenterAll", ]
+        detectionTask = self._setup_detection(nSkySources=0, badSourceFlags=badSourceFlags)
+        excludeMaskPlanes = ["EDGE", "BAD", "SAT"]
         nBad = len(excludeMaskPlanes)
         self.assertGreater(nBad, 0)
         kwargs["seed"] = transientSeed
@@ -563,10 +571,8 @@ def _detection_wrapper(setFlags=True):
                     srcBbox = lsst.geom.Box2I(lsst.geom.Point2I(srcX - radius, srcY - radius),
                                               lsst.geom.Extent2I(2*radius + 1, 2*radius + 1))
                     difference[srcBbox].mask.array |= lsst.afw.image.Mask.getPlaneBitMask(badMask)
-            if setFlags:
                 with self.assertRaises(AlgorithmError):
                     output = detectionTask.run(science, matchedTemplate, difference, sources)
-                return
             else:
                 output = detectionTask.run(science, matchedTemplate, difference, sources)
                 refIds = []
@@ -881,9 +887,7 @@ def test_detection_xy0(self):
 
         self.assertImagesEqual(subtractedMeasuredExposure.image, difference.image)
 
-        # Not all of the sources will be detected: preconvolution results in
-        # a larger edge mask, so we miss an edge source.
-        self.assertEqual(len(output.diaSources), len(sources)-1)
+        self.assertEqual(len(output.diaSources), len(sources))
         # no sources should be flagged as negative
         self.assertEqual(len(~output.diaSources["is_negative"]), len(output.diaSources))
         # TODO DM-41496: restore this block once we handle detections on edge
@@ -1103,13 +1107,19 @@ def test_exclude_mask_detections(self):
         kwargs = {"seed": staticSeed, "psfSize": 2.4, "fluxLevel": fluxLevel}
         science, sources = makeTestImage(noiseLevel=noiseLevel, noiseSeed=6, **kwargs)
         science.getInfo().setVisitInfo(makeVisitInfo())
+        detector = DetectorWrapper(numAmps=1).detector
+        science.setDetector(detector)
         matchedTemplate, _ = makeTestImage(noiseLevel=noiseLevel/4, noiseSeed=7, **kwargs)
 
         subtractTask = subtractImages.AlardLuptonPreconvolveSubtractTask()
         scienceKernel = science.psf.getKernel()
         # Configure the detection Task
-        detectionTask = self._setup_detection()
-        excludeMaskPlanes = detectionTask.config.detection.excludeMaskPlanes
+        badSourceFlags = ["base_PixelFlags_flag_offimage",
+                          "base_PixelFlags_flag_edgeCenterAll",
+                          "base_PixelFlags_flag_badCenterAll",
+                          "base_PixelFlags_flag_saturatedCenterAll", ]
+        detectionTask = self._setup_detection(nSkySources=0, badSourceFlags=badSourceFlags)
+        excludeMaskPlanes = ["EDGE", "BAD", "SAT"]
         nBad = len(excludeMaskPlanes)
         self.assertGreater(nBad, 0)
         kwargs["seed"] = transientSeed
@@ -1133,7 +1143,6 @@ def _detection_wrapper(setFlags=True):
             if setFlags:
                 with self.assertRaises(AlgorithmError):
                     output = detectionTask.run(science, matchedTemplate, difference, score, sources)
-                return
             else:
                 output = detectionTask.run(science, matchedTemplate, difference, score, sources)
                 refIds = []