DM-49210: Set is_negative flag for negative-detected diffim sources

python/lsst/ip/diffim/detectAndMeasure.py

@@ -324,6 +324,10 @@ def __init__(self, **kwargs):
             self.makeSubtask("maskStreaks")
             self.makeSubtask("streakDetection")
 
+        # To get the "merge_*" fields in the schema; have to re-initialize
+        # this later, once we have a peak schema post-detection.
+        lsst.afw.detection.FootprintMergeList(self.schema, ["positive", "negative"])
+
         # Check that the schema and config are consistent
         for flag in self.config.badSourceFlags:
             if flag not in self.schema:
@@ -474,20 +478,24 @@ def processResults(self, science, matchedTemplate, difference, sources, idFactor
         if self.config.doMerge:
             # preserve peak schema, if there are any footprints
             if len(positives) > 0:
-                schema = positives[0].getFootprint().peaks.schema
+                peakSchema = positives[0].getFootprint().peaks.schema
             elif len(negatives) > 0:
-                schema = negatives[0].getFootprint().peaks.schema
+                peakSchema = negatives[0].getFootprint().peaks.schema
             else:
-                schema = afwDetection.PeakTable.makeMinimalSchema()
-            mergeList = lsst.afw.detection.FootprintMergeList(afwTable.SourceTable.makeMinimalSchema(),
-                                                              ["positive", "negative"], schema)
-            initialDiaSources = lsst.afw.table.SourceCatalog(positives.schema)
+                peakSchema = afwDetection.PeakTable.makeMinimalSchema()
+            mergeList = afwDetection.FootprintMergeList(self.schema,
+                                                        ["positive", "negative"], peakSchema)
+            initialDiaSources = afwTable.SourceCatalog(self.schema)
             # Start with positive, as FootprintMergeList will self-merge the
             # subsequent added catalogs, and we want to try to preserve
             # deblended positive sources.
             mergeList.addCatalog(initialDiaSources.table, positives, "positive", minNewPeakDist=0)
             mergeList.addCatalog(initialDiaSources.table, negatives, "negative", minNewPeakDist=0)
             mergeList.getFinalSources(initialDiaSources)
+            # Flag as negative those sources that *only* came from the negative
+            # footprint set.
+            initialDiaSources["is_negative"] = initialDiaSources["merge_footprint_negative"] & \
+                ~initialDiaSources["merge_footprint_positive"]
             self.log.info("Merging detections into %d sources", len(initialDiaSources))
         else:
             initialDiaSources = sources
@@ -592,6 +600,8 @@ def deblend(footprints, negative=False):
         sources.reserve(len(positives) + len(negatives))
         sources.extend(positives, deep=True)
         sources.extend(negatives, deep=True)
+        if len(negatives) > 0:
+            sources[-len(negatives):]["is_negative"] = True
         return sources, positives, negatives
 
     def _removeBadSources(self, diaSources):

tests/test_detectAndMeasure.py

@@ -168,6 +168,8 @@ def test_detection_xy0(self):
 
         # all of the sources should have been detected
         self.assertEqual(len(output.diaSources), len(sources))
+        # no sources should be flagged as negative
+        self.assertEqual(len(~output.diaSources["is_negative"]), len(output.diaSources))
         refIds = []
         for source in sources:
             self._check_diaSource(output.diaSources, source, refIds=refIds)
@@ -383,6 +385,8 @@ def test_detect_dipoles(self):
         detectionTask = self._setup_detection(doMerge=True)
         output = detectionTask.run(science, matchedTemplate, difference)
         self.assertEqual(len(output.diaSources), len(sources))
+        # no sources should be flagged as negative
+        self.assertEqual(len(~output.diaSources["is_negative"]), len(output.diaSources))
         refIds = []
         for diaSource in output.diaSources:
             if diaSource[dipoleFlag]:
@@ -668,6 +672,8 @@ def test_detection_xy0(self):
         # 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)
+        # 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
         # pixels better; at least one of these sources currently has a bad
         # centroid because most of the source is rejected as EDGE.
@@ -977,32 +983,47 @@ def testDeblendNegatives(self):
         self.assertFloatsAlmostEqual(negatives[1].getFootprint().getImageArray().min(),
                                      -catalog[2].getFootprint().getImageArray().max(), rtol=1e-4)
 
+        self.assertEqual(sources["is_negative"].sum(), 2)
+
     def testMergeFootprints(self):
         """Test that merging footprints does not cause negative ones to
         disappear (e.g. get merged into non-connected footprints).
+
+        As implemented, the diffim will have three positive sources (one a
+        blended pair), and 7 negative sources (two a blended pair).
         """
-        # Make a science image with multiple blends, designed to trigger the
+        # Make a template image with multiple blends, designed to trigger the
         # negative-footprint-related bug in `FootprintSet.merge`.
         bbox = lsst.geom.Box2I(lsst.geom.Point2I(0, 0), lsst.geom.Point2I(400, 200))
         dataset = lsst.meas.base.tests.TestDataset(bbox)
+        # Because these sources are on the template, they will be negative on
+        # the diffim, unless the pixels are explicitly made negative via
+        # `template.image.subset` below.
+        # positive isolated source on diffim
         dataset.addSource(instFlux=.5E5, centroid=lsst.geom.Point2D(25, 26))
+        # negative isolated source on diffim
         dataset.addSource(instFlux=.7E5, centroid=lsst.geom.Point2D(75, 24),
                           shape=lsst.afw.geom.Quadrupole(12, 7, 2))
         delta = 10
+        # negative blended pair on diffim
         with dataset.addBlend() as family:
             family.addChild(instFlux=1E5, centroid=lsst.geom.Point2D(150, 72))
             family.addChild(instFlux=1.5E5, centroid=lsst.geom.Point2D(150+delta, 74))
+        # positive blended pair on diffim
         with dataset.addBlend() as family:
             family.addChild(instFlux=2E5, centroid=lsst.geom.Point2D(250, 72))
             family.addChild(instFlux=2.5E5, centroid=lsst.geom.Point2D(250+delta, 74))
+        # negative blended pair on diffim
         with dataset.addBlend() as family:
             family.addChild(instFlux=3E5, centroid=lsst.geom.Point2D(350, 72))
             family.addChild(instFlux=3.5E5, centroid=lsst.geom.Point2D(350+delta, 74))
+        # negative isolated source on diffim
         dataset.addSource(instFlux=4E5, centroid=lsst.geom.Point2D(75, 124))
+        # negative isolated source on diffim
         dataset.addSource(instFlux=5E5, centroid=lsst.geom.Point2D(175, 124))
         template, catalog = dataset.realize(noise=100.0,
                                             schema=lsst.meas.base.tests.TestDataset.makeMinimalSchema())
-        # These will be positive sources on the diffim.
+        # These will be positive sources on the diffim (as noted above).
         template.image.subset(lsst.geom.Box2I(lsst.geom.Point2I(15, 15),
                                               lsst.geom.Point2I(35, 35))).array *= -1
         template.image.subset(lsst.geom.Box2I(lsst.geom.Point2I(230, 60),
@@ -1046,6 +1067,11 @@ def testMergeFootprints(self):
         unique_peaks, counts = np.unique(peaks, axis=0, return_counts=True)
         self.assertEqual(np.sum(counts > 1), 0)
 
+        # There are three positive sources that should not have `is_negative`
+        # set, independent of how deblending/merging of negative sources is
+        # handled.
+        self.assertEqual((~result.diaSources["is_negative"]).sum(), 3)
+
 
 def setup_module(module):
     lsst.utils.tests.init()