Skip to content

Commit

Permalink
fixup reprocess
Browse files Browse the repository at this point in the history
got finalized catalog merging to work!
  • Loading branch information
parejkoj committed Sep 10, 2024
1 parent 38de6df commit 7769a5d
Show file tree
Hide file tree
Showing 3 changed files with 168 additions and 25 deletions.
89 changes: 78 additions & 11 deletions python/lsst/drp/tasks/reprocess_visit_image.py
Original file line number Diff line number Diff line change
Expand Up @@ -21,6 +21,8 @@

__all__ = ["ReprocessVisitImageTask", "ReprocessVisitImageConfig", "combine_backgrounds"]

import smatch

import lsst.afw.image as afwImage
import lsst.afw.table as afwTable
import lsst.meas.algorithms
Expand Down Expand Up @@ -49,6 +51,13 @@ class ReprocessVisitImageConnections(
storageClass="ExposureCatalog",
dimensions=["instrument", "visit"],
)
initial_photo_calib = connectionTypes.Input(
doc="Photometric calibration that was applied to exposure during the measurement of background_1."
" Used to uncalibrate the background before subtracting it from the input exposure.",
name="initial_photoCalib_detector",
storageClass="PhotoCalib",
dimensions=("instrument", "visit", "detector"),
)
background_1 = connectionTypes.Input(
doc="Background models estimated during calibration.",
name="initial_pvi_background",
Expand All @@ -61,9 +70,13 @@ class ReprocessVisitImageConnections(
dimensions=("instrument", "visit", "detector"),
storageClass="Background",
)

# TODO: we'll want the finalized catalog to get calib_psf_* from.
# TODO: eventually we'll want to pull in the STREAK mask from the diffim.
calib_sources = connectionTypes.Input(
doc="Per-visit catalog of measurements to get 'calib_*' flags from.",
name="finalized_src_table",
storageClass="ArrowAstropy",
dimensions=["instrument", "visit"],
)
# TODO DM-45980: pull in the STREAK mask from the diffim.

# outputs
source_schema = connectionTypes.InitOutput(
Expand Down Expand Up @@ -151,11 +164,15 @@ class ReprocessVisitImageConfig(
target=lsst.meas.base.CatalogCalculationTask,
doc="Task to compute catalog values using only the catalog entries.",
)
# TODO: hopefully we can get rid of this, once Jim sorts out SDM details.
post_calculations = pexConfig.ConfigurableField(
target=lsst.meas.base.SingleFrameMeasurementTask,
doc="Task to compute catalog values after all other calculations have been done.",
)
calib_match_radius = pexConfig.Field(
dtype=float,
default=0.2,
doc="Radius in arcseconds to match calib_sources to the output catalog."
)

def setDefaults(self):
super().setDefaults()
Expand Down Expand Up @@ -252,6 +269,7 @@ def __init__(self, source_schema=None, **kwargs):
type="Flag",
doc="Set if source was used in the PSF determination by FinalizeCharacterizationTask.",
)
self.psf_fields = ("calib_psf_candidate", "calib_psf_used", "calib_psf_reserved")

# These fields are only here to satisfy the SDM schema, and will
# be removed from there as they are misleading (because we don't
Expand Down Expand Up @@ -286,6 +304,8 @@ def runQuantum(self, butlerQC, inputRefs, outputRefs):
detector = outputRefs.exposure.dataId["detector"]
exposures = inputs.pop("exposures")
visit_summary = inputs.pop("visit_summary")
calib_sources = inputs.pop("calib_sources")
initial_photo_calib = inputs.pop("initial_photo_calib")
background_1 = inputs.pop("background_1")
if self.config.do_use_sky_corr:
background_2 = inputs.pop("background_2")
Expand All @@ -297,6 +317,10 @@ def runQuantum(self, butlerQC, inputRefs, outputRefs):
assert not inputs, "runQuantum got more inputs than expected"

detector_summary = visit_summary.find(detector)
if detector_summary is None:
msg = f"Detector {detector} not found in visit summary table; not reprocessing this exposure."
raise pipeBase.NoWorkFound(msg)

# Specify the fields that `annotate` needs below, to ensure they
# exist, even as None.
result = pipeBase.Struct(
Expand All @@ -306,12 +330,14 @@ def runQuantum(self, butlerQC, inputRefs, outputRefs):
try:
self.run(
exposures=exposures,
initial_photo_calib=initial_photo_calib,
psf=detector_summary.psf,
background=background,
ap_corr=detector_summary.apCorrMap,
photo_calib=detector_summary.photoCalib,
wcs=detector_summary.wcs,
summary_stats=afwImage.ExposureSummaryStats.from_record(detector_summary),
calib_sources=calib_sources,
result=result,
id_generator=id_generator,
)
Expand All @@ -328,12 +354,14 @@ def run(
self,
*,
exposures,
initial_photo_calib,
psf,
background,
ap_corr,
photo_calib,
wcs,
summary_stats,
calib_sources,
id_generator=None,
result=None
):
Expand All @@ -349,6 +377,9 @@ def run(
Modified in-place during processing if only one is passed.
If two exposures are passed, treat them as snaps and combine
before doing further processing.
initial_photo_calib : `lsst.afw.image.PhotoCalib`
Photometric calibration that was applied to exposure during the
measurement of the background.
psf : `lsst.afw.detection.Psf`
PSF model for this exposure.
background : `lsst.afw.math.BackgroundList`
Expand All @@ -362,6 +393,8 @@ def run(
World Coordinate System model for this exposure.
summary_stats : `lsst.afw.image.ExposureSummaryStats`
Summary statistics measured on this exposure.
calib_sources : `astropy.table.Table`
Per-visit catalog of measurements to get 'calib_*' flags from.
id_generator : `lsst.meas.base.IdGenerator`, optional
Object that generates source IDs and provides random seeds.
result : `lsst.pipe.base.Struct`, optional
Expand Down Expand Up @@ -393,11 +426,11 @@ def run(
result.exposure = self.snap_combine.run(exposures).exposure

# Calibrate the image, so it's on the same units as the background.
result.exposure.maskedImage = photo_calib.calibrateImage(result.exposure.maskedImage)
result.exposure.maskedImage = initial_photo_calib.calibrateImage(result.exposure.maskedImage)
result.exposure.maskedImage -= background.getImage()
# Uncalibrate so that we do the measurements in instFlux, because
# we don't have a way to identify measurements as being in nJy (TODO!).
result.exposure.maskedImage /= photo_calib.getCalibrationMean()
result.exposure.maskedImage /= initial_photo_calib.getCalibrationMean()

result.exposure.setPsf(psf)
result.exposure.setApCorrMap(ap_corr)
Expand All @@ -407,9 +440,10 @@ def run(
# to calibrate the image and catalog, and thus the image will have a
# photoCalib of exactly 1.

result.source_footprints = self._find_sources(result.exposure, background, id_generator)
result.source_footprints["visit"] = result.exposure.visitInfo.id
result.source_footprints["detector"] = result.exposure.info.getDetector().getId()
result.source_footprints = self._find_sources(result.exposure,
background,
calib_sources,
id_generator)
result.background = background
self._apply_photo_calib(result.exposure, result.source_footprints, photo_calib)

Expand All @@ -419,16 +453,18 @@ def run(

return result

def _find_sources(self, exposure, background, id_generator):
def _find_sources(self, exposure, background, calib_sources, id_generator):
"""Detect and measure sources on the exposure.
Parameters
----------
exposures : `lsst.afw.image.Exposure`
exposure : `lsst.afw.image.Exposure`
Exposure to detect and measure sources on; must have a valid PSF.
background : `lsst.afw.math.BackgroundList`
Background that was fit to the exposure during detection;
modified in-place during subsequent detection.
calib_sources : `astropy.table.Table`
Per-visit catalog of measurements to get 'calib_*' flags from.
id_generator : `lsst.meas.base.IdGenerator`
Object that generates source IDs and provides random seeds.
Expand Down Expand Up @@ -456,8 +492,39 @@ def _find_sources(self, exposure, background, id_generator):
self.catalog_calculation.run(source)
self.set_primary_flags.run(source)

source["visit"] = exposure.visitInfo.id
source["detector"] = exposure.info.getDetector().getId()

self._match_calib_sources(source, calib_sources, exposure.info.getDetector().getId())

return source

def _match_calib_sources(self, source, calib_sources, detector):
"""Match with calib_sources to set `calib_*` flags in the output
catalog.
Parameters
----------
source : `lsst.afw.table.SourceCatalog`
Catalog that was detected and measured on the exposure. Modified
in place to set the psf_fields.
calib_sources : `astropy.table.Table`
Per-visit catalog of measurements to get 'calib_*' flags from.
detector : `int`
Id of detector for this exposure, to get the correct sources from
calib_sources for cross-matching.
"""
use = calib_sources["detector"] == detector
sky = source["sky_source"]
with smatch.Matcher(source[~sky]["coord_ra"], source[~sky]["coord_dec"]) as matcher:
_, i1, i2, _ = matcher.query_knn(calib_sources[use]["coord_ra"],
calib_sources[use]["coord_dec"],
1,
self.config.calib_match_radius/3600.0,
return_indices=True)
for field in self.psf_fields:
source[~sky][field] = calib_sources[use][i1][field]

def _apply_photo_calib(self, exposure, source_footprints, photo_calib):
"""Photometrically calibrate the exposure and catalog with the
supplied PhotoCalib, and set the exposure's PhotoCalib to 1.
Expand Down
Loading

0 comments on commit 7769a5d

Please sign in to comment.