Coverage for python / lsst / analysis / tools / tasks / photometricCatalogMatch.py: 23%
99 statements
« prev ^ index » next coverage.py v7.13.5, created at 2026-04-21 10:55 +0000
1# This file is part of analysis_tools.
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/>.
22__all__ = ("PhotometricCatalogMatchConfig", "PhotometricCatalogMatchTask")
25import numpy as np
26from astropy.table import Table
28import lsst.geom
29import lsst.pipe.base as pipeBase
30from lsst.pipe.tasks.loadReferenceCatalog import LoadReferenceCatalogTask
32from ..actions.vector import StarSelector, VisitPlotFlagSelector
33from ..tasks.catalogMatch import CatalogMatchConfig, CatalogMatchConnections, CatalogMatchTask
36class PhotometricCatalogMatchConnections(
37 CatalogMatchConnections,
38 dimensions=("tract", "skymap"),
39 defaultTemplates={"targetCatalog": "objectTable_tract", "refCatalog": "ps1_pv3_3pi_20170110"},
40):
41 matchedCatalog = pipeBase.connectionTypes.Output(
42 doc="Catalog with matched target and reference objects with separations",
43 name="{targetCatalog}_{refCatalog}_match_photom",
44 storageClass="ArrowAstropy",
45 dimensions=("tract", "skymap"),
46 )
49class PhotometricCatalogMatchConfig(
50 CatalogMatchConfig,
51 pipelineConnections=PhotometricCatalogMatchConnections,
52):
53 def setDefaults(self):
54 super().setDefaults()
55 self.matchesRefCat = True
56 self.referenceCatalogLoader.doReferenceSelection = False
57 self.referenceCatalogLoader.doApplyColorTerms = True
60class PhotometricCatalogMatchTask(CatalogMatchTask):
61 """A wrapper task to provide the information that
62 is specific to the photometric reference catalog.
63 """
65 ConfigClass = PhotometricCatalogMatchConfig
66 _DefaultName = "analysisToolsPhotometricCatalogMatch"
68 def runQuantum(self, butlerQC, inputRefs, outputRefs):
69 """Run the matching to the photometric reference catalog.
71 Parameters
72 ----------
73 butlerQC : `lsst.pipe.base.QuantumContext`
74 inputRefs : `lsst.pipe.base.InputQuantizedConnection`
75 outputRefs : `lsst.pipe.base.OutputQuantizedConnection`
76 """
78 inputs = butlerQC.get(inputRefs)
79 bands = []
80 for filterName in self.config.filterNames:
81 bands.append(self.config.referenceCatalogLoader.refObjLoader.filterMap[filterName])
83 # For some reason the imsim filterMaps don't work the same way
84 # as the HSC ones do, this is a bit hacky but fixes this.
85 # This code makes the assumption that filterMap is a dict
86 # mapping observed filter names: band, but filterMap is
87 # currently defined as a mapping of:
88 # observed filter names -> reference catalog filter names.
89 # Therefore, when the reference catalog filter name is
90 # not the band name we need this work around.
91 # TODO: workaround for DM-46728
92 if bands[0].startswith("lsst") or "sim" in bands[0] or "smeared" in bands[0]:
93 bands = self.config.filterNames
94 elif bands[0].startswith("monster"):
95 # for the_monster_* the reference catalog filter name is
96 # "monster_{system}_{band}" the last character is the band
97 bands = [band[-1] for band in bands]
99 columns = self.prepColumns(bands)
100 table = inputs["catalog"].get(parameters={"columns": columns})
102 tract = butlerQC.quantum.dataId["tract"]
104 loaderTask = LoadReferenceCatalogTask(
105 config=self.config.referenceCatalogLoader,
106 dataIds=[ref.dataId for ref in inputRefs.refCat],
107 name=inputs["refCat"][0].ref.datasetType.name,
108 refCats=inputs["refCat"],
109 )
111 skymap = inputs.pop("skymap")
112 loadedRefCat = self._loadRefCat(loaderTask, skymap[tract])
113 outputs = self.run(targetCatalog=table, refCatalog=loadedRefCat, bands=bands)
115 butlerQC.put(outputs, outputRefs)
118class PhotometricCatalogMatchVisitConnections(
119 pipeBase.PipelineTaskConnections,
120 dimensions=("visit",),
121 defaultTemplates={"targetCatalog": "sourceTable_visit", "refCatalog": "ps1_pv3_3pi_20170110"},
122):
123 catalog = pipeBase.connectionTypes.Input(
124 doc="The visit-wide catalog to make plots from.",
125 storageClass="ArrowAstropy",
126 name="{targetCatalog}",
127 dimensions=("visit",),
128 deferLoad=True,
129 )
131 refCat = pipeBase.connectionTypes.PrerequisiteInput(
132 doc="The photometric reference catalog to match to.",
133 name="{refCatalog}",
134 storageClass="SimpleCatalog",
135 dimensions=("skypix",),
136 deferLoad=True,
137 multiple=True,
138 )
140 visitSummaryTable = pipeBase.connectionTypes.Input(
141 doc="A summary table of the ccds in the visit",
142 storageClass="ExposureCatalog",
143 name="finalVisitSummary",
144 dimensions=("visit",),
145 )
147 matchedCatalog = pipeBase.connectionTypes.Output(
148 doc="Catalog with matched target and reference objects with separations",
149 name="{targetCatalog}_{refCatalog}_match_photom",
150 storageClass="ArrowAstropy",
151 dimensions=("visit",),
152 )
155class PhotometricCatalogMatchVisitConfig(
156 PhotometricCatalogMatchConfig, pipelineConnections=PhotometricCatalogMatchVisitConnections
157):
158 def setDefaults(self):
159 super().setDefaults()
160 self.idColumn = "sourceId"
161 self.matchesRefCat = True
162 self.filterNames = []
163 self.extraPerBandColumns = []
164 self.patchColumn = ""
165 self.selectorBands = []
166 self.extraColumns = ["x", "y", "ap09Flux", "ap09FluxErr"]
167 self.selectorActions.flagSelector = VisitPlotFlagSelector
168 self.sourceSelectorActions.sourceSelector = StarSelector()
169 self.sourceSelectorActions.sourceSelector.vectorKey = "extendedness"
170 self.extraColumnSelectors.selector1.fluxType = "psfFlux"
171 self.extraColumnSelectors.selector2.vectorKey = "extendedness"
172 self.extraColumnSelectors.selector3.vectorKey = "extendedness"
173 self.extraColumnSelectors.selector4 = VisitPlotFlagSelector
176class PhotometricCatalogMatchVisitTask(PhotometricCatalogMatchTask):
177 """A wrapper task to provide the information that
178 is specific to the photometric reference catalog.
179 """
181 ConfigClass = PhotometricCatalogMatchVisitConfig
182 _DefaultName = "analysisToolsPhotometricCatalogMatchVisit"
184 def runQuantum(self, butlerQC, inputRefs, outputRefs):
185 """Run the matching to the photometric reference catalog.
187 Parameters
188 ----------
189 butlerQC : `lsst.pipe.base.QuantumContext`
190 inputRefs : `lsst.pipe.base.InputQuantizedConnection`
191 outputRefs : `lsst.pipe.base.OutputQuantizedConnection`
192 """
194 inputs = butlerQC.get(inputRefs)
195 physicalFilter = inputs["catalog"].dataId["physical_filter"]
197 # For some reason the imsim filterMaps don't work the same way as
198 # the HSC ones do, this is a bit hacky but fixes this
199 # This code makes the assumption that filterMap is a dict
200 # mapping observed filter names: band, but filterMap is
201 # currently defined as a mapping of:
202 # observed filter names -> reference catalog filter names.
203 # Therefore, when the reference catalog filter name is
204 # not the band name we need this work around.
205 # TODO: workaround for DM-46728
206 if "sim" in physicalFilter:
207 physicalFilter = physicalFilter[0]
208 bands = [physicalFilter]
209 elif self.config.referenceCatalogLoader.refObjLoader.filterMap[physicalFilter].startswith("monster"):
210 # filtermap gives the column name for the refcat
211 # which is monster_{system}_{band}
212 bands = [self.config.referenceCatalogLoader.refObjLoader.filterMap[physicalFilter][-1]]
213 else:
214 bands = [self.config.referenceCatalogLoader.refObjLoader.filterMap[physicalFilter]]
216 # No bands needed for visit tables
217 # but we do need them later for the matching
218 columns = [self.config.idColumn, "coord_ra", "coord_dec", "detector"]
219 columns.extend(self.config.extraColumns)
220 for selectorAction in [
221 self.config.selectorActions,
222 self.config.sourceSelectorActions,
223 self.config.extraColumnSelectors,
224 ]:
225 for selector in selectorAction:
226 selectorSchema = selector.getFormattedInputSchema()
227 columns += [s[0] for s in selectorSchema]
229 table = inputs["catalog"].get(parameters={"columns": columns})
231 loaderTask = LoadReferenceCatalogTask(
232 config=self.config.referenceCatalogLoader,
233 dataIds=[ref.dataId for ref in inputRefs.refCat],
234 name=inputs["refCat"][0].ref.datasetType.name,
235 refCats=inputs["refCat"],
236 )
238 visitSummaryTable = inputs.pop("visitSummaryTable")
239 loadedRefCat = self._loadRefCat(loaderTask, visitSummaryTable, physicalFilter)
240 outputs = self.run(targetCatalog=table, refCatalog=loadedRefCat, bands=bands)
242 # The matcher adds the band to the front of the columns
243 # but the visit plots aren't expecting it
244 cols = list(outputs.matchedCatalog.columns)
245 for col in cols:
246 if col[:2] == bands[0] + "_":
247 outputs.matchedCatalog.rename_column(col, col[2:])
249 butlerQC.put(outputs, outputRefs)
251 def _loadRefCat(self, loaderTask, visitSummaryTable, physicalFilter):
252 """Make a reference catalog with coordinates in degrees
254 Parameters
255 ----------
256 visitSummaryTable : `lsst.afw.table.ExposureCatalog`
257 The table of visit information
258 """
259 # Get convex hull around the detectors, then get its center and radius
260 corners = []
261 for visSum in visitSummaryTable:
262 for ra, dec in zip(visSum["raCorners"], visSum["decCorners"]):
263 # If the coordinates are nan then don't keep going
264 # because it crashes later
265 if not np.isfinite(ra) or not np.isfinite(dec):
266 raise pipeBase.NoWorkFound("Visit summary corners not finite")
267 corners.append(lsst.geom.SpherePoint(ra, dec, units=lsst.geom.degrees).getVector())
268 visitBoundingCircle = lsst.sphgeom.ConvexPolygon.convexHull(corners).getBoundingCircle()
269 center = lsst.geom.SpherePoint(visitBoundingCircle.getCenter())
270 radius = visitBoundingCircle.getOpeningAngle()
272 # Get the observation date of the visit
273 epoch = visSum.getVisitInfo().getDate().toAstropy()
275 # Load the reference catalog in the skyCircle of the detectors, then
276 # convert the coordinates to degrees and convert the catalog to a
277 # dataframe
279 loadedRefCat = loaderTask.getSkyCircleCatalog(center, radius, [physicalFilter], epoch=epoch)
281 return Table(loadedRefCat)