LSST Applications 26.0.0,g0265f82a02+6660c170cc,g07994bdeae+30b05a742e,g0a0026dc87+17526d298f,g0a60f58ba1+17526d298f,g0e4bf8285c+96dd2c2ea9,g0ecae5effc+c266a536c8,g1e7d6db67d+6f7cb1f4bb,g26482f50c6+6346c0633c,g2bbee38e9b+6660c170cc,g2cc88a2952+0a4e78cd49,g3273194fdb+f6908454ef,g337abbeb29+6660c170cc,g337c41fc51+9a8f8f0815,g37c6e7c3d5+7bbafe9d37,g44018dc512+6660c170cc,g4a941329ef+4f7594a38e,g4c90b7bd52+5145c320d2,g58be5f913a+bea990ba40,g635b316a6c+8d6b3a3e56,g67924a670a+bfead8c487,g6ae5381d9b+81bc2a20b4,g93c4d6e787+26b17396bd,g98cecbdb62+ed2cb6d659,g98ffbb4407+81bc2a20b4,g9ddcbc5298+7f7571301f,ga1e77700b3+99e9273977,gae46bcf261+6660c170cc,gb2715bf1a1+17526d298f,gc86a011abf+17526d298f,gcf0d15dbbd+96dd2c2ea9,gdaeeff99f8+0d8dbea60f,gdb4ec4c597+6660c170cc,ge23793e450+96dd2c2ea9,gf041782ebf+171108ac67
LSST Data Management Base Package
Loading...
Searching...
No Matches
forcedPhotCcd.py
Go to the documentation of this file.
1# This file is part of meas_base.
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
22import warnings
23
24import pandas as pd
25import numpy as np
26
27import lsst.pex.config
29import lsst.pipe.base
30import lsst.geom
32import lsst.afw.geom
33import lsst.afw.image
34import lsst.afw.table
35import lsst.sphgeom
36
37from lsst.utils.introspection import find_outside_stacklevel
38from lsst.pipe.base import PipelineTaskConnections
39import lsst.pipe.base.connectionTypes as cT
40
41import lsst.pipe.base as pipeBase
42from lsst.skymap import BaseSkyMap
43
44from .forcedMeasurement import ForcedMeasurementTask
45from .applyApCorr import ApplyApCorrTask
46from .catalogCalculation import CatalogCalculationTask
47from ._id_generator import DetectorVisitIdGeneratorConfig
48
49__all__ = ("ForcedPhotCcdConfig", "ForcedPhotCcdTask",
50 "ForcedPhotCcdFromDataFrameTask", "ForcedPhotCcdFromDataFrameConfig")
51
52
53class ForcedPhotCcdConnections(PipelineTaskConnections,
54 dimensions=("instrument", "visit", "detector", "skymap", "tract"),
55 defaultTemplates={"inputCoaddName": "deep",
56 "inputName": "calexp",
57 "skyWcsName": "gbdesAstrometricFit",
58 "photoCalibName": "fgcm"},
59 # TODO: remove on DM-39854
60 deprecatedTemplates={"skyWcsName": "Deprecated; will be removed after v26.",
61 "photoCalibName": "Deprecated; will be removed after v26."
62 }):
63 inputSchema = cT.InitInput(
64 doc="Schema for the input measurement catalogs.",
65 name="{inputCoaddName}Coadd_ref_schema",
66 storageClass="SourceCatalog",
67 )
68 outputSchema = cT.InitOutput(
69 doc="Schema for the output forced measurement catalogs.",
70 name="forced_src_schema",
71 storageClass="SourceCatalog",
72 )
73 exposure = cT.Input(
74 doc="Input exposure to perform photometry on.",
75 name="{inputName}",
76 storageClass="ExposureF",
77 dimensions=["instrument", "visit", "detector"],
78 )
79 refCat = cT.Input(
80 doc="Catalog of shapes and positions at which to force photometry.",
81 name="{inputCoaddName}Coadd_ref",
82 storageClass="SourceCatalog",
83 dimensions=["skymap", "tract", "patch"],
84 multiple=True,
85 deferLoad=True,
86 )
87 skyMap = cT.Input(
88 doc="SkyMap dataset that defines the coordinate system of the reference catalog.",
89 name=BaseSkyMap.SKYMAP_DATASET_TYPE_NAME,
90 storageClass="SkyMap",
91 dimensions=["skymap"],
92 )
93 skyCorr = cT.Input(
94 doc="Input Sky Correction to be subtracted from the calexp if doApplySkyCorr=True",
95 name="skyCorr",
96 storageClass="Background",
97 dimensions=("instrument", "visit", "detector"),
98 )
99 visitSummary = cT.Input(
100 doc="Input visit-summary catalog with updated calibration objects.",
101 name="finalVisitSummary",
102 storageClass="ExposureCatalog",
103 dimensions=("instrument", "visit"),
104 )
105 externalSkyWcsTractCatalog = cT.Input(
106 doc=("Per-tract, per-visit wcs calibrations. These catalogs use the detector "
107 "id for the catalog id, sorted on id for fast lookup."),
108 name="{skyWcsName}SkyWcsCatalog",
109 storageClass="ExposureCatalog",
110 dimensions=["instrument", "visit", "tract"],
111 # TODO: remove on DM-39854
112 deprecated="Deprecated in favor of 'visitSummary'; will be removed after v26."
113 )
114 externalSkyWcsGlobalCatalog = cT.Input(
115 doc=("Per-visit wcs calibrations computed globally (with no tract information). "
116 "These catalogs use the detector id for the catalog id, sorted on id for "
117 "fast lookup."),
118 name="finalVisitSummary",
119 storageClass="ExposureCatalog",
120 dimensions=["instrument", "visit"],
121 # TODO: remove on DM-39854
122 deprecated="Deprecated in favor of 'visitSummary'; will be removed after v26."
123 )
124 externalPhotoCalibTractCatalog = cT.Input(
125 doc=("Per-tract, per-visit photometric calibrations. These catalogs use the "
126 "detector id for the catalog id, sorted on id for fast lookup."),
127 name="{photoCalibName}PhotoCalibCatalog",
128 storageClass="ExposureCatalog",
129 dimensions=["instrument", "visit", "tract"],
130 # TODO: remove on DM-39854
131 deprecated="Deprecated in favor of 'visitSummary'; will be removed after v26."
132 )
133 externalPhotoCalibGlobalCatalog = cT.Input(
134 doc=("Per-visit photometric calibrations computed globally (with no tract "
135 "information). These catalogs use the detector id for the catalog id, "
136 "sorted on id for fast lookup."),
137 name="finalVisitSummary",
138 storageClass="ExposureCatalog",
139 dimensions=["instrument", "visit"],
140 # TODO: remove on DM-39854
141 deprecated="Deprecated in favor of 'visitSummary'; will be removed after v26."
142 )
143 finalizedPsfApCorrCatalog = cT.Input(
144 doc=("Per-visit finalized psf models and aperture correction maps. "
145 "These catalogs use the detector id for the catalog id, "
146 "sorted on id for fast lookup."),
147 name="finalized_psf_ap_corr_catalog",
148 storageClass="ExposureCatalog",
149 dimensions=["instrument", "visit"],
150 # TODO: remove on DM-39854
151 deprecated="Deprecated in favor of 'visitSummary'; will be removed after v26."
152 )
153 measCat = cT.Output(
154 doc="Output forced photometry catalog.",
155 name="forced_src",
156 storageClass="SourceCatalog",
157 dimensions=["instrument", "visit", "detector", "skymap", "tract"],
158 )
159
160 def __init__(self, *, config=None):
161 super().__init__(config=config)
162 if not config.doApplySkyCorr:
163 self.inputs.remove("skyCorr")
164 if config.doApplyExternalSkyWcs:
165 if config.useGlobalExternalSkyWcs:
166 self.inputs.remove("externalSkyWcsTractCatalog")
167 else:
168 self.inputs.remove("externalSkyWcsGlobalCatalog")
169 else:
170 self.inputs.remove("externalSkyWcsTractCatalog")
171 self.inputs.remove("externalSkyWcsGlobalCatalog")
172 if config.doApplyExternalPhotoCalib:
173 if config.useGlobalExternalPhotoCalib:
174 self.inputs.remove("externalPhotoCalibTractCatalog")
175 else:
176 self.inputs.remove("externalPhotoCalibGlobalCatalog")
177 else:
178 self.inputs.remove("externalPhotoCalibTractCatalog")
179 self.inputs.remove("externalPhotoCalibGlobalCatalog")
180 if not config.doApplyFinalizedPsf:
181 self.inputs.remove("finalizedPsfApCorrCatalog")
182
183
184class ForcedPhotCcdConfig(pipeBase.PipelineTaskConfig,
185 pipelineConnections=ForcedPhotCcdConnections):
186 """Config class for forced measurement driver task."""
188 target=ForcedMeasurementTask,
189 doc="subtask to do forced measurement"
190 )
191 coaddName = lsst.pex.config.Field(
192 doc="coadd name: typically one of deep or goodSeeing",
193 dtype=str,
194 default="deep",
195 )
196 doApCorr = lsst.pex.config.Field(
197 dtype=bool,
198 default=True,
199 doc="Run subtask to apply aperture corrections"
200 )
202 target=ApplyApCorrTask,
203 doc="Subtask to apply aperture corrections"
204 )
205 catalogCalculation = lsst.pex.config.ConfigurableField(
206 target=CatalogCalculationTask,
207 doc="Subtask to run catalogCalculation plugins on catalog"
208 )
209 doApplyExternalPhotoCalib = lsst.pex.config.Field(
210 dtype=bool,
211 default=False,
212 doc=("Whether to apply external photometric calibration via an "
213 "`lsst.afw.image.PhotoCalib` object."),
214 # TODO: remove on DM-39854
215 deprecated="Removed in favor of the 'visitSummary' connection; will be removed after v26.",
216 )
217 useGlobalExternalPhotoCalib = lsst.pex.config.Field(
218 dtype=bool,
219 default=False,
220 doc=("When using doApplyExternalPhotoCalib, use 'global' calibrations "
221 "that are not run per-tract. When False, use per-tract photometric "
222 "calibration files."),
223 # TODO: remove on DM-39854
224 deprecated="Removed in favor of the 'visitSummary' connection; will be removed after v26.",
225 )
226 doApplyExternalSkyWcs = lsst.pex.config.Field(
227 dtype=bool,
228 default=False,
229 doc=("Whether to apply external astrometric calibration via an "
230 "`lsst.afw.geom.SkyWcs` object."),
231 # TODO: remove on DM-39854
232 deprecated="Removed in favor of the 'visitSummary' connection; will be removed after v26.",
233 )
234 useGlobalExternalSkyWcs = lsst.pex.config.Field(
235 dtype=bool,
236 default=False,
237 doc=("When using doApplyExternalSkyWcs, use 'global' calibrations "
238 "that are not run per-tract. When False, use per-tract wcs "
239 "files."),
240 # TODO: remove on DM-39854
241 deprecated="Removed in favor of the 'visitSummary' connection; will be removed after v26.",
242 )
243 doApplyFinalizedPsf = lsst.pex.config.Field(
244 dtype=bool,
245 default=False,
246 doc="Whether to apply finalized psf models and aperture correction map.",
247 # TODO: remove on DM-39854
248 deprecated="Removed in favor of the 'visitSummary' connection; will be removed after v26.",
249 )
250 doApplySkyCorr = lsst.pex.config.Field(
251 dtype=bool,
252 default=False,
253 doc="Apply sky correction?",
254 )
255 includePhotoCalibVar = lsst.pex.config.Field(
256 dtype=bool,
257 default=False,
258 doc="Add photometric calibration variance to warp variance plane?",
259 )
260 footprintSource = lsst.pex.config.ChoiceField(
261 dtype=str,
262 doc="Where to obtain footprints to install in the measurement catalog, prior to measurement.",
263 allowed={
264 "transformed": "Transform footprints from the reference catalog (downgrades HeavyFootprints).",
265 "psf": ("Use the scaled shape of the PSF at the position of each source (does not generate "
266 "HeavyFootprints)."),
267 },
268 optional=True,
269 default="transformed",
270 )
271 psfFootprintScaling = lsst.pex.config.Field(
272 dtype=float,
273 doc="Scaling factor to apply to the PSF shape when footprintSource='psf' (ignored otherwise).",
274 default=3.0,
275 )
276 idGenerator = DetectorVisitIdGeneratorConfig.make_field()
277
278 def setDefaults(self):
279 # Docstring inherited.
280 super().setDefaults()
281 # Footprints here will not be entirely correct, so don't try to make
282 # a biased correction for blended neighbors.
283 self.measurement.doReplaceWithNoise = False
284 # Only run a minimal set of plugins, as these measurements are only
285 # needed for PSF-like sources.
286 self.measurement.plugins.names = ["base_PixelFlags",
287 "base_TransformedCentroid",
288 "base_PsfFlux",
289 "base_LocalBackground",
290 "base_LocalPhotoCalib",
291 "base_LocalWcs",
292 ]
293 self.measurement.slots.shape = None
294 # Make catalogCalculation a no-op by default as no modelFlux is setup
295 # by default in ForcedMeasurementTask.
296 self.catalogCalculation.plugins.names = []
297
298
299class ForcedPhotCcdTask(pipeBase.PipelineTask):
300 """A pipeline task for performing forced measurement on CCD images.
301
302 Parameters
303 ----------
304 refSchema : `lsst.afw.table.Schema`, optional
305 The schema of the reference catalog, passed to the constructor of the
306 references subtask. Optional, but must be specified if ``initInputs``
307 is not; if both are specified, ``initInputs`` takes precedence.
308 initInputs : `dict`
309 Dictionary that can contain a key ``inputSchema`` containing the
310 schema. If present will override the value of ``refSchema``.
311 **kwargs
312 Keyword arguments are passed to the supertask constructor.
313 """
314
315 ConfigClass = ForcedPhotCcdConfig
316 _DefaultName = "forcedPhotCcd"
317 dataPrefix = ""
318
319 def __init__(self, refSchema=None, initInputs=None, **kwargs):
320 super().__init__(**kwargs)
321
322 if initInputs is not None:
323 refSchema = initInputs['inputSchema'].schema
324
325 if refSchema is None:
326 raise ValueError("No reference schema provided.")
327
328 self.makeSubtask("measurement", refSchema=refSchema)
329 # It is necessary to get the schema internal to the forced measurement
330 # task until such a time that the schema is not owned by the
331 # measurement task, but is passed in by an external caller.
332 if self.config.doApCorr:
333 self.makeSubtask("applyApCorr", schema=self.measurement.schema)
334 self.makeSubtask('catalogCalculation', schema=self.measurement.schema)
335 self.outputSchema = lsst.afw.table.SourceCatalog(self.measurement.schema)
336
337 def runQuantum(self, butlerQC, inputRefs, outputRefs):
338 inputs = butlerQC.get(inputRefs)
339
340 tract = butlerQC.quantum.dataId['tract']
341 skyMap = inputs.pop('skyMap')
342 inputs['refWcs'] = skyMap[tract].getWcs()
343
344 # Connections only exist if they are configured to be used.
345 skyCorr = inputs.pop('skyCorr', None)
346 if self.config.useGlobalExternalSkyWcs:
347 externalSkyWcsCatalog = inputs.pop('externalSkyWcsGlobalCatalog', None)
348 else:
349 externalSkyWcsCatalog = inputs.pop('externalSkyWcsTractCatalog', None)
350 if self.config.useGlobalExternalPhotoCalib:
351 externalPhotoCalibCatalog = inputs.pop('externalPhotoCalibGlobalCatalog', None)
352 else:
353 externalPhotoCalibCatalog = inputs.pop('externalPhotoCalibTractCatalog', None)
354 finalizedPsfApCorrCatalog = inputs.pop('finalizedPsfApCorrCatalog', None)
355
356 inputs['exposure'] = self.prepareCalibratedExposure(
357 inputs['exposure'],
358 skyCorr=skyCorr,
359 externalSkyWcsCatalog=externalSkyWcsCatalog,
360 externalPhotoCalibCatalog=externalPhotoCalibCatalog,
361 finalizedPsfApCorrCatalog=finalizedPsfApCorrCatalog,
362 visitSummary=inputs.pop("visitSummary"),
363 )
364
365 inputs['refCat'] = self.mergeAndFilterReferences(inputs['exposure'], inputs['refCat'],
366 inputs['refWcs'])
367
368 if inputs['refCat'] is None:
369 self.log.info("No WCS for exposure %s. No %s catalog will be written.",
370 butlerQC.quantum.dataId, outputRefs.measCat.datasetType.name)
371 else:
372 inputs['measCat'], inputs['exposureId'] = self.generateMeasCat(inputRefs.exposure.dataId,
373 inputs['exposure'],
374 inputs['refCat'], inputs['refWcs'])
375 self.attachFootprints(inputs['measCat'], inputs['refCat'], inputs['exposure'], inputs['refWcs'])
376 outputs = self.run(**inputs)
377 butlerQC.put(outputs, outputRefs)
378
379 def prepareCalibratedExposure(self, exposure, skyCorr=None, externalSkyWcsCatalog=None,
380 externalPhotoCalibCatalog=None, finalizedPsfApCorrCatalog=None,
381 visitSummary=None):
382 """Prepare a calibrated exposure and apply external calibrations
383 and sky corrections if so configured.
384
385 Parameters
386 ----------
387 exposure : `lsst.afw.image.exposure.Exposure`
388 Input exposure to adjust calibrations.
389 skyCorr : `lsst.afw.math.backgroundList`, optional
390 Sky correction frame to apply if doApplySkyCorr=True.
391 externalSkyWcsCatalog : `lsst.afw.table.ExposureCatalog`, optional
392 Exposure catalog with external skyWcs to be applied if
393 config.doApplyExternalSkyWcs=True. Catalog uses the detector id
394 for the catalog id, sorted on id for fast lookup. Deprecated in
395 favor of ``visitSummary`` and will be removed after v26.
396 externalPhotoCalibCatalog : `lsst.afw.table.ExposureCatalog`, optional
397 Exposure catalog with external photoCalib to be applied if
398 config.doApplyExternalPhotoCalib=True. Catalog uses the detector
399 id for the catalog id, sorted on id for fast lookup. Deprecated in
400 favor of ``visitSummary`` and will be removed after v26.
401 finalizedPsfApCorrCatalog : `lsst.afw.table.ExposureCatalog`, optional
402 Exposure catalog with finalized psf models and aperture correction
403 maps to be applied if config.doApplyFinalizedPsf=True. Catalog
404 uses the detector id for the catalog id, sorted on id for fast
405 lookup. Deprecated in favor of ``visitSummary`` and will be removed
406 after v26.
407 visitSummary : `lsst.afw.table.ExposureCatalog`, optional
408 Exposure catalog with update calibrations; any not-None calibration
409 objects attached will be used. These are applied first and may be
410 overridden by other arguments.
411
412 Returns
413 -------
414 exposure : `lsst.afw.image.exposure.Exposure`
415 Exposure with adjusted calibrations.
416 """
417 detectorId = exposure.getInfo().getDetector().getId()
418
419 if visitSummary is not None:
420 row = visitSummary.find(detectorId)
421 if row is None:
422 raise RuntimeError(f"Detector id {detectorId} not found in visitSummary.")
423 if (photoCalib := row.getPhotoCalib()) is not None:
424 exposure.setPhotoCalib(photoCalib)
425 if (skyWcs := row.getWcs()) is not None:
426 exposure.setWcs(skyWcs)
427 if (psf := row.getPsf()) is not None:
428 exposure.setPsf(psf)
429 if (apCorrMap := row.getApCorrMap()) is not None:
430 exposure.info.setApCorrMap(apCorrMap)
431
432 if externalPhotoCalibCatalog is not None:
433 # TODO: remove on DM-39854
434 warnings.warn(
435 "The 'externalPhotoCalibCatalog' argument is deprecated in favor of 'visitSummary' and will "
436 "be removed after v26.",
437 FutureWarning,
438 stacklevel=find_outside_stacklevel("lsst.meas.base"),
439 )
440 row = externalPhotoCalibCatalog.find(detectorId)
441 if row is None:
442 self.log.warning("Detector id %s not found in externalPhotoCalibCatalog; "
443 "Using original photoCalib.", detectorId)
444 else:
445 photoCalib = row.getPhotoCalib()
446 if photoCalib is None:
447 self.log.warning("Detector id %s has None for photoCalib in externalPhotoCalibCatalog; "
448 "Using original photoCalib.", detectorId)
449 else:
450 exposure.setPhotoCalib(photoCalib)
451
452 if externalSkyWcsCatalog is not None:
453 # TODO: remove on DM-39854
454 warnings.warn(
455 "The 'externalSkyWcsCatalog' argument is deprecated in favor of 'visitSummary' and will "
456 "be removed after v26.",
457 FutureWarning,
458 stacklevel=find_outside_stacklevel("lsst.meas.base"),
459 )
460 row = externalSkyWcsCatalog.find(detectorId)
461 if row is None:
462 self.log.warning("Detector id %s not found in externalSkyWcsCatalog; "
463 "Using original skyWcs.", detectorId)
464 else:
465 skyWcs = row.getWcs()
466 if skyWcs is None:
467 self.log.warning("Detector id %s has None for skyWcs in externalSkyWcsCatalog; "
468 "Using original skyWcs.", detectorId)
469 else:
470 exposure.setWcs(skyWcs)
471
472 if finalizedPsfApCorrCatalog is not None:
473 # TODO: remove on DM-39854
474 warnings.warn(
475 "The 'finalizedPsfApCorrCatalog' argument is deprecated in favor of 'visitSummary' and will "
476 "be removed after v26.",
477 FutureWarning,
478 stacklevel=find_outside_stacklevel("lsst.meas.base"),
479 )
480 row = finalizedPsfApCorrCatalog.find(detectorId)
481 if row is None:
482 self.log.warning("Detector id %s not found in finalizedPsfApCorrCatalog; "
483 "Using original psf.", detectorId)
484 else:
485 psf = row.getPsf()
486 apCorrMap = row.getApCorrMap()
487 if psf is None or apCorrMap is None:
488 self.log.warning("Detector id %s has None for psf/apCorrMap in "
489 "finalizedPsfApCorrCatalog; Using original psf.", detectorId)
490 else:
491 exposure.setPsf(psf)
492 exposure.setApCorrMap(apCorrMap)
493
494 if skyCorr is not None:
495 exposure.maskedImage -= skyCorr.getImage()
496
497 return exposure
498
499 def mergeAndFilterReferences(self, exposure, refCats, refWcs):
500 """Filter reference catalog so that all sources are within the
501 boundaries of the exposure.
502
503 Parameters
504 ----------
505 exposure : `lsst.afw.image.exposure.Exposure`
506 Exposure to generate the catalog for.
507 refCats : sequence of `lsst.daf.butler.DeferredDatasetHandle`
508 Handles for catalogs of shapes and positions at which to force
509 photometry.
510 refWcs : `lsst.afw.image.SkyWcs`
511 Reference world coordinate system.
512
513 Returns
514 -------
516 Filtered catalog of forced sources to measure.
517
518 Notes
519 -----
520 The majority of this code is based on the methods of
522
523 """
524 mergedRefCat = None
525
526 # Step 1: Determine bounds of the exposure photometry will
527 # be performed on.
528 expWcs = exposure.getWcs()
529 if expWcs is None:
530 self.log.info("Exposure has no WCS. Returning None for mergedRefCat.")
531 else:
532 expRegion = exposure.getBBox(lsst.afw.image.PARENT)
533 expBBox = lsst.geom.Box2D(expRegion)
534 expBoxCorners = expBBox.getCorners()
535 expSkyCorners = [expWcs.pixelToSky(corner).getVector() for
536 corner in expBoxCorners]
537 expPolygon = lsst.sphgeom.ConvexPolygon(expSkyCorners)
538
539 # Step 2: Filter out reference catalog sources that are
540 # not contained within the exposure boundaries, or whose
541 # parents are not within the exposure boundaries. Note
542 # that within a single input refCat, the parents always
543 # appear before the children.
544 for refCat in refCats:
545 refCat = refCat.get()
546 if mergedRefCat is None:
547 mergedRefCat = lsst.afw.table.SourceCatalog(refCat.table)
548 containedIds = {0} # zero as a parent ID means "this is a parent"
549 for record in refCat:
550 if (expPolygon.contains(record.getCoord().getVector()) and record.getParent()
551 in containedIds):
552 record.setFootprint(record.getFootprint())
553 mergedRefCat.append(record)
554 containedIds.add(record.getId())
555 if mergedRefCat is None:
556 raise RuntimeError("No reference objects for forced photometry.")
557 mergedRefCat.sort(lsst.afw.table.SourceTable.getParentKey())
558 return mergedRefCat
559
560 def generateMeasCat(self, dataId, exposure, refCat, refWcs):
561 """Generate a measurement catalog.
562
563 Parameters
564 ----------
565 dataId : `lsst.daf.butler.DataCoordinate`
566 Butler data ID for this image, with ``{visit, detector}`` keys.
567 exposure : `lsst.afw.image.exposure.Exposure`
568 Exposure to generate the catalog for.
570 Catalog of shapes and positions at which to force photometry.
571 refWcs : `lsst.afw.image.SkyWcs`
572 Reference world coordinate system.
573 This parameter is not currently used.
574
575 Returns
576 -------
578 Catalog of forced sources to measure.
579 expId : `int`
580 Unique binary id associated with the input exposure
581 """
582 id_generator = self.config.idGenerator.apply(dataId)
583 measCat = self.measurement.generateMeasCat(exposure, refCat, refWcs,
584 idFactory=id_generator.make_table_id_factory())
585 return measCat, id_generator.catalog_id
586
587 def run(self, measCat, exposure, refCat, refWcs, exposureId=None):
588 """Perform forced measurement on a single exposure.
589
590 Parameters
591 ----------
593 The measurement catalog, based on the sources listed in the
594 reference catalog.
595 exposure : `lsst.afw.image.Exposure`
596 The measurement image upon which to perform forced detection.
598 The reference catalog of sources to measure.
599 refWcs : `lsst.afw.image.SkyWcs`
600 The WCS for the references.
601 exposureId : `int`
602 Optional unique exposureId used for random seed in measurement
603 task.
604
605 Returns
606 -------
607 result : `lsst.pipe.base.Struct`
608 Structure with fields:
609
610 ``measCat``
611 Catalog of forced measurement results
613 """
614 self.measurement.run(measCat, exposure, refCat, refWcs, exposureId=exposureId)
615 if self.config.doApCorr:
616 apCorrMap = exposure.getInfo().getApCorrMap()
617 if apCorrMap is None:
618 self.log.warning("Forced exposure image does not have valid aperture correction; skipping.")
619 else:
620 self.applyApCorr.run(
621 catalog=measCat,
622 apCorrMap=apCorrMap,
623 )
624 self.catalogCalculation.run(measCat)
625
626 return pipeBase.Struct(measCat=measCat)
627
628 def attachFootprints(self, sources, refCat, exposure, refWcs):
629 """Attach footprints to blank sources prior to measurements.
630
631 Notes
632 -----
633 `~lsst.afw.detection.Footprint` objects for forced photometry must
634 be in the pixel coordinate system of the image being measured, while
635 the actual detections may start out in a different coordinate system.
636
637 Subclasses of this class may implement this method to define how
638 those `~lsst.afw.detection.Footprint` objects should be generated.
639
640 This default implementation transforms depends on the
641 ``footprintSource`` configuration parameter.
642 """
643 if self.config.footprintSource == "transformed":
644 return self.measurement.attachTransformedFootprints(sources, refCat, exposure, refWcs)
645 elif self.config.footprintSource == "psf":
646 return self.measurement.attachPsfShapeFootprints(sources, exposure,
647 scaling=self.config.psfFootprintScaling)
648
649
650class ForcedPhotCcdFromDataFrameConnections(PipelineTaskConnections,
651 dimensions=("instrument", "visit", "detector", "skymap", "tract"),
652 defaultTemplates={"inputCoaddName": "goodSeeing",
653 "inputName": "calexp",
654 "skyWcsName": "gbdesAstrometricFit",
655 "photoCalibName": "fgcm"},
656 # TODO: remove on DM-39854
657 deprecatedTemplates={
658 "skyWcsName": "Deprecated; will be removed after v26.",
659 "photoCalibName": "Deprecated; will be removed after v26."
660 }):
661 refCat = cT.Input(
662 doc="Catalog of positions at which to force photometry.",
663 name="{inputCoaddName}Diff_fullDiaObjTable",
664 storageClass="DataFrame",
665 dimensions=["skymap", "tract", "patch"],
666 multiple=True,
667 deferLoad=True,
668 )
669 exposure = cT.Input(
670 doc="Input exposure to perform photometry on.",
671 name="{inputName}",
672 storageClass="ExposureF",
673 dimensions=["instrument", "visit", "detector"],
674 )
675 skyCorr = cT.Input(
676 doc="Input Sky Correction to be subtracted from the calexp if doApplySkyCorr=True",
677 name="skyCorr",
678 storageClass="Background",
679 dimensions=("instrument", "visit", "detector"),
680 )
681 visitSummary = cT.Input(
682 doc="Input visit-summary catalog with updated calibration objects.",
683 name="finalVisitSummary",
684 storageClass="ExposureCatalog",
685 dimensions=("instrument", "visit"),
686 )
687 externalSkyWcsTractCatalog = cT.Input(
688 doc=("Per-tract, per-visit wcs calibrations. These catalogs use the detector "
689 "id for the catalog id, sorted on id for fast lookup."),
690 name="{skyWcsName}SkyWcsCatalog",
691 storageClass="ExposureCatalog",
692 dimensions=["instrument", "visit", "tract"],
693 # TODO: remove on DM-39854
694 deprecated="Deprecated in favor of 'visitSummary'; will be removed after v26."
695 )
696 externalSkyWcsGlobalCatalog = cT.Input(
697 doc=("Per-visit wcs calibrations computed globally (with no tract information). "
698 "These catalogs use the detector id for the catalog id, sorted on id for "
699 "fast lookup."),
700 name="{skyWcsName}SkyWcsCatalog",
701 storageClass="ExposureCatalog",
702 dimensions=["instrument", "visit"],
703 # TODO: remove on DM-39854
704 deprecated="Deprecated in favor of 'visitSummary'; will be removed after v26."
705 )
706 externalPhotoCalibTractCatalog = cT.Input(
707 doc=("Per-tract, per-visit photometric calibrations. These catalogs use the "
708 "detector id for the catalog id, sorted on id for fast lookup."),
709 name="{photoCalibName}PhotoCalibCatalog",
710 storageClass="ExposureCatalog",
711 dimensions=["instrument", "visit", "tract"],
712 # TODO: remove on DM-39854
713 deprecated="Deprecated in favor of 'visitSummary'; will be removed after v26."
714 )
715 externalPhotoCalibGlobalCatalog = cT.Input(
716 doc=("Per-visit photometric calibrations computed globally (with no tract "
717 "information). These catalogs use the detector id for the catalog id, "
718 "sorted on id for fast lookup."),
719 name="{photoCalibName}PhotoCalibCatalog",
720 storageClass="ExposureCatalog",
721 dimensions=["instrument", "visit"],
722 # TODO: remove on DM-39854
723 deprecated="Deprecated in favor of 'visitSummary'; will be removed after v26."
724 )
725 finalizedPsfApCorrCatalog = cT.Input(
726 doc=("Per-visit finalized psf models and aperture correction maps. "
727 "These catalogs use the detector id for the catalog id, "
728 "sorted on id for fast lookup."),
729 name="finalized_psf_ap_corr_catalog",
730 storageClass="ExposureCatalog",
731 dimensions=["instrument", "visit"],
732 # TODO: remove on DM-39854
733 deprecated="Deprecated in favor of 'visitSummary'; will be removed after v26."
734 )
735 measCat = cT.Output(
736 doc="Output forced photometry catalog.",
737 name="forced_src_diaObject",
738 storageClass="SourceCatalog",
739 dimensions=["instrument", "visit", "detector", "skymap", "tract"],
740 )
741 outputSchema = cT.InitOutput(
742 doc="Schema for the output forced measurement catalogs.",
743 name="forced_src_diaObject_schema",
744 storageClass="SourceCatalog",
745 )
746
747 def __init__(self, *, config=None):
748 super().__init__(config=config)
749 if not config.doApplySkyCorr:
750 self.inputs.remove("skyCorr")
751 if config.doApplyExternalSkyWcs:
752 if config.useGlobalExternalSkyWcs:
753 self.inputs.remove("externalSkyWcsTractCatalog")
754 else:
755 self.inputs.remove("externalSkyWcsGlobalCatalog")
756 else:
757 self.inputs.remove("externalSkyWcsTractCatalog")
758 self.inputs.remove("externalSkyWcsGlobalCatalog")
759 if config.doApplyExternalPhotoCalib:
760 if config.useGlobalExternalPhotoCalib:
761 self.inputs.remove("externalPhotoCalibTractCatalog")
762 else:
763 self.inputs.remove("externalPhotoCalibGlobalCatalog")
764 else:
765 self.inputs.remove("externalPhotoCalibTractCatalog")
766 self.inputs.remove("externalPhotoCalibGlobalCatalog")
767 if not config.doApplyFinalizedPsf:
768 self.inputs.remove("finalizedPsfApCorrCatalog")
769
770
771class ForcedPhotCcdFromDataFrameConfig(ForcedPhotCcdConfig,
772 pipelineConnections=ForcedPhotCcdFromDataFrameConnections):
773 def setDefaults(self):
774 super().setDefaults()
775 self.footprintSource = "psf"
776 self.measurement.doReplaceWithNoise = False
777 # Only run a minimal set of plugins, as these measurements are only
778 # needed for PSF-like sources.
779 self.measurement.plugins.names = ["base_PixelFlags",
780 "base_TransformedCentroidFromCoord",
781 "base_PsfFlux",
782 "base_LocalBackground",
783 "base_LocalPhotoCalib",
784 "base_LocalWcs",
785 ]
786 self.measurement.slots.shape = None
787 # Make catalogCalculation a no-op by default as no modelFlux is setup
788 # by default in ForcedMeasurementTask.
789 self.catalogCalculation.plugins.names = []
790
791 self.measurement.copyColumns = {'id': 'diaObjectId', 'coord_ra': 'coord_ra', 'coord_dec': 'coord_dec'}
792 self.measurement.slots.centroid = "base_TransformedCentroidFromCoord"
793 self.measurement.slots.psfFlux = "base_PsfFlux"
794
795 def validate(self):
796 super().validate()
797 if self.footprintSource == "transformed":
798 raise ValueError("Cannot transform footprints from reference catalog, "
799 "because DataFrames can't hold footprints.")
800
801
802class ForcedPhotCcdFromDataFrameTask(ForcedPhotCcdTask):
803 """Force Photometry on a per-detector exposure with coords from a DataFrame
804
805 Uses input from a DataFrame instead of SourceCatalog
806 like the base class ForcedPhotCcd does.
807 Writes out a SourceCatalog so that the downstream
808 WriteForcedSourceTableTask can be reused with output from this Task.
809 """
810 _DefaultName = "forcedPhotCcdFromDataFrame"
811 ConfigClass = ForcedPhotCcdFromDataFrameConfig
812
813 def __init__(self, refSchema=None, initInputs=None, **kwargs):
814 # Parent's init assumes that we have a reference schema; Cannot reuse
815 pipeBase.PipelineTask.__init__(self, **kwargs)
816
817 self.makeSubtask("measurement", refSchema=lsst.afw.table.SourceTable.makeMinimalSchema())
818
819 if self.config.doApCorr:
820 self.makeSubtask("applyApCorr", schema=self.measurement.schema)
821 self.makeSubtask('catalogCalculation', schema=self.measurement.schema)
822 self.outputSchema = lsst.afw.table.SourceCatalog(self.measurement.schema)
823
824 def runQuantum(self, butlerQC, inputRefs, outputRefs):
825 inputs = butlerQC.get(inputRefs)
826
827 # When run with dataframes, we do not need a reference wcs.
828 inputs['refWcs'] = None
829
830 # Connections only exist if they are configured to be used.
831 skyCorr = inputs.pop('skyCorr', None)
832 if self.config.useGlobalExternalSkyWcs:
833 externalSkyWcsCatalog = inputs.pop('externalSkyWcsGlobalCatalog', None)
834 else:
835 externalSkyWcsCatalog = inputs.pop('externalSkyWcsTractCatalog', None)
836 if self.config.useGlobalExternalPhotoCalib:
837 externalPhotoCalibCatalog = inputs.pop('externalPhotoCalibGlobalCatalog', None)
838 else:
839 externalPhotoCalibCatalog = inputs.pop('externalPhotoCalibTractCatalog', None)
840 finalizedPsfApCorrCatalog = inputs.pop('finalizedPsfApCorrCatalog', None)
841
842 inputs['exposure'] = self.prepareCalibratedExposure(
843 inputs['exposure'],
844 skyCorr=skyCorr,
845 externalSkyWcsCatalog=externalSkyWcsCatalog,
846 externalPhotoCalibCatalog=externalPhotoCalibCatalog,
847 finalizedPsfApCorrCatalog=finalizedPsfApCorrCatalog,
848 visitSummary=inputs.pop("visitSummary"),
849 )
850
851 self.log.info("Filtering ref cats: %s", ','.join([str(i.dataId) for i in inputs['refCat']]))
852 if inputs["exposure"].getWcs() is not None:
853 refCat = self.df2RefCat([i.get(parameters={"columns": ['diaObjectId', 'ra', 'dec']})
854 for i in inputs['refCat']],
855 inputs['exposure'].getBBox(), inputs['exposure'].getWcs())
856 inputs['refCat'] = refCat
857 # generateMeasCat does not use the refWcs.
858 inputs['measCat'], inputs['exposureId'] = self.generateMeasCat(
859 inputRefs.exposure.dataId, inputs['exposure'], inputs['refCat'], inputs['refWcs']
860 )
861 # attachFootprints only uses refWcs in ``transformed`` mode, which is not
862 # supported in the DataFrame-backed task.
863 self.attachFootprints(inputs["measCat"], inputs["refCat"], inputs["exposure"], inputs["refWcs"])
864 outputs = self.run(**inputs)
865
866 butlerQC.put(outputs, outputRefs)
867 else:
868 self.log.info("No WCS for %s. Skipping and no %s catalog will be written.",
869 butlerQC.quantum.dataId, outputRefs.measCat.datasetType.name)
870
871 def df2RefCat(self, dfList, exposureBBox, exposureWcs):
872 """Convert list of DataFrames to reference catalog
873
874 Concatenate list of DataFrames presumably from multiple patches and
875 downselect rows that overlap the exposureBBox using the exposureWcs.
876
877 Parameters
878 ----------
879 dfList : `list` of `pandas.DataFrame`
880 Each element containst diaObjects with ra/dec position in degrees
881 Columns 'diaObjectId', 'ra', 'dec' are expected
882 exposureBBox : `lsst.geom.Box2I`
883 Bounding box on which to select rows that overlap
884 exposureWcs : `lsst.afw.geom.SkyWcs`
885 World coordinate system to convert sky coords in ref cat to
886 pixel coords with which to compare with exposureBBox
887
888 Returns
889 -------
891 Source Catalog with minimal schema that overlaps exposureBBox
892 """
893 df = pd.concat(dfList)
894 # translate ra/dec coords in dataframe to detector pixel coords
895 # to down select rows that overlap the detector bbox
896 mapping = exposureWcs.getTransform().getMapping()
897 x, y = mapping.applyInverse(np.array(df[['ra', 'dec']].values*2*np.pi/360).T)
898 inBBox = lsst.geom.Box2D(exposureBBox).contains(x, y)
899 refCat = self.df2SourceCat(df[inBBox])
900 return refCat
901
902 def df2SourceCat(self, df):
903 """Create minimal schema SourceCatalog from a pandas DataFrame.
904
905 The forced measurement subtask expects this as input.
906
907 Parameters
908 ----------
909 df : `pandas.DataFrame`
910 DiaObjects with locations and ids.
911
912 Returns
913 -------
914 outputCatalog : `lsst.afw.table.SourceTable`
915 Output catalog with minimal schema.
916 """
918 outputCatalog = lsst.afw.table.SourceCatalog(schema)
919 outputCatalog.reserve(len(df))
920
921 for diaObjectId, ra, dec in df[['ra', 'dec']].itertuples():
922 outputRecord = outputCatalog.addNew()
923 outputRecord.setId(diaObjectId)
924 outputRecord.setCoord(lsst.geom.SpherePoint(ra, dec, lsst.geom.degrees))
925 return outputCatalog
table::Key< int > to
Class to describe the properties of a detected object from an image.
Definition Footprint.h:63
A 2-dimensional celestial WCS that transform pixels to ICRS RA/Dec, using the LSST standard for pixel...
Definition SkyWcs.h:117
A class to contain the data, WCS, and other information needed to describe an image of the sky.
Definition Exposure.h:72
Custom catalog class for ExposureRecord/Table.
Definition Exposure.h:311
Defines the fields and offsets for a table.
Definition Schema.h:51
Table class that contains measurements made on a single exposure.
Definition Source.h:217
static Schema makeMinimalSchema()
Return a minimal schema for Source tables and records.
Definition Source.h:258
static Key< RecordId > getParentKey()
Key for the parent ID.
Definition Source.h:273
A floating-point coordinate rectangle geometry.
Definition Box.h:413
An integer coordinate rectangle.
Definition Box.h:55
Point in an unspecified spherical coordinate system.
Definition SpherePoint.h:57
ConvexPolygon is a closed convex polygon on the unit sphere.
df2RefCat(self, dfList, exposureBBox, exposureWcs)
runQuantum(self, butlerQC, inputRefs, outputRefs)