Loading [MathJax]/extensions/tex2jax.js
LSST Applications g0fba68d861+83433b07ee,g16d25e1f1b+23bc9e47ac,g1ec0fe41b4+3ea9d11450,g1fd858c14a+9be2b0f3b9,g2440f9efcc+8c5ae1fdc5,g35bb328faa+8c5ae1fdc5,g4a4af6cd76+d25431c27e,g4d2262a081+c74e83464e,g53246c7159+8c5ae1fdc5,g55585698de+1e04e59700,g56a49b3a55+92a7603e7a,g60b5630c4e+1e04e59700,g67b6fd64d1+3fc8cb0b9e,g78460c75b0+7e33a9eb6d,g786e29fd12+668abc6043,g8352419a5c+8c5ae1fdc5,g8852436030+60e38ee5ff,g89139ef638+3fc8cb0b9e,g94187f82dc+1e04e59700,g989de1cb63+3fc8cb0b9e,g9d31334357+1e04e59700,g9f33ca652e+0a83e03614,gabe3b4be73+8856018cbb,gabf8522325+977d9fabaf,gb1101e3267+8b4b9c8ed7,gb89ab40317+3fc8cb0b9e,gc0af124501+57ccba3ad1,gcf25f946ba+60e38ee5ff,gd6cbbdb0b4+1cc2750d2e,gd794735e4e+7be992507c,gdb1c4ca869+be65c9c1d7,gde0f65d7ad+c7f52e58fe,ge278dab8ac+6b863515ed,ge410e46f29+3fc8cb0b9e,gf35d7ec915+97dd712d81,gf5e32f922b+8c5ae1fdc5,gf618743f1b+747388abfa,gf67bdafdda+3fc8cb0b9e,w.2025.18
LSST Data Management Base Package
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
calibrateImage.py
Go to the documentation of this file.
1# This file is part of pipe_tasks.
2#
3# Developed for the LSST Data Management System.
4# This product includes software developed by the LSST Project
5# (https://www.lsst.org).
6# See the COPYRIGHT file at the top-level directory of this distribution
7# for details of code ownership.
8#
9# This program is free software: you can redistribute it and/or modify
10# it under the terms of the GNU General Public License as published by
11# the Free Software Foundation, either version 3 of the License, or
12# (at your option) any later version.
13#
14# This program is distributed in the hope that it will be useful,
15# but WITHOUT ANY WARRANTY; without even the implied warranty of
16# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17# GNU General Public License for more details.
18#
19# You should have received a copy of the GNU General Public License
20# along with this program. If not, see <https://www.gnu.org/licenses/>.
21
22__all__ = ["CalibrateImageTask", "CalibrateImageConfig", "NoPsfStarsToStarsMatchError",
23 "AllCentroidsFlaggedError"]
24
25import numpy as np
26
27import lsst.afw.table as afwTable
28import lsst.afw.image as afwImage
29from lsst.ip.diffim.utils import evaluateMaskFraction
34import lsst.meas.base
37import lsst.meas.extensions.shapeHSM
38import lsst.pex.config as pexConfig
39import lsst.pipe.base as pipeBase
40from lsst.pipe.base import connectionTypes
41from lsst.utils.timer import timeMethod
42
43from . import measurePsf, repair, photoCal, computeExposureSummaryStats, snapCombine
44
45
46class AllCentroidsFlaggedError(pipeBase.AlgorithmError):
47 """Raised when there are no valid centroids during psf fitting.
48 """
49 def __init__(self, n_sources, psf_shape_ixx, psf_shape_iyy, psf_shape_ixy, psf_size):
50 msg = (f"All source centroids (out of {n_sources}) flagged during PSF fitting. "
51 "Original image PSF is likely unuseable; best-fit PSF shape parameters: "
52 f"Ixx={psf_shape_ixx}, Iyy={psf_shape_iyy}, Ixy={psf_shape_ixy}, size={psf_size}"
53 )
54 super().__init__(msg)
55 self.n_sources = n_sources
56 self.psf_shape_ixx = psf_shape_ixx
57 self.psf_shape_iyy = psf_shape_iyy
58 self.psf_shape_ixy = psf_shape_ixy
59 self.psf_size = psf_size
60
61 @property
62 def metadata(self):
63 return {"n_sources": self.n_sources,
64 "psf_shape_ixx": self.psf_shape_ixx,
65 "psf_shape_iyy": self.psf_shape_iyy,
66 "psf_shape_ixy": self.psf_shape_ixy,
67 "psf_size": self.psf_size
68 }
69
70
71class NoPsfStarsToStarsMatchError(pipeBase.AlgorithmError):
72 """Raised when there are no matches between the psf_stars and stars
73 catalogs.
74 """
75 def __init__(self, *, n_psf_stars, n_stars):
76 msg = (f"No psf stars out of {n_psf_stars} matched {n_stars} calib stars."
77 " Downstream processes probably won't have useful {calib_type} stars in this case."
78 " Is `star_selector` too strict or is this a bad image?")
79 super().__init__(msg)
80 self.n_psf_stars = n_psf_stars
81 self.n_stars = n_stars
82
83 @property
84 def metadata(self):
85 return {"n_psf_stars": self.n_psf_stars,
86 "n_stars": self.n_stars
87 }
88
89
90class CalibrateImageConnections(pipeBase.PipelineTaskConnections,
91 dimensions=("instrument", "visit", "detector")):
92
93 astrometry_ref_cat = connectionTypes.PrerequisiteInput(
94 doc="Reference catalog to use for astrometric calibration.",
95 name="gaia_dr3_20230707",
96 storageClass="SimpleCatalog",
97 dimensions=("skypix",),
98 deferLoad=True,
99 multiple=True,
100 )
101 photometry_ref_cat = connectionTypes.PrerequisiteInput(
102 doc="Reference catalog to use for photometric calibration.",
103 name="ps1_pv3_3pi_20170110",
104 storageClass="SimpleCatalog",
105 dimensions=("skypix",),
106 deferLoad=True,
107 multiple=True
108 )
109
110 exposures = connectionTypes.Input(
111 doc="Exposure (or two snaps) to be calibrated, and detected and measured on.",
112 name="postISRCCD",
113 storageClass="Exposure",
114 multiple=True, # to handle 1 exposure or 2 snaps
115 dimensions=["instrument", "exposure", "detector"],
116 )
117
118 background_flat = connectionTypes.PrerequisiteInput(
119 name="flat",
120 doc="Flat calibration frame used for background correction.",
121 storageClass="ExposureF",
122 dimensions=["instrument", "detector", "physical_filter"],
123 isCalibration=True,
124 )
125 illumination_correction = connectionTypes.PrerequisiteInput(
126 name="illuminationCorrection",
127 doc="Illumination correction frame.",
128 storageClass="Exposure",
129 dimensions=["instrument", "detector", "physical_filter"],
130 isCalibration=True,
131 )
132
133 # outputs
134 initial_stars_schema = connectionTypes.InitOutput(
135 doc="Schema of the output initial stars catalog.",
136 name="initial_stars_schema",
137 storageClass="SourceCatalog",
138 )
139
140 # TODO DM-38732: We want some kind of flag on Exposures/Catalogs to make
141 # it obvious which components had failed to be computed/persisted.
142 exposure = connectionTypes.Output(
143 doc="Photometrically calibrated, background-subtracted exposure with fitted calibrations and "
144 "summary statistics. To recover the original exposure, first add the background "
145 "(`initial_pvi_background`), and then uncalibrate (divide by `initial_photoCalib_detector`).",
146 name="initial_pvi",
147 storageClass="ExposureF",
148 dimensions=("instrument", "visit", "detector"),
149 )
150 stars = connectionTypes.Output(
151 doc="Catalog of unresolved sources detected on the calibrated exposure.",
152 name="initial_stars_detector",
153 storageClass="ArrowAstropy",
154 dimensions=["instrument", "visit", "detector"],
155 )
156 stars_footprints = connectionTypes.Output(
157 doc="Catalog of unresolved sources detected on the calibrated exposure; "
158 "includes source footprints.",
159 name="initial_stars_footprints_detector",
160 storageClass="SourceCatalog",
161 dimensions=["instrument", "visit", "detector"],
162 )
163 applied_photo_calib = connectionTypes.Output(
164 doc=(
165 "Photometric calibration that was applied to exposure's pixels. "
166 "This connection is disabled when do_calibrate_pixels=False."
167 ),
168 name="initial_photoCalib_detector",
169 storageClass="PhotoCalib",
170 dimensions=("instrument", "visit", "detector"),
171 )
172 background = connectionTypes.Output(
173 doc="Background models estimated during calibration task; calibrated to be in nJy units.",
174 name="initial_pvi_background",
175 storageClass="Background",
176 dimensions=("instrument", "visit", "detector"),
177 )
178 background_to_photometric_ratio = connectionTypes.Output(
179 doc="Ratio of a background-flattened image to a photometric-flattened image. Only persisted "
180 "if do_illumination_correction is True.",
181 name="background_to_photometric_ratio",
182 storageClass="Image",
183 dimensions=("instrument", "visit", "detector"),
184 )
185
186 # Optional outputs
187 psf_stars_footprints = connectionTypes.Output(
188 doc="Catalog of bright unresolved sources detected on the exposure used for PSF determination; "
189 "includes source footprints.",
190 name="initial_psf_stars_footprints_detector",
191 storageClass="SourceCatalog",
192 dimensions=["instrument", "visit", "detector"],
193 )
194 psf_stars = connectionTypes.Output(
195 doc="Catalog of bright unresolved sources detected on the exposure used for PSF determination.",
196 name="initial_psf_stars_detector",
197 storageClass="ArrowAstropy",
198 dimensions=["instrument", "visit", "detector"],
199 )
200 astrometry_matches = connectionTypes.Output(
201 doc="Source to reference catalog matches from the astrometry solver.",
202 name="initial_astrometry_match_detector",
203 storageClass="Catalog",
204 dimensions=("instrument", "visit", "detector"),
205 )
206 photometry_matches = connectionTypes.Output(
207 doc="Source to reference catalog matches from the photometry solver.",
208 name="initial_photometry_match_detector",
209 storageClass="Catalog",
210 dimensions=("instrument", "visit", "detector"),
211 )
212
213 def __init__(self, *, config=None):
214 super().__init__(config=config)
215 if "psf_stars" not in config.optional_outputs:
216 del self.psf_stars
217 if "psf_stars_footprints" not in config.optional_outputs:
218 del self.psf_stars_footprints
219 if "astrometry_matches" not in config.optional_outputs:
220 del self.astrometry_matches
221 if "photometry_matches" not in config.optional_outputs:
222 del self.photometry_matches
223 if not config.do_calibrate_pixels:
224 del self.applied_photo_calib
225 if not config.do_illumination_correction:
226 del self.background_flat
229
230
231class CalibrateImageConfig(pipeBase.PipelineTaskConfig, pipelineConnections=CalibrateImageConnections):
232 optional_outputs = pexConfig.ListField(
233 doc="Which optional outputs to save (as their connection name)?",
234 dtype=str,
235 # TODO: note somewhere to disable this for benchmarking, but should
236 # we always have it on for production runs?
237 default=["psf_stars", "psf_stars_footprints", "astrometry_matches", "photometry_matches"],
238 optional=False
239 )
240
241 # To generate catalog ids consistently across subtasks.
242 id_generator = lsst.meas.base.DetectorVisitIdGeneratorConfig.make_field()
243
244 snap_combine = pexConfig.ConfigurableField(
246 doc="Task to combine two snaps to make one exposure.",
247 )
248
249 # subtasks used during psf characterization
250 install_simple_psf = pexConfig.ConfigurableField(
252 doc="Task to install a simple PSF model into the input exposure to use "
253 "when detecting bright sources for PSF estimation.",
254 )
255 psf_repair = pexConfig.ConfigurableField(
256 target=repair.RepairTask,
257 doc="Task to repair cosmic rays on the exposure before PSF determination.",
258 )
259 psf_subtract_background = pexConfig.ConfigurableField(
261 doc="Task to perform intial background subtraction, before first detection pass.",
262 )
263 psf_detection = pexConfig.ConfigurableField(
265 doc="Task to detect sources for PSF determination."
266 )
267 psf_source_measurement = pexConfig.ConfigurableField(
269 doc="Task to measure sources to be used for psf estimation."
270 )
271 psf_measure_psf = pexConfig.ConfigurableField(
273 doc="Task to measure the psf on bright sources."
274 )
275 psf_normalized_calibration_flux = pexConfig.ConfigurableField(
277 doc="Task to normalize the calibration flux (e.g. compensated tophats) "
278 "for the bright stars used for psf estimation.",
279 )
280
281 # TODO DM-39203: we can remove aperture correction from this task once we are
282 # using the shape-based star/galaxy code.
283 measure_aperture_correction = pexConfig.ConfigurableField(
285 doc="Task to compute the aperture correction from the bright stars."
286 )
287
288 # subtasks used during star measurement
289 star_detection = pexConfig.ConfigurableField(
291 doc="Task to detect stars to return in the output catalog."
292 )
293 star_sky_sources = pexConfig.ConfigurableField(
295 doc="Task to generate sky sources ('empty' regions where there are no detections).",
296 )
297 star_deblend = pexConfig.ConfigurableField(
299 doc="Split blended sources into their components."
300 )
301 star_measurement = pexConfig.ConfigurableField(
303 doc="Task to measure stars to return in the output catalog."
304 )
305 star_normalized_calibration_flux = pexConfig.ConfigurableField(
307 doc="Task to apply the normalization for calibration fluxes (e.g. compensated tophats) "
308 "for the final output star catalog.",
309 )
310 star_apply_aperture_correction = pexConfig.ConfigurableField(
312 doc="Task to apply aperture corrections to the selected stars."
313 )
314 star_catalog_calculation = pexConfig.ConfigurableField(
316 doc="Task to compute extendedness values on the star catalog, "
317 "for the star selector to remove extended sources."
318 )
319 star_set_primary_flags = pexConfig.ConfigurableField(
321 doc="Task to add isPrimary to the catalog."
322 )
323 star_selector = lsst.meas.algorithms.sourceSelectorRegistry.makeField(
324 default="science",
325 doc="Task to select reliable stars to use for calibration."
326 )
327
328 # final calibrations and statistics
329 astrometry = pexConfig.ConfigurableField(
331 doc="Task to perform astrometric calibration to fit a WCS.",
332 )
333 astrometry_ref_loader = pexConfig.ConfigField(
335 doc="Configuration of reference object loader for astrometric fit.",
336 )
337 photometry = pexConfig.ConfigurableField(
339 doc="Task to perform photometric calibration to fit a PhotoCalib.",
340 )
341 photometry_ref_loader = pexConfig.ConfigField(
343 doc="Configuration of reference object loader for photometric fit.",
344 )
345
346 compute_summary_stats = pexConfig.ConfigurableField(
348 doc="Task to to compute summary statistics on the calibrated exposure."
349 )
350
351 do_illumination_correction = pexConfig.Field(
352 dtype=bool,
353 default=False,
354 doc="If True, apply the illumination correction. This assumes that the "
355 "input image has already been flat-fielded such that it is suitable "
356 "for background subtraction.",
357 )
358
359 do_calibrate_pixels = pexConfig.Field(
360 dtype=bool,
361 default=True,
362 doc=(
363 "If True, apply the photometric calibration to the image pixels "
364 "and background model, and attach an identity PhotoCalib to "
365 "the output image to reflect this. If False`, leave the image "
366 "and background uncalibrated and attach the PhotoCalib that maps "
367 "them to physical units."
368 )
369 )
370
371 def setDefaults(self):
372 super().setDefaults()
373
374 # Use a very broad PSF here, to throughly reject CRs.
375 # TODO investigation: a large initial psf guess may make stars look
376 # like CRs for very good seeing images.
377 self.install_simple_psf.fwhm = 4
378
379 # S/N>=50 sources for PSF determination, but detection to S/N=10.
380 # The thresholdValue sets the minimum flux in a pixel to be included in the
381 # footprint, while peaks are only detected when they are above
382 # thresholdValue * includeThresholdMultiplier. The low thresholdValue
383 # ensures that the footprints are large enough for the noise replacer
384 # to mask out faint undetected neighbors that are not to be measured.
385 self.psf_detection.thresholdValue = 10.0
386 self.psf_detection.includeThresholdMultiplier = 5.0
387 # TODO investigation: Probably want False here, but that may require
388 # tweaking the background spatial scale, to make it small enough to
389 # prevent extra peaks in the wings of bright objects.
390 self.psf_detection.doTempLocalBackground = False
391 # NOTE: we do want reEstimateBackground=True in psf_detection, so that
392 # each measurement step is done with the best background available.
393
394 # Minimal measurement plugins for PSF determination.
395 # TODO DM-39203: We can drop GaussianFlux and PsfFlux, if we use
396 # shapeHSM/moments for star/galaxy separation.
397 # TODO DM-39203: we can remove aperture correction from this task once
398 # we are using the shape-based star/galaxy code.
399 self.psf_source_measurement.plugins = ["base_PixelFlags",
400 "base_SdssCentroid",
401 "ext_shapeHSM_HsmSourceMoments",
402 "base_CircularApertureFlux",
403 "base_GaussianFlux",
404 "base_PsfFlux",
405 "base_CompensatedTophatFlux",
406 ]
407 self.psf_source_measurement.slots.shape = "ext_shapeHSM_HsmSourceMoments"
408 # Only measure apertures we need for PSF measurement.
409 self.psf_source_measurement.plugins["base_CircularApertureFlux"].radii = [12.0]
410 self.psf_source_measurement.plugins["base_CompensatedTophatFlux"].apertures = [12]
411 # TODO DM-40843: Remove this line once this is the psfex default.
412 self.psf_measure_psf.psfDeterminer["psfex"].photometricFluxField = \
413 "base_CircularApertureFlux_12_0_instFlux"
414
415 # No extendeness information available: we need the aperture
416 # corrections to determine that.
417 self.measure_aperture_correction.sourceSelector["science"].doUnresolved = False
418 self.measure_aperture_correction.sourceSelector["science"].flags.good = ["calib_psf_used"]
419 self.measure_aperture_correction.sourceSelector["science"].flags.bad = []
420
421 # Detection for good S/N for astrometry/photometry and other
422 # downstream tasks; detection mask to S/N>=5, but S/N>=10 peaks.
423 self.star_detection.thresholdValue = 5.0
424 self.star_detection.includeThresholdMultiplier = 2.0
425 self.star_measurement.plugins = ["base_PixelFlags",
426 "base_SdssCentroid",
427 "ext_shapeHSM_HsmSourceMoments",
428 "ext_shapeHSM_HsmPsfMoments",
429 "base_GaussianFlux",
430 "base_PsfFlux",
431 "base_CircularApertureFlux",
432 "base_ClassificationSizeExtendedness",
433 "base_CompensatedTophatFlux",
434 ]
435 self.star_measurement.slots.psfShape = "ext_shapeHSM_HsmPsfMoments"
436 self.star_measurement.slots.shape = "ext_shapeHSM_HsmSourceMoments"
437 # Only measure the apertures we need for star selection.
438 self.star_measurement.plugins["base_CircularApertureFlux"].radii = [12.0]
439 self.star_measurement.plugins["base_CompensatedTophatFlux"].apertures = [12]
440
441 # We measure and apply the normalization aperture correction with the
442 # psf_normalized_calibration_flux task, and we only apply the normalization
443 # aperture correction for the full list of stars.
444 self.star_normalized_calibration_flux.do_measure_ap_corr = False
445
446 # Select stars with reliable measurements and no bad flags.
447 self.star_selector["science"].doFlags = True
448 self.star_selector["science"].doUnresolved = True
449 self.star_selector["science"].doSignalToNoise = True
450 self.star_selector["science"].signalToNoise.minimum = 10.0
451 # Keep sky sources in the output catalog, even though they aren't
452 # wanted for calibration.
453 self.star_selector["science"].doSkySources = True
454 # Set the flux and error fields
455 self.star_selector["science"].signalToNoise.fluxField = "slot_CalibFlux_instFlux"
456 self.star_selector["science"].signalToNoise.errField = "slot_CalibFlux_instFluxErr"
457
458 # Use the affine WCS fitter (assumes we have a good camera geometry).
459 self.astrometry.wcsFitter.retarget(lsst.meas.astrom.FitAffineWcsTask)
460 # phot_g_mean is the primary Gaia band for all input bands.
461 self.astrometry_ref_loader.anyFilterMapsToThis = "phot_g_mean"
462
463 # Only reject sky sources; we already selected good stars.
464 self.astrometry.sourceSelector["science"].doFlags = True
465 self.astrometry.sourceSelector["science"].flags.good = ["calib_psf_candidate"]
466 self.astrometry.sourceSelector["science"].flags.bad = []
467 self.astrometry.sourceSelector["science"].doUnresolved = False
468 self.astrometry.sourceSelector["science"].doIsolated = False
469 self.astrometry.sourceSelector["science"].doRequirePrimary = False
470 self.photometry.match.sourceSelection.doFlags = True
471 self.photometry.match.sourceSelection.flags.bad = ["sky_source"]
472 # Unset the (otherwise reasonable, but we've already made the
473 # selections we want above) selection settings in PhotoCalTask.
474 self.photometry.match.sourceSelection.doRequirePrimary = False
475 self.photometry.match.sourceSelection.doUnresolved = False
476
477 def validate(self):
478 super().validate()
479
480 # Ensure that the normalization calibration flux tasks
481 # are configured correctly.
482 if not self.psf_normalized_calibration_flux.do_measure_ap_corr:
483 msg = ("psf_normalized_calibration_flux task must be configured with do_measure_ap_corr=True "
484 "or else the normalization and calibration flux will not be properly measured.")
485 raise pexConfig.FieldValidationError(
486 CalibrateImageConfig.psf_normalized_calibration_flux, self, msg,
487 )
488 if self.star_normalized_calibration_flux.do_measure_ap_corr:
489 msg = ("star_normalized_calibration_flux task must be configured with do_measure_ap_corr=False "
490 "to apply the previously measured normalization to the full catalog of calibration "
491 "fluxes.")
492 raise pexConfig.FieldValidationError(
493 CalibrateImageConfig.star_normalized_calibration_flux, self, msg,
494 )
495
496 # Ensure base_LocalPhotoCalib and base_LocalWcs plugins are not run,
497 # because they'd be running too early to pick up the fitted PhotoCalib
498 # and WCS.
499 if "base_LocalWcs" in self.psf_source_measurement.plugins.names:
500 raise pexConfig.FieldValidationError(
501 CalibrateImageConfig.psf_source_measurement,
502 self,
503 "base_LocalWcs cannot run CalibrateImageTask, as it would be run before the astrometry fit."
504 )
505 if "base_LocalWcs" in self.star_measurement.plugins.names:
506 raise pexConfig.FieldValidationError(
507 CalibrateImageConfig.star_measurement,
508 self,
509 "base_LocalWcs cannot run CalibrateImageTask, as it would be run before the astrometry fit."
510 )
511 if "base_LocalPhotoCalib" in self.psf_source_measurement.plugins.names:
512 raise pexConfig.FieldValidationError(
513 CalibrateImageConfig.psf_source_measurement,
514 self,
515 "base_LocalPhotoCalib cannot run CalibrateImageTask, "
516 "as it would be run before the photometry fit."
517 )
518 if "base_LocalPhotoCalib" in self.star_measurement.plugins.names:
519 raise pexConfig.FieldValidationError(
520 CalibrateImageConfig.star_measurement,
521 self,
522 "base_LocalPhotoCalib cannot run CalibrateImageTask, "
523 "as it would be run before the photometry fit."
524 )
525
526 # Check for illumination correction and background consistency.
528 if not self.psf_subtract_background.doApplyFlatBackgroundRatio:
529 raise pexConfig.FieldValidationError(
530 CalibrateImageConfig.psf_subtract_background,
531 self,
532 "CalibrateImageTask.psf_subtract_background must be configured with "
533 "doApplyFlatBackgroundRatio=True if do_illumination_correction=True."
534 )
535 if self.psf_detection.reEstimateBackground:
536 if not self.psf_detection.doApplyFlatBackgroundRatio:
537 raise pexConfig.FieldValidationError(
538 CalibrateImageConfig.psf_detection,
539 self,
540 "CalibrateImageTask.psf_detection background must be configured with "
541 "doApplyFlatBackgroundRatio=True if do_illumination_correction=True."
542 )
543 if self.star_detection.reEstimateBackground:
544 if not self.star_detection.doApplyFlatBackgroundRatio:
545 raise pexConfig.FieldValidationError(
546 CalibrateImageConfig.star_detection,
547 self,
548 "CalibrateImageTask.star_detection background must be configured with "
549 "doApplyFlatBackgroundRatio=True if do_illumination_correction=True."
550 )
551
552
553class CalibrateImageTask(pipeBase.PipelineTask):
554 """Compute the PSF, aperture corrections, astrometric and photometric
555 calibrations, and summary statistics for a single science exposure, and
556 produce a catalog of brighter stars that were used to calibrate it.
557
558 Parameters
559 ----------
560 initial_stars_schema : `lsst.afw.table.Schema`
561 Schema of the initial_stars output catalog.
562 """
563 _DefaultName = "calibrateImage"
564 ConfigClass = CalibrateImageConfig
565
566 def __init__(self, initial_stars_schema=None, **kwargs):
567 super().__init__(**kwargs)
568 self.makeSubtask("snap_combine")
569
570 # PSF determination subtasks
571 self.makeSubtask("install_simple_psf")
572 self.makeSubtask("psf_repair")
573 self.makeSubtask("psf_subtract_background")
574 self.psf_schema = afwTable.SourceTable.makeMinimalSchema()
575 afwTable.CoordKey.addErrorFields(self.psf_schema)
576 self.makeSubtask("psf_detection", schema=self.psf_schema)
577 self.makeSubtask("psf_source_measurement", schema=self.psf_schema)
578 self.makeSubtask("psf_measure_psf", schema=self.psf_schema)
579 self.makeSubtask("psf_normalized_calibration_flux", schema=self.psf_schema)
580
581 self.makeSubtask("measure_aperture_correction", schema=self.psf_schema)
582 self.makeSubtask("astrometry", schema=self.psf_schema)
583
584 # star measurement subtasks
585 if initial_stars_schema is None:
586 initial_stars_schema = afwTable.SourceTable.makeMinimalSchema()
587
588 # These fields let us track which sources were used for psf modeling,
589 # astrometric fitting, and aperture correction calculations.
590 self.psf_fields = ("calib_psf_candidate", "calib_psf_used", "calib_psf_reserved",
591 "calib_astrometry_used",
592 # TODO DM-39203: these can be removed once apcorr is gone.
593 "apcorr_slot_CalibFlux_used", "apcorr_base_GaussianFlux_used",
594 "apcorr_base_PsfFlux_used")
595 for field in self.psf_fields:
596 item = self.psf_schema.find(field)
597 initial_stars_schema.addField(item.getField())
598 id_type = self.psf_schema["id"].asField().getTypeString()
599 initial_stars_schema.addField("psf_id",
600 type=id_type,
601 doc="id of this source in psf_stars; 0 if there is no match.")
602
603 afwTable.CoordKey.addErrorFields(initial_stars_schema)
604 self.makeSubtask("star_detection", schema=initial_stars_schema)
605 self.makeSubtask("star_sky_sources", schema=initial_stars_schema)
606 self.makeSubtask("star_deblend", schema=initial_stars_schema)
607 self.makeSubtask("star_measurement", schema=initial_stars_schema)
608 self.makeSubtask("star_normalized_calibration_flux", schema=initial_stars_schema)
609
610 self.makeSubtask("star_apply_aperture_correction", schema=initial_stars_schema)
611 self.makeSubtask("star_catalog_calculation", schema=initial_stars_schema)
612 self.makeSubtask("star_set_primary_flags", schema=initial_stars_schema, isSingleFrame=True)
613 self.makeSubtask("star_selector")
614 self.makeSubtask("photometry", schema=initial_stars_schema)
615 self.makeSubtask("compute_summary_stats")
616
617 # The final catalog will have calibrated flux columns, which we add to
618 # the init-output schema by calibrating our zero-length catalog with an
619 # arbitrary dummy PhotoCalib. We also use this schema to initialze
620 # the stars catalog in order to ensure it's the same even when we hit
621 # an error (and write partial outputs) before calibrating the catalog
622 # - note that calibrateCatalog will happily reuse existing output
623 # columns.
624 dummy_photo_calib = afwImage.PhotoCalib(1.0, 0, bbox=lsst.geom.Box2I())
625 self.initial_stars_schema = dummy_photo_calib.calibrateCatalog(
626 afwTable.SourceCatalog(initial_stars_schema)
627 )
628
629 def runQuantum(self, butlerQC, inputRefs, outputRefs):
630 inputs = butlerQC.get(inputRefs)
631 exposures = inputs.pop("exposures")
632
633 id_generator = self.config.id_generator.apply(butlerQC.quantum.dataId)
634
636 dataIds=[ref.datasetRef.dataId for ref in inputRefs.astrometry_ref_cat],
637 refCats=inputs.pop("astrometry_ref_cat"),
638 name=self.config.connections.astrometry_ref_cat,
639 config=self.config.astrometry_ref_loader, log=self.log)
640 self.astrometry.setRefObjLoader(astrometry_loader)
641
643 dataIds=[ref.datasetRef.dataId for ref in inputRefs.photometry_ref_cat],
644 refCats=inputs.pop("photometry_ref_cat"),
645 name=self.config.connections.photometry_ref_cat,
646 config=self.config.photometry_ref_loader, log=self.log)
647 self.photometry.match.setRefObjLoader(photometry_loader)
648
649 if self.config.do_illumination_correction:
650 background_flat = inputs.pop("background_flat")
651 illumination_correction = inputs.pop("illumination_correction")
652 else:
653 background_flat = None
654 illumination_correction = None
655
656 # This should not happen with a properly configured execution context.
657 assert not inputs, "runQuantum got more inputs than expected"
658
659 # Specify the fields that `annotate` needs below, to ensure they
660 # exist, even as None.
661 result = pipeBase.Struct(
662 exposure=None,
663 stars_footprints=None,
664 psf_stars_footprints=None,
665 background_to_photometric_ratio=None,
666 )
667 try:
668 self.run(
669 exposures=exposures,
670 result=result,
671 id_generator=id_generator,
672 background_flat=background_flat,
673 illumination_correction=illumination_correction,
674 )
675 except pipeBase.AlgorithmError as e:
676 error = pipeBase.AnnotatedPartialOutputsError.annotate(
677 e,
678 self,
679 result.exposure,
680 result.psf_stars_footprints,
681 result.stars_footprints,
682 log=self.log
683 )
684 butlerQC.put(result, outputRefs)
685 raise error from e
686
687 butlerQC.put(result, outputRefs)
688
689 @timeMethod
690 def run(
691 self,
692 *,
693 exposures,
694 id_generator=None,
695 result=None,
696 background_flat=None,
697 illumination_correction=None,
698 ):
699 """Find stars and perform psf measurement, then do a deeper detection
700 and measurement and calibrate astrometry and photometry from that.
701
702 Parameters
703 ----------
704 exposures : `lsst.afw.image.Exposure` or `list` [`lsst.afw.image.Exposure`]
705 Post-ISR exposure(s), with an initial WCS, VisitInfo, and Filter.
706 Modified in-place during processing if only one is passed.
707 If two exposures are passed, treat them as snaps and combine
708 before doing further processing.
709 id_generator : `lsst.meas.base.IdGenerator`, optional
710 Object that generates source IDs and provides random seeds.
711 result : `lsst.pipe.base.Struct`, optional
712 Result struct that is modified to allow saving of partial outputs
713 for some failure conditions. If the task completes successfully,
714 this is also returned.
715 background_flat : `lsst.afw.image.Exposure`, optional
716 Background flat-field image.
717 illumination_correction : `lsst.afw.image.Exposure`, optional
718 Illumination correction image.
719
720 Returns
721 -------
722 result : `lsst.pipe.base.Struct`
723 Results as a struct with attributes:
724
725 ``exposure``
726 Calibrated exposure, with pixels in nJy units.
727 (`lsst.afw.image.Exposure`)
728 ``stars``
729 Stars that were used to calibrate the exposure, with
730 calibrated fluxes and magnitudes.
731 (`astropy.table.Table`)
732 ``stars_footprints``
733 Footprints of stars that were used to calibrate the exposure.
734 (`lsst.afw.table.SourceCatalog`)
735 ``psf_stars``
736 Stars that were used to determine the image PSF.
737 (`astropy.table.Table`)
738 ``psf_stars_footprints``
739 Footprints of stars that were used to determine the image PSF.
740 (`lsst.afw.table.SourceCatalog`)
741 ``background``
742 Background that was fit to the exposure when detecting
743 ``stars``. (`lsst.afw.math.BackgroundList`)
744 ``applied_photo_calib``
745 Photometric calibration that was fit to the star catalog and
746 applied to the exposure. (`lsst.afw.image.PhotoCalib`)
747 This is `None` if ``config.do_calibrate_pixels`` is `False`.
748 ``astrometry_matches``
749 Reference catalog stars matches used in the astrometric fit.
750 (`list` [`lsst.afw.table.ReferenceMatch`] or `lsst.afw.table.BaseCatalog`)
751 ``photometry_matches``
752 Reference catalog stars matches used in the photometric fit.
753 (`list` [`lsst.afw.table.ReferenceMatch`] or `lsst.afw.table.BaseCatalog`)
754 """
755 if result is None:
756 result = pipeBase.Struct()
757 if id_generator is None:
758 id_generator = lsst.meas.base.IdGenerator()
759
760 result.exposure = self.snap_combine.run(exposures).exposure
761 self._recordMaskedPixelFractions(result.exposure)
762 self.log.info("Initial PhotoCalib: %s", result.exposure.getPhotoCalib())
763
764 result.exposure.metadata["LSST CALIB ILLUMCORR APPLIED"] = False
765
766 # Check input image processing.
767 if self.config.do_illumination_correction:
768 if not result.exposure.metadata.get("LSST ISR FLAT APPLIED", False):
769 raise pipeBase.InvalidQuantumError(
770 "Cannot use do_illumination_correction with an image that has not had a flat applied",
771 )
772
773 result.background = None
774 summary_stat_catalog = None
775 # Some exposure components are set to initial placeholder objects
776 # while we try to bootstrap them. If we fail before we fit for them,
777 # we want to reset those components to None so the placeholders don't
778 # masquerade as the real thing.
779 have_fit_psf = False
780 have_fit_astrometry = False
781 have_fit_photometry = False
782 try:
783 result.background_to_photometric_ratio = self._apply_illumination_correction(
784 result.exposure,
785 background_flat,
786 illumination_correction,
787 )
788
789 result.psf_stars_footprints, result.background, _ = self._compute_psf(
790 result.exposure,
791 id_generator,
792 background_to_photometric_ratio=result.background_to_photometric_ratio,
793 )
794 have_fit_psf = True
795 self._measure_aperture_correction(result.exposure, result.psf_stars_footprints)
796 result.psf_stars = result.psf_stars_footprints.asAstropy()
797 # Run astrometry using PSF candidate stars
798 astrometry_matches, astrometry_meta = self._fit_astrometry(
799 result.exposure, result.psf_stars_footprints
800 )
801 self.metadata["astrometry_matches_count"] = len(astrometry_matches)
802 if "astrometry_matches" in self.config.optional_outputs:
803 result.astrometry_matches = lsst.meas.astrom.denormalizeMatches(astrometry_matches,
804 astrometry_meta)
805 result.psf_stars = result.psf_stars_footprints.asAstropy()
806
807 # Run the stars_detection subtask for the photometric calibration.
808 result.stars_footprints = self._find_stars(
809 result.exposure,
810 result.background,
811 id_generator,
812 background_to_photometric_ratio=result.background_to_photometric_ratio,
813 )
814 self._match_psf_stars(result.psf_stars_footprints, result.stars_footprints)
815
816 # Update the source cooordinates with the current wcs.
817 afwTable.updateSourceCoords(result.exposure.wcs, sourceList=result.stars_footprints)
818
819 summary_stat_catalog = result.stars_footprints
820 result.stars = result.stars_footprints.asAstropy()
821 self.metadata["star_count"] = np.sum(~result.stars["sky_source"])
822
823 # Validate the astrometric fit. Send in the stars_footprints
824 # catalog so that its coords get set to NaN if the fit is deemed
825 # a failure.
826 self.astrometry.check(result.exposure, result.stars_footprints, len(astrometry_matches))
827 result.stars = result.stars_footprints.asAstropy()
828 have_fit_astrometry = True
829
830 result.stars_footprints, photometry_matches, \
831 photometry_meta, photo_calib = self._fit_photometry(result.exposure, result.stars_footprints)
832 have_fit_photometry = True
833 self.metadata["photometry_matches_count"] = len(photometry_matches)
834 # fit_photometry returns a new catalog, so we need a new astropy table view.
835 result.stars = result.stars_footprints.asAstropy()
836 # summary stats don't make use of the calibrated fluxes, but we
837 # might as well use the best catalog we've got in case that
838 # changes, and help the old one get garbage-collected.
839 summary_stat_catalog = result.stars_footprints
840 if "photometry_matches" in self.config.optional_outputs:
841 result.photometry_matches = lsst.meas.astrom.denormalizeMatches(photometry_matches,
842 photometry_meta)
843 except pipeBase.AlgorithmError:
844 if not have_fit_psf:
845 result.exposure.setPsf(None)
846 if not have_fit_astrometry:
847 result.exposure.setWcs(None)
848 if not have_fit_photometry:
849 result.exposure.setPhotoCalib(None)
850 # Summary stat calculations can handle missing components gracefully,
851 # but we want to run them as late as possible (but still before we
852 # calibrate pixels, if we do that at all).
853 # So we run them after we succeed or if we get an AlgorithmError. We
854 # intentionally don't use 'finally' here because we don't want to run
855 # them if we get some other kind of error.
856 self._summarize(result.exposure, summary_stat_catalog, result.background)
857 raise
858 else:
859 self._summarize(result.exposure, summary_stat_catalog, result.background)
860
861 if self.config.do_calibrate_pixels:
863 result.exposure,
864 result.background,
865 background_to_photometric_ratio=result.background_to_photometric_ratio,
866 )
867 result.applied_photo_calib = photo_calib
868 else:
869 result.applied_photo_calib = None
870 return result
871
872 def _apply_illumination_correction(self, exposure, background_flat, illumination_correction):
873 """Apply the illumination correction to a background-flattened image.
874
875 Parameters
876 ----------
877 exposure : `lsst.afw.image.Exposure`
878 Exposure to convert to a photometric-flattened image.
879 background_flat : `lsst.afw.image.Exposure`
880 Flat image that had previously been applied to exposure.
881 illumination_correction : `lsst.afw.image.Exposure`
882 Illumination correction image to convert to photometric-flattened image.
883
884 Returns
885 -------
886 background_to_photometric_ratio : `lsst.afw.image.Image`
887 Ratio image to convert a photometric-flattened image to/from
888 a background-flattened image. Will be None if task not
889 configured to use the illumination correction.
890 """
891 if not self.config.do_illumination_correction:
892 return None
893
894 # From a raw image to a background-flattened image, we have:
895 # bfi = image / background_flat
896 # From a raw image to a photometric-flattened image, we have:
897 # pfi = image / reference_flux_flat
898 # pfi = image / (dome_flat * illumination_correction),
899 # where the illumination correction contains the jacobian
900 # of the wcs, converting to fluence units.
901 # Currently background_flat == dome_flat, so we have for the
902 # "background_to_photometric_ratio", the ratio of the background-
903 # flattened image to the photometric-flattened image:
904 # bfi / pfi = illumination_correction.
905
906 background_to_photometric_ratio = illumination_correction.image.clone()
907
908 # Dividing the ratio will convert a background-flattened image to
909 # a photometric-flattened image.
910 exposure.maskedImage /= background_to_photometric_ratio
911
912 exposure.metadata["LSST CALIB ILLUMCORR APPLIED"] = True
913
914 return background_to_photometric_ratio
915
916 def _compute_psf(self, exposure, id_generator, background_to_photometric_ratio=None):
917 """Find bright sources detected on an exposure and fit a PSF model to
918 them, repairing likely cosmic rays before detection.
919
920 Repair, detect, measure, and compute PSF twice, to ensure the PSF
921 model does not include contributions from cosmic rays.
922
923 Parameters
924 ----------
925 exposure : `lsst.afw.image.Exposure`
926 Exposure to detect and measure bright stars on.
927 id_generator : `lsst.meas.base.IdGenerator`
928 Object that generates source IDs and provides random seeds.
929 background_to_photometric_ratio : `lsst.afw.image.Image`, optional
930 Image to convert photometric-flattened image to
931 background-flattened image.
932
933 Returns
934 -------
935 sources : `lsst.afw.table.SourceCatalog`
936 Catalog of detected bright sources.
937 background : `lsst.afw.math.BackgroundList`
938 Background that was fit to the exposure during detection.
939 cell_set : `lsst.afw.math.SpatialCellSet`
940 PSF candidates returned by the psf determiner.
941 """
942 def log_psf(msg, addToMetadata=False):
943 """Log the parameters of the psf and background, with a prepended
944 message. There is also the option to add the PSF sigma to the task
945 metadata.
946
947 Parameters
948 ----------
949 msg : `str`
950 Message to prepend the log info with.
951 addToMetadata : `bool`, optional
952 Whether to add the final psf sigma value to the task metadata
953 (the default is False).
954 """
955 position = exposure.psf.getAveragePosition()
956 sigma = exposure.psf.computeShape(position).getDeterminantRadius()
957 dimensions = exposure.psf.computeImage(position).getDimensions()
958 median_background = np.median(background.getImage().array)
959 self.log.info("%s sigma=%0.4f, dimensions=%s; median background=%0.2f",
960 msg, sigma, dimensions, median_background)
961 if addToMetadata:
962 self.metadata["final_psf_sigma"] = sigma
963
964 self.log.info("First pass detection with Guassian PSF FWHM=%s pixels",
965 self.config.install_simple_psf.fwhm)
966 self.install_simple_psf.run(exposure=exposure)
967
968 background = self.psf_subtract_background.run(
969 exposure=exposure,
970 backgroundToPhotometricRatio=background_to_photometric_ratio,
971 ).background
972 log_psf("Initial PSF:")
973 self.psf_repair.run(exposure=exposure, keepCRs=True)
974
975 table = afwTable.SourceTable.make(self.psf_schema, id_generator.make_table_id_factory())
976 # Re-estimate the background during this detection step, so that
977 # measurement uses the most accurate background-subtraction.
978 detections = self.psf_detection.run(
979 table=table,
980 exposure=exposure,
981 background=background,
982 backgroundToPhotometricRatio=background_to_photometric_ratio,
983 )
984 self.metadata["initial_psf_positive_footprint_count"] = detections.numPos
985 self.metadata["initial_psf_negative_footprint_count"] = detections.numNeg
986 self.metadata["initial_psf_positive_peak_count"] = detections.numPosPeaks
987 self.metadata["initial_psf_negative_peak_count"] = detections.numNegPeaks
988 self.psf_source_measurement.run(detections.sources, exposure)
989 psf_result = self.psf_measure_psf.run(exposure=exposure, sources=detections.sources)
990 # Replace the initial PSF with something simpler for the second
991 # repair/detect/measure/measure_psf step: this can help it converge.
992 self.install_simple_psf.run(exposure=exposure)
993
994 log_psf("Rerunning with simple PSF:")
995 # TODO investigation: Should we only re-run repair here, to use the
996 # new PSF? Maybe we *do* need to re-run measurement with PsfFlux, to
997 # use the fitted PSF?
998 # TODO investigation: do we need a separate measurement task here
999 # for the post-psf_measure_psf step, since we only want to do PsfFlux
1000 # and GaussianFlux *after* we have a PSF? Maybe that's not relevant
1001 # once DM-39203 is merged?
1002 self.psf_repair.run(exposure=exposure, keepCRs=True)
1003 # Re-estimate the background during this detection step, so that
1004 # measurement uses the most accurate background-subtraction.
1005 detections = self.psf_detection.run(
1006 table=table,
1007 exposure=exposure,
1008 background=background,
1009 backgroundToPhotometricRatio=background_to_photometric_ratio,
1010 )
1011 self.metadata["simple_psf_positive_footprint_count"] = detections.numPos
1012 self.metadata["simple_psf_negative_footprint_count"] = detections.numNeg
1013 self.metadata["simple_psf_positive_peak_count"] = detections.numPosPeaks
1014 self.metadata["simple_psf_negative_peak_count"] = detections.numNegPeaks
1015 self.psf_source_measurement.run(detections.sources, exposure)
1016 psf_result = self.psf_measure_psf.run(exposure=exposure, sources=detections.sources)
1017
1018 log_psf("Final PSF:", addToMetadata=True)
1019
1020 # Final repair with final PSF, removing cosmic rays this time.
1021 self.psf_repair.run(exposure=exposure)
1022 # Final measurement with the CRs removed.
1023 self.psf_source_measurement.run(detections.sources, exposure)
1024
1025 if detections.sources["slot_Centroid_flag"].all():
1026 psf_shape = exposure.psf.computeShape(exposure.psf.getAveragePosition())
1027 raise AllCentroidsFlaggedError(n_sources=len(detections.sources),
1028 psf_shape_ixx=psf_shape.getIxx(),
1029 psf_shape_iyy=psf_shape.getIyy(),
1030 psf_shape_ixy=psf_shape.getIxy(),
1031 psf_size=psf_shape.getDeterminantRadius()
1032 )
1033
1034 # PSF is set on exposure; candidates are returned to use for
1035 # calibration flux normalization and aperture corrections.
1036 return detections.sources, background, psf_result.cellSet
1037
1038 def _measure_aperture_correction(self, exposure, bright_sources):
1039 """Measure and set the ApCorrMap on the Exposure, using
1040 previously-measured bright sources.
1041
1042 This function first normalizes the calibration flux and then
1043 the full set of aperture corrections are measured relative
1044 to this normalized calibration flux.
1045
1046 Parameters
1047 ----------
1048 exposure : `lsst.afw.image.Exposure`
1049 Exposure to set the ApCorrMap on.
1050 bright_sources : `lsst.afw.table.SourceCatalog`
1051 Catalog of detected bright sources; modified to include columns
1052 necessary for point source determination for the aperture correction
1053 calculation.
1054 """
1055 norm_ap_corr_map = self.psf_normalized_calibration_flux.run(
1056 exposure=exposure,
1057 catalog=bright_sources,
1058 ).ap_corr_map
1059
1060 ap_corr_map = self.measure_aperture_correction.run(exposure, bright_sources).apCorrMap
1061
1062 # Need to merge the aperture correction map from the normalization.
1063 for key in norm_ap_corr_map:
1064 ap_corr_map[key] = norm_ap_corr_map[key]
1065
1066 exposure.info.setApCorrMap(ap_corr_map)
1067
1068 def _find_stars(self, exposure, background, id_generator, background_to_photometric_ratio=None):
1069 """Detect stars on an exposure that has a PSF model, and measure their
1070 PSF, circular aperture, compensated gaussian fluxes.
1071
1072 Parameters
1073 ----------
1074 exposure : `lsst.afw.image.Exposure`
1075 Exposure to detect and measure stars on.
1076 background : `lsst.afw.math.BackgroundList`
1077 Background that was fit to the exposure during detection;
1078 modified in-place during subsequent detection.
1079 id_generator : `lsst.meas.base.IdGenerator`
1080 Object that generates source IDs and provides random seeds.
1081 background_to_photometric_ratio : `lsst.afw.image.Image`, optional
1082 Image to convert photometric-flattened image to
1083 background-flattened image.
1084
1085 Returns
1086 -------
1087 stars : `SourceCatalog`
1088 Sources that are very likely to be stars, with a limited set of
1089 measurements performed on them.
1090 """
1091 table = afwTable.SourceTable.make(self.initial_stars_schema.schema,
1092 id_generator.make_table_id_factory())
1093 # Re-estimate the background during this detection step, so that
1094 # measurement uses the most accurate background-subtraction.
1095 detections = self.star_detection.run(
1096 table=table,
1097 exposure=exposure,
1098 background=background,
1099 backgroundToPhotometricRatio=background_to_photometric_ratio,
1100 )
1101 sources = detections.sources
1102 self.star_sky_sources.run(exposure.mask, id_generator.catalog_id, sources)
1103
1104 # TODO investigation: Could this deblender throw away blends of non-PSF sources?
1105 self.star_deblend.run(exposure=exposure, sources=sources)
1106 # The deblender may not produce a contiguous catalog; ensure
1107 # contiguity for subsequent tasks.
1108 if not sources.isContiguous():
1109 sources = sources.copy(deep=True)
1110
1111 # Measure everything, and use those results to select only stars.
1112 self.star_measurement.run(sources, exposure)
1113 self.metadata["post_deblend_source_count"] = np.sum(~sources["sky_source"])
1114 self.metadata["saturated_source_count"] = np.sum(sources["base_PixelFlags_flag_saturated"])
1115 self.metadata["bad_source_count"] = np.sum(sources["base_PixelFlags_flag_bad"])
1116
1117 # Run the normalization calibration flux task to apply the
1118 # normalization correction to create normalized
1119 # calibration fluxes.
1120 self.star_normalized_calibration_flux.run(exposure=exposure, catalog=sources)
1121 self.star_apply_aperture_correction.run(sources, exposure.apCorrMap)
1122 self.star_catalog_calculation.run(sources)
1123 self.star_set_primary_flags.run(sources)
1124
1125 result = self.star_selector.run(sources)
1126 # The star selector may not produce a contiguous catalog.
1127 if not result.sourceCat.isContiguous():
1128 return result.sourceCat.copy(deep=True)
1129 else:
1130 return result.sourceCat
1131
1132 def _match_psf_stars(self, psf_stars, stars):
1133 """Match calibration stars to psf stars, to identify which were psf
1134 candidates, and which were used or reserved during psf measurement
1135 and the astrometric fit.
1136
1137 Parameters
1138 ----------
1139 psf_stars : `lsst.afw.table.SourceCatalog`
1140 PSF candidate stars that were sent to the psf determiner and
1141 used in the astrometric fit. Used to populate psf and astrometry
1142 related flag fields.
1143 stars : `lsst.afw.table.SourceCatalog`
1144 Stars that will be used for calibration; psf-related fields will
1145 be updated in-place.
1146
1147 Notes
1148 -----
1149 This code was adapted from CalibrateTask.copyIcSourceFields().
1150 """
1151 control = afwTable.MatchControl()
1152 # Return all matched objects, to separate blends.
1153 control.findOnlyClosest = False
1154 matches = afwTable.matchXy(psf_stars, stars, 3.0, control)
1155 deblend_key = stars.schema["deblend_nChild"].asKey()
1156 matches = [m for m in matches if m[1].get(deblend_key) == 0]
1157
1158 # Because we had to allow multiple matches to handle parents, we now
1159 # need to prune to the best (closest) matches.
1160 # Closest matches is a dict of psf_stars source ID to Match record
1161 # (psf_stars source, sourceCat source, distance in pixels).
1162 best = {}
1163 for match_psf, match_stars, d in matches:
1164 match = best.get(match_psf.getId())
1165 if match is None or d <= match[2]:
1166 best[match_psf.getId()] = (match_psf, match_stars, d)
1167 matches = list(best.values())
1168 # We'll use this to construct index arrays into each catalog.
1169 ids = np.array([(match_psf.getId(), match_stars.getId()) for match_psf, match_stars, d in matches]).T
1170
1171 if (n_matches := len(matches)) == 0:
1172 raise NoPsfStarsToStarsMatchError(n_psf_stars=len(psf_stars), n_stars=len(stars))
1173
1174 self.log.info("%d psf/astrometry stars out of %d matched %d calib stars",
1175 n_matches, len(psf_stars), len(stars))
1176 self.metadata["matched_psf_star_count"] = n_matches
1177
1178 # Check that no stars sources are listed twice; we already know
1179 # that each match has a unique psf_stars id, due to using as the key
1180 # in best above.
1181 n_unique = len(set(m[1].getId() for m in matches))
1182 if n_unique != n_matches:
1183 self.log.warning("%d psf_stars matched only %d stars", n_matches, n_unique)
1184
1185 # The indices of the IDs, so we can update the flag fields as arrays.
1186 idx_psf_stars = np.searchsorted(psf_stars["id"], ids[0])
1187 idx_stars = np.searchsorted(stars["id"], ids[1])
1188 for field in self.psf_fields:
1189 result = np.zeros(len(stars), dtype=bool)
1190 result[idx_stars] = psf_stars[field][idx_psf_stars]
1191 stars[field] = result
1192 stars['psf_id'][idx_stars] = psf_stars['id'][idx_psf_stars]
1193
1194 def _fit_astrometry(self, exposure, stars):
1195 """Fit an astrometric model to the data and return the reference
1196 matches used in the fit, and the fitted WCS.
1197
1198 Parameters
1199 ----------
1200 exposure : `lsst.afw.image.Exposure`
1201 Exposure that is being fit, to get PSF and other metadata from.
1202 Modified to add the fitted skyWcs.
1203 stars : `SourceCatalog`
1204 Good stars selected for use in calibration, with RA/Dec coordinates
1205 computed from the pixel positions and fitted WCS.
1206
1207 Returns
1208 -------
1209 matches : `list` [`lsst.afw.table.ReferenceMatch`]
1210 Reference/stars matches used in the fit.
1211 """
1212 result = self.astrometry.run(stars, exposure)
1213 return result.matches, result.matchMeta
1214
1215 def _fit_photometry(self, exposure, stars):
1216 """Fit a photometric model to the data and return the reference
1217 matches used in the fit, and the fitted PhotoCalib.
1218
1219 Parameters
1220 ----------
1221 exposure : `lsst.afw.image.Exposure`
1222 Exposure that is being fit, to get PSF and other metadata from.
1223 Has the fit `lsst.afw.image.PhotoCalib` attached, with pixel values
1224 unchanged.
1225 stars : `lsst.afw.table.SourceCatalog`
1226 Good stars selected for use in calibration.
1227 background : `lsst.afw.math.BackgroundList`
1228 Background that was fit to the exposure during detection of the
1229 above stars.
1230
1231 Returns
1232 -------
1233 calibrated_stars : `lsst.afw.table.SourceCatalog`
1234 Star catalog with flux/magnitude columns computed from the fitted
1235 photoCalib (instFlux columns are retained as well).
1236 matches : `list` [`lsst.afw.table.ReferenceMatch`]
1237 Reference/stars matches used in the fit.
1238 matchMeta : `lsst.daf.base.PropertyList`
1239 Metadata needed to unpersist matches, as returned by the matcher.
1240 photo_calib : `lsst.afw.image.PhotoCalib`
1241 Photometric calibration that was fit to the star catalog.
1242 """
1243 result = self.photometry.run(exposure, stars)
1244 calibrated_stars = result.photoCalib.calibrateCatalog(stars)
1245 exposure.setPhotoCalib(result.photoCalib)
1246 return calibrated_stars, result.matches, result.matchMeta, result.photoCalib
1247
1248 def _apply_photometry(self, exposure, background, background_to_photometric_ratio=None):
1249 """Apply the photometric model attached to the exposure to the
1250 exposure's pixels and an associated background model.
1251
1252 Parameters
1253 ----------
1254 exposure : `lsst.afw.image.Exposure`
1255 Exposure with the target `lsst.afw.image.PhotoCalib` attached.
1256 On return, pixel values will be calibrated and an identity
1257 photometric transform will be attached.
1258 background : `lsst.afw.math.BackgroundList`
1259 Background model to convert to nanojansky units in place.
1260 background_to_photometric_ratio : `lsst.afw.image.Image`, optional
1261 Image to convert photometric-flattened image to
1262 background-flattened image.
1263 """
1264 photo_calib = exposure.getPhotoCalib()
1265 exposure.maskedImage = photo_calib.calibrateImage(exposure.maskedImage)
1266 identity = afwImage.PhotoCalib(1.0,
1267 photo_calib.getCalibrationErr(),
1268 bbox=exposure.getBBox())
1269 exposure.setPhotoCalib(identity)
1270 exposure.metadata["BUNIT"] = "nJy"
1271
1272 assert photo_calib._isConstant, \
1273 "Background calibration assumes a constant PhotoCalib; PhotoCalTask should always return that."
1274
1275 for bg in background:
1276 # The statsImage is a view, but we can't assign to a function call in python.
1277 binned_image = bg[0].getStatsImage()
1278 binned_image *= photo_calib.getCalibrationMean()
1279
1280 def _summarize(self, exposure, stars, background):
1281 """Compute summary statistics on the exposure and update in-place the
1282 calibrations attached to it.
1283
1284 Parameters
1285 ----------
1286 exposure : `lsst.afw.image.Exposure`
1287 Exposure that was calibrated, to get PSF and other metadata from.
1288 Should be in instrumental units with the photometric calibration
1289 attached.
1290 Modified to contain the computed summary statistics.
1291 stars : `SourceCatalog`
1292 Good stars selected used in calibration.
1293 background : `lsst.afw.math.BackgroundList`
1294 Background that was fit to the exposure during detection of the
1295 above stars. Should be in instrumental units.
1296 """
1297 summary = self.compute_summary_stats.run(exposure, stars, background)
1298 exposure.info.setSummaryStats(summary)
1299
1300 def _recordMaskedPixelFractions(self, exposure):
1301 """Record the fraction of all the pixels in an exposure
1302 that are masked with a given flag. Each fraction is
1303 recorded in the task metadata. One record per flag type.
1304
1305 Parameters
1306 ----------
1307 exposure : `lsst.afw.image.ExposureF`
1308 The target exposure to calculate masked pixel fractions for.
1309 """
1310
1311 mask = exposure.mask
1312 maskPlanes = list(mask.getMaskPlaneDict().keys())
1313 for maskPlane in maskPlanes:
1314 self.metadata[f"{maskPlane.lower()}_mask_fraction"] = (
1315 evaluateMaskFraction(mask, maskPlane)
1316 )
The photometric calibration of an exposure.
Definition PhotoCalib.h:114
Pass parameters to algorithms that match list of sources.
Definition Match.h:45
An integer coordinate rectangle.
Definition Box.h:55
__init__(self, n_sources, psf_shape_ixx, psf_shape_iyy, psf_shape_ixy, psf_size)
runQuantum(self, butlerQC, inputRefs, outputRefs)
_find_stars(self, exposure, background, id_generator, background_to_photometric_ratio=None)
_apply_illumination_correction(self, exposure, background_flat, illumination_correction)
_measure_aperture_correction(self, exposure, bright_sources)
_apply_photometry(self, exposure, background, background_to_photometric_ratio=None)
__init__(self, initial_stars_schema=None, **kwargs)
_compute_psf(self, exposure, id_generator, background_to_photometric_ratio=None)
_summarize(self, exposure, stars, background)
run(self, *, exposures, id_generator=None, result=None, background_flat=None, illumination_correction=None)
SourceMatchVector matchXy(SourceCatalog const &cat1, SourceCatalog const &cat2, double radius, MatchControl const &mc=MatchControl())
Compute all tuples (s1,s2,d) where s1 belings to cat1, s2 belongs to cat2 and d, the distance between...
Definition Match.cc:305
void updateSourceCoords(geom::SkyWcs const &wcs, SourceCollection &sourceList, bool include_covariance=true)
Update sky coordinates in a collection of source objects.
Definition wcsUtils.cc:125