Merge pull request #1018 from lsst/tickets/DM-47738

DM-47738: Update for scarlet lite changes

Author: fred3m

python/lsst/pipe/tasks/deblendCoaddSourcesPipeline.py

@@ -134,6 +134,13 @@ class DeblendCoaddSourcesMultiConnections(PipelineTaskConnections,
         multiple=True,
         dimensions=("tract", "patch", "band", "skymap")
     )
+    deconvolvedCoadds = cT.Input(
+        doc="Deconvolved coadds",
+        name="deconvolved_{inputCoaddName}_coadd",
+        storageClass="ExposureF",
+        multiple=True,
+        dimensions=("tract", "patch", "band", "skymap")
+    )
     outputSchema = cT.InitOutput(
         doc="Output of the schema used in deblending task",
         name="{outputCoaddName}Coadd_deblendedFlux_schema",
@@ -169,6 +176,12 @@ class DeblendCoaddSourcesMultiConnections(PipelineTaskConnections,
         storageClass="ScarletModelData",
         dimensions=("tract", "patch", "skymap"),
     )
+    objectParents = cT.Output(
+        doc="Parents of the deblended objects",
+        name="object_parents",
+        storageClass="SourceCatalog",
+        dimensions=("tract", "patch", "skymap"),
+    )
 
     def __init__(self, *, config=None):
         super().__init__(config=config)
@@ -286,21 +299,33 @@ def runQuantum(self, butlerQC, inputRefs, outputRefs):
             coadds = exposures
         else:
             coadds = inputs.pop("coadds")
+
+        # Ensure that the coadd bands and deconvolved coadd bands match
+        deconvBands = [dRef.dataId["band"] for dRef in inputRefs.deconvolvedCoadds]
+        if bands != deconvBands:
+            self.log.error("Coadd bands %s != deconvolved coadd bands %s", bands, deconvBands)
+            raise RuntimeError("Number of coadd bands and deconvolved coadd bands do not match")
+
+        deconvolvedCoadds = inputs.pop("deconvolvedCoadds")
+
+        # Check that all inputs have been extracted correctly.
         assert not inputs, "runQuantum got extra inputs"
+
         outputs = self.run(
             coadds=coadds,
             bands=bands,
             mergedDetections=mergedDetections,
             idFactory=self.config.idGenerator.apply(butlerQC.quantum.dataId).make_table_id_factory(),
+            deconvolvedCoadds=deconvolvedCoadds,
         )
         butlerQC.put(outputs, outputRefs)
 
-    def run(self, coadds, bands, mergedDetections, idFactory):
+    def run(self, coadds, bands, mergedDetections, deconvolvedCoadds, idFactory):
         sources = self._makeSourceCatalog(mergedDetections, idFactory)
         multiExposure = afwImage.MultibandExposure.fromExposures(bands, coadds)
-        catalog, modelData = self.multibandDeblend.run(multiExposure, sources)
-        retStruct = Struct(deblendedCatalog=catalog, scarletModelData=modelData)
-        return retStruct
+        mDeconvolved = afwImage.MultibandExposure.fromExposures(bands, deconvolvedCoadds)
+        result = self.multibandDeblend.run(multiExposure, mDeconvolved, sources)
+        return result
 
     def _makeSourceCatalog(self, mergedDetections, idFactory):
         # There may be gaps in the mergeDet catalog, which will cause the

schemas/Object.yaml

@@ -1099,6 +1099,9 @@ refFlags:
     - deblend_peak_center_x
     - deblend_peak_center_y
     - deblend_logL
+    - deblend_chi2
+    - deblend_blendId
+    - deblend_blendNChild
     - slot_Shape_flag
 flags:
     # flags are columns taken and exploded per-band from the meas tables

schemas/initial_stars_detector_standardized.yaml

@@ -174,6 +174,10 @@ flags:
    - deblend_rampedTemplate
    - deblend_skipped
    - deblend_tooManyPeaks
+   - deblend_logL
+   - deblend_chi2
+   - deblend_blendId
+   - deblend_blendNChild
    - sky_source
    - detect_isPrimary
 
@@ -189,4 +193,3 @@ flag_rename_rules:
     - ['ext_shapeHSM_Hsm', 'hsm']
     - ['base_', '']
     - ['slot_', '']
-

tests/test_isPrimaryFlag.py

@@ -33,6 +33,7 @@
 from lsst.meas.algorithms import SourceDetectionTask, SkyObjectsTask, SetPrimaryFlagsTask
 import lsst.meas.extensions.scarlet as mes
 from lsst.meas.extensions.scarlet.scarletDeblendTask import ScarletDeblendTask
+from lsst.meas.extensions.scarlet.deconvolveExposureTask import DeconvolveExposureTask
 from lsst.meas.base import SingleFrameMeasurementTask
 from lsst.afw.table import SourceCatalog
 
@@ -214,6 +215,10 @@ def testIsScarletPrimaryFlag(self):
         skySourcesTask = SkyObjectsTask(name="skySources", config=skyConfig)
         schema.addField("merge_peak_sky", type="Flag")
 
+        # Initialize the deconvolution task
+        deconvolveConfig = DeconvolveExposureTask.ConfigClass()
+        deconvolveTask = DeconvolveExposureTask(config=deconvolveConfig)
+
         # Initialize the deblender task
         scarletConfig = ScarletDeblendTask.ConfigClass()
         scarletConfig.maxIter = 20
@@ -243,8 +248,18 @@ def testIsScarletPrimaryFlag(self):
             src = catalog.addNew()
             src.setFootprint(foot)
             src.set("merge_peak_sky", True)
+        # deconvolve the images
+        deconvolved = deconvolveTask.run(coadds["test"], catalog).deconvolved
+        mDeconvolved = afwImage.MultibandExposure.fromExposures(["test"], [deconvolved])
         # deblend
-        catalog, modelData = deblendTask.run(coadds, catalog)
+        # This is a hack because the variance is not calibrated properly
+        # (it is 3 orders of magnitude too high), which causes the deblender
+        # to improperly deblend most sources due to the sparsity constraint.
+        coadds.variance.array[:] = 2e-1
+        mDeconvolved.variance.array[:] = 2e-1
+        result = deblendTask.run(coadds, mDeconvolved, catalog)
+        modelData = result.scarletModelData
+        catalog = result.deblendedCatalog
         # Attach footprints to the catalog
         mes.io.updateCatalogFootprints(
             modelData=modelData,
@@ -254,6 +269,7 @@ def testIsScarletPrimaryFlag(self):
             removeScarletData=True,
             updateFluxColumns=True,
         )
+
         # measure
         measureTask.run(catalog, self.exposure)
         outputCat = catalog