Merge pull request #81 from lsst/tickets/DM-47535

DM-47535: Actually save calibrated catalog in ReprocessVisitImage.

python/lsst/drp/tasks/reprocess_visit_image.py

@@ -26,6 +26,7 @@
 
 import lsst.afw.image as afwImage
 import lsst.afw.table as afwTable
+import lsst.geom
 import lsst.meas.algorithms
 import lsst.meas.deblender
 import lsst.meas.extensions.photometryKron
@@ -238,32 +239,32 @@ class ReprocessVisitImageTask(pipeBase.PipelineTask):
     ConfigClass = ReprocessVisitImageConfig
     _DefaultName = "reprocessVisitImage"
 
-    def __init__(self, sources_schema=None, **kwargs):
+    def __init__(self, schema=None, **kwargs):
         super().__init__(**kwargs)
 
-        if sources_schema is None:
-            sources_schema = afwTable.SourceTable.makeMinimalSchema()
+        if schema is None:
+            schema = afwTable.SourceTable.makeMinimalSchema()
 
-        afwTable.CoordKey.addErrorFields(sources_schema)
+        afwTable.CoordKey.addErrorFields(schema)
         self.makeSubtask("snap_combine")
         self.makeSubtask("repair")
-        self.makeSubtask("detection", schema=sources_schema)
-        self.makeSubtask("sky_sources", schema=sources_schema)
-        self.makeSubtask("deblend", schema=sources_schema)
-        self.makeSubtask("measurement", schema=sources_schema)
-        self.makeSubtask("normalized_calibration_flux", schema=sources_schema)
-        self.makeSubtask("apply_aperture_correction", schema=sources_schema)
-        self.makeSubtask("catalog_calculation", schema=sources_schema)
-        self.makeSubtask("set_primary_flags", schema=sources_schema, isSingleFrame=True)
-        self.makeSubtask("post_calculations", schema=sources_schema)
+        self.makeSubtask("detection", schema=schema)
+        self.makeSubtask("sky_sources", schema=schema)
+        self.makeSubtask("deblend", schema=schema)
+        self.makeSubtask("measurement", schema=schema)
+        self.makeSubtask("normalized_calibration_flux", schema=schema)
+        self.makeSubtask("apply_aperture_correction", schema=schema)
+        self.makeSubtask("catalog_calculation", schema=schema)
+        self.makeSubtask("set_primary_flags", schema=schema, isSingleFrame=True)
+        self.makeSubtask("post_calculations", schema=schema)
         self.makeSubtask("compute_summary_stats")
 
-        sources_schema.addField(
+        schema.addField(
             "visit",
             type="I",
             doc="Visit this source appeared on.",
         )
-        sources_schema.addField(
+        schema.addField(
             "detector",
             type="U",
             doc="Detector this source appeared on.",
@@ -272,18 +273,18 @@ def __init__(self, sources_schema=None, **kwargs):
         # These fields will be propagated from finalizeCharacterization.
         # It might be better to get them from the finalized catalog instead
         # (if it output a schema), so the docstrings exactly match.
-        sources_schema.addField(
+        schema.addField(
             "calib_psf_candidate",
             type="Flag",
             doc="Set if the source was a candidate for PSF determination, "
             "as determined from FinalizeCharacterizationTask.",
         )
-        sources_schema.addField(
+        schema.addField(
             "calib_psf_reserved",
             type="Flag",
             doc="set if source was reserved from PSF determination by FinalizeCharacterizationTask.",
         )
-        sources_schema.addField(
+        schema.addField(
             "calib_psf_used",
             type="Flag",
             doc="Set if source was used in the PSF determination by FinalizeCharacterizationTask.",
@@ -294,28 +295,33 @@ def __init__(self, sources_schema=None, **kwargs):
         # These fields are only here to satisfy the SDM schema, and will
         # be removed from there as they are misleading (because we don't
         # propagate this information from gbdes/fgcmcal).
-        sources_schema.addField(
+        schema.addField(
             "calib_photometry_used",
             type="Flag",
             doc="Unused; placeholder for SDM schemas.",
         )
-        sources_schema.addField(
+        schema.addField(
             "calib_photometry_reserved",
             type="Flag",
             doc="Unused; placeholder for SDM schemas.",
         )
-        sources_schema.addField(
+        schema.addField(
             "calib_astrometry_used",
             type="Flag",
             doc="Unused; placeholder for SDM schemas.",
         )
-        sources_schema.addField(
+        schema.addField(
             "calib_astrometry_reserved",
             type="Flag",
             doc="Unused; placeholder for SDM schemas.",
         )
-
-        self.sources_schema = afwTable.SourceCatalog(sources_schema)
+        # This pre-calibration schema is the one that most methods should use.
+        self.schema = schema
+        # The final catalog will have calibrated flux columns, which we add to
+        # the init-output schema by calibrating our zero-length catalog with an
+        # arbitrary dummy PhotoCalib.
+        dummy_photo_calib = afwImage.PhotoCalib(1.0, 0, bbox=lsst.geom.Box2I())
+        self.sources_schema = dummy_photo_calib.calibrateCatalog(afwTable.SourceCatalog(schema))
 
     def runQuantum(self, butlerQC, inputRefs, outputRefs):
         inputs = butlerQC.get(inputRefs)
@@ -491,7 +497,11 @@ def run(
         result.exposure.info.setSummaryStats(
             self.compute_summary_stats.run(result.exposure, result.sources_footprints, background)
         )
-        self._apply_photo_calib(result.exposure, result.sources_footprints, photo_calib)
+        result.sources_footprints = self._apply_photo_calib(
+            result.exposure,
+            result.sources_footprints,
+            photo_calib,
+        )
         result.sources = result.sources_footprints.asAstropy()
 
         return result
@@ -516,7 +526,7 @@ def _find_sources(self, exposure, background, calib_sources, id_generator):
         sources
             Catalog that was detected and measured on the exposure.
         """
-        table = afwTable.SourceTable.make(self.sources_schema.schema, id_generator.make_table_id_factory())
+        table = afwTable.SourceTable.make(self.schema, id_generator.make_table_id_factory())
 
         self.repair.run(exposure=exposure)
         detections = self.detection.run(table=table, exposure=exposure, background=background)
@@ -592,11 +602,11 @@ def _apply_photo_calib(self, exposure, sources_footprints, photo_calib):
             Star catalog with flux/magnitude columns computed from the
             supplied PhotoCalib.
         """
-        sources_footprints = photo_calib.calibrateCatalog(sources_footprints)
+        calibrated_sources_footprints = photo_calib.calibrateCatalog(sources_footprints)
         exposure.maskedImage = photo_calib.calibrateImage(exposure.maskedImage)
         identity = afwImage.PhotoCalib(1.0, photo_calib.getCalibrationErr(), bbox=exposure.getBBox())
         exposure.setPhotoCalib(identity)
-        return sources_footprints
+        return calibrated_sources_footprints
 
 
 def combine_backgrounds(initial_pvi_background, sky_corr):

tests/test_reprocess_visit_image.py

@@ -71,7 +71,7 @@ def make_visit_summary(summary=None, psf=None, wcs=None, photo_calib=None, detec
     return summary
 
 
-class CalibrateImageTaskTests(lsst.utils.tests.TestCase):
+class ReprocessVisitImageTaskTests(lsst.utils.tests.TestCase):
     def setUp(self):
         # Different x/y dimensions so they're easy to distinguish in a plot,
         # and non-zero minimum, to help catch xy0 errors.
@@ -211,6 +211,13 @@ def test_run(self):
         # Faintest non-sky source should be marked as used.
         flux_sorted = result.sources[result.sources.argsort("slot_CalibFlux_instFlux")]
         self.assertTrue(flux_sorted[~flux_sorted["sky_source"]]["calib_psf_used"][0])
+        # Test that the schema init-output agrees with the catalog output.
+        self.assertEqual(task.sources_schema.schema, result.sources_footprints.schema)
+        # The flux/instFlux ratio should be the LocalPhotoCalib value.
+        for record in result.sources_footprints:
+            self.assertAlmostEqual(
+                record["base_PsfFlux_flux"] / record["base_PsfFlux_instFlux"], record["base_LocalPhotoCalib"]
+            )
 
 
 class ReprocessVisitImageTaskRunQuantumTests(lsst.utils.tests.TestCase):