LSST Applications 27.0.0,g0265f82a02+469cd937ee,g02d81e74bb+21ad69e7e1,g1470d8bcf6+cbe83ee85a,g2079a07aa2+e67c6346a6,g212a7c68fe+04a9158687,g2305ad1205+94392ce272,g295015adf3+81dd352a9d,g2bbee38e9b+469cd937ee,g337abbeb29+469cd937ee,g3939d97d7f+72a9f7b576,g487adcacf7+71499e7cba,g50ff169b8f+5929b3527e,g52b1c1532d+a6fc98d2e7,g591dd9f2cf+df404f777f,g5a732f18d5+be83d3ecdb,g64a986408d+21ad69e7e1,g858d7b2824+21ad69e7e1,g8a8a8dda67+a6fc98d2e7,g99cad8db69+f62e5b0af5,g9ddcbc5298+d4bad12328,ga1e77700b3+9c366c4306,ga8c6da7877+71e4819109,gb0e22166c9+25ba2f69a1,gb6a65358fc+469cd937ee,gbb8dafda3b+69d3c0e320,gc07e1c2157+a98bf949bb,gc120e1dc64+615ec43309,gc28159a63d+469cd937ee,gcf0d15dbbd+72a9f7b576,gdaeeff99f8+a38ce5ea23,ge6526c86ff+3a7c1ac5f1,ge79ae78c31+469cd937ee,gee10cc3b42+a6fc98d2e7,gf1cff7945b+21ad69e7e1,gfbcc870c63+9a11dc8c8f
LSST Data Management Base Package
Loading...
Searching...
No Matches
matchFakes.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__ = ["MatchFakesTask",
23 "MatchFakesConfig",
24 "MatchVariableFakesConfig",
25 "MatchVariableFakesTask"]
26
27import astropy.units as u
28import numpy as np
29import pandas as pd
30from scipy.spatial import cKDTree
31
32from lsst.geom import Box2D
33import lsst.pex.config as pexConfig
34from lsst.pipe.base import PipelineTask, PipelineTaskConnections, Struct
35import lsst.pipe.base.connectionTypes as connTypes
36from lsst.skymap import BaseSkyMap
37
38from lsst.pipe.tasks.insertFakes import InsertFakesConfig
39
40
41class MatchFakesConnections(PipelineTaskConnections,
42 defaultTemplates={"coaddName": "deep",
43 "fakesType": "fakes_"},
44 dimensions=("instrument",
45 "visit",
46 "detector")):
47 skyMap = connTypes.Input(
48 doc="Input definition of geometry/bbox and projection/wcs for "
49 "template exposures. Needed to test which tract to generate ",
50 name=BaseSkyMap.SKYMAP_DATASET_TYPE_NAME,
51 dimensions=("skymap",),
52 storageClass="SkyMap",
53 )
54 fakeCats = connTypes.Input(
55 doc="Catalog of fake sources inserted into an image.",
56 name="{fakesType}fakeSourceCat",
57 storageClass="DataFrame",
58 dimensions=("tract", "skymap"),
59 deferLoad=True,
60 multiple=True
61 )
62 diffIm = connTypes.Input(
63 doc="Difference image on which the DiaSources were detected.",
64 name="{fakesType}{coaddName}Diff_differenceExp",
65 storageClass="ExposureF",
66 dimensions=("instrument", "visit", "detector"),
67 )
68 associatedDiaSources = connTypes.Input(
69 doc="A DiaSource catalog to match against fakeCat. Assumed "
70 "to be SDMified.",
71 name="{fakesType}{coaddName}Diff_assocDiaSrc",
72 storageClass="DataFrame",
73 dimensions=("instrument", "visit", "detector"),
74 )
75 matchedDiaSources = connTypes.Output(
76 doc="A catalog of those fakeCat sources that have a match in "
77 "associatedDiaSources. The schema is the union of the schemas for "
78 "``fakeCat`` and ``associatedDiaSources``.",
79 name="{fakesType}{coaddName}Diff_matchDiaSrc",
80 storageClass="DataFrame",
81 dimensions=("instrument", "visit", "detector"),
82 )
83
84
85class MatchFakesConfig(
86 InsertFakesConfig,
87 pipelineConnections=MatchFakesConnections):
88 """Config for MatchFakesTask.
89 """
90 matchDistanceArcseconds = pexConfig.RangeField(
91 doc="Distance in arcseconds to match within.",
92 dtype=float,
93 default=0.5,
94 min=0,
95 max=10,
96 )
97
98 doMatchVisit = pexConfig.Field(
99 dtype=bool,
100 default=False,
101 doc="Match visit to trim the fakeCat"
102 )
103
104 trimBuffer = pexConfig.Field(
105 doc="Size of the pixel buffer surrounding the image. Only those fake sources with a centroid"
106 "falling within the image+buffer region will be considered matches.",
107 dtype=int,
108 default=100,
109 )
110
111
112class MatchFakesTask(PipelineTask):
113 """Match a pre-existing catalog of fakes to a catalog of detections on
114 a difference image.
115
116 This task is generally for injected sources that cannot be easily
117 identified by their footprints such as in the case of detector sources
118 post image differencing.
119 """
120
121 _DefaultName = "matchFakes"
122 ConfigClass = MatchFakesConfig
123
124 def run(self, fakeCats, skyMap, diffIm, associatedDiaSources):
125 """Compose fakes into a single catalog and match fakes to detected
126 diaSources within a difference image bound.
127
128 Parameters
129 ----------
130 fakeCats : `pandas.DataFrame`
131 List of catalog of fakes to match to detected diaSources.
132 skyMap : `lsst.skymap.SkyMap`
133 SkyMap defining the tracts and patches the fakes are stored over.
134 diffIm : `lsst.afw.image.Exposure`
135 Difference image where ``associatedDiaSources`` were detected.
136 associatedDiaSources : `pandas.DataFrame`
137 Catalog of difference image sources detected in ``diffIm``.
138
139 Returns
140 -------
141 result : `lsst.pipe.base.Struct`
142 Results struct with components.
143
144 - ``matchedDiaSources`` : Fakes matched to input diaSources. Has
145 length of ``fakeCat``. (`pandas.DataFrame`)
146 """
147 fakeCat = self.composeFakeCat(fakeCats, skyMap)
148
149 if self.config.doMatchVisit:
150 fakeCat = self.getVisitMatchedFakeCat(fakeCat, diffIm)
151
152 return self._processFakes(fakeCat, diffIm, associatedDiaSources)
153
154 def _processFakes(self, fakeCat, diffIm, associatedDiaSources):
155 """Match fakes to detected diaSources within a difference image bound.
156
157 Parameters
158 ----------
159 fakeCat : `pandas.DataFrame`
160 Catalog of fakes to match to detected diaSources.
161 diffIm : `lsst.afw.image.Exposure`
162 Difference image where ``associatedDiaSources`` were detected.
163 associatedDiaSources : `pandas.DataFrame`
164 Catalog of difference image sources detected in ``diffIm``.
165
166 Returns
167 -------
168 result : `lsst.pipe.base.Struct`
169 Results struct with components.
170
171 - ``matchedDiaSources`` : Fakes matched to input diaSources. Has
172 length of ``fakeCat``. (`pandas.DataFrame`)
173 """
174 trimmedFakes = self._trimFakeCat(fakeCat, diffIm)
175 nPossibleFakes = len(trimmedFakes)
176
177 fakeVects = self._getVectors(trimmedFakes[self.config.ra_col],
178 trimmedFakes[self.config.dec_col])
179 diaSrcVects = self._getVectors(
180 np.radians(associatedDiaSources.loc[:, "ra"]),
181 np.radians(associatedDiaSources.loc[:, "dec"]))
182
183 diaSrcTree = cKDTree(diaSrcVects)
184 dist, idxs = diaSrcTree.query(
185 fakeVects,
186 distance_upper_bound=np.radians(self.config.matchDistanceArcseconds / 3600))
187 nFakesFound = np.isfinite(dist).sum()
188
189 self.log.info("Found %d out of %d possible.", nFakesFound, nPossibleFakes)
190 diaSrcIds = associatedDiaSources.iloc[np.where(np.isfinite(dist), idxs, 0)]["diaSourceId"].to_numpy()
191 matchedFakes = trimmedFakes.assign(diaSourceId=np.where(np.isfinite(dist), diaSrcIds, 0))
192
193 return Struct(
194 matchedDiaSources=matchedFakes.merge(
195 associatedDiaSources.reset_index(drop=True), on="diaSourceId", how="left")
196 )
197
198 def composeFakeCat(self, fakeCats, skyMap):
199 """Concatenate the fakeCats from tracts that may cover the exposure.
200
201 Parameters
202 ----------
203 fakeCats : `list` of `lsst.daf.butler.DeferredDatasetHandle`
204 Set of fake cats to concatenate.
205 skyMap : `lsst.skymap.SkyMap`
206 SkyMap defining the geometry of the tracts and patches.
207
208 Returns
209 -------
210 combinedFakeCat : `pandas.DataFrame`
211 All fakes that cover the inner polygon of the tracts in this
212 quantum.
213 """
214 if len(fakeCats) == 1:
215 return fakeCats[0].get()
216 outputCat = []
217 for fakeCatRef in fakeCats:
218 cat = fakeCatRef.get()
219 tractId = fakeCatRef.dataId["tract"]
220 # Make sure all data is within the inner part of the tract.
221 outputCat.append(cat[
222 skyMap.findTractIdArray(cat[self.config.ra_col],
223 cat[self.config.dec_col],
224 degrees=False)
225 == tractId])
226
227 return pd.concat(outputCat)
228
229 def getVisitMatchedFakeCat(self, fakeCat, exposure):
230 """Trim the fakeCat to select particular visit
231
232 Parameters
233 ----------
234 fakeCat : `pandas.core.frame.DataFrame`
235 The catalog of fake sources to add to the exposure
236 exposure : `lsst.afw.image.exposure.exposure.ExposureF`
237 The exposure to add the fake sources to
238
239 Returns
240 -------
241 movingFakeCat : `pandas.DataFrame`
242 All fakes that belong to the visit
243 """
244 selected = exposure.getInfo().getVisitInfo().getId() == fakeCat["visit"]
245
246 return fakeCat[selected]
247
248 def _addPixCoords(self, fakeCat, image):
249
250 """Add pixel coordinates to the catalog of fakes.
251
252 Parameters
253 ----------
254 fakeCat : `pandas.core.frame.DataFrame`
255 The catalog of fake sources to be input
256 image : `lsst.afw.image.exposure.exposure.ExposureF`
257 The image into which the fake sources should be added
258 Returns
259 -------
260 fakeCat : `pandas.core.frame.DataFrame`
261 """
262 wcs = image.getWcs()
263 ras = fakeCat[self.config.ra_col].values
264 decs = fakeCat[self.config.dec_col].values
265 xs, ys = wcs.skyToPixelArray(ras, decs)
266 fakeCat["x"] = xs
267 fakeCat["y"] = ys
268
269 return fakeCat
270
271 def _trimFakeCat(self, fakeCat, image):
272 """Trim the fake cat to the exact size of the input image.
273
274 Parameters
275 ----------
276 fakeCat : `pandas.core.frame.DataFrame`
277 The catalog of fake sources that was input
278 image : `lsst.afw.image.exposure.exposure.ExposureF`
279 The image into which the fake sources were added
280 Returns
281 -------
282 fakeCat : `pandas.core.frame.DataFrame`
283 The original fakeCat trimmed to the area of the image
284 """
285
286 # fakeCat must be processed with _addPixCoords before trimming
287 if ('x' not in fakeCat.columns) or ('y' not in fakeCat.columns):
288 fakeCat = self._addPixCoords(fakeCat, image)
289
290 # Prefilter in ra/dec to avoid cases where the wcs incorrectly maps
291 # input fakes which are really off the chip onto it.
292 ras = fakeCat[self.config.ra_col].values * u.rad
293 decs = fakeCat[self.config.dec_col].values * u.rad
294
295 isContainedRaDec = image.containsSkyCoords(ras, decs, padding=0)
296
297 # now use the exact pixel BBox to filter to only fakes that were inserted
298 xs = fakeCat["x"].values
299 ys = fakeCat["y"].values
300
301 bbox = Box2D(image.getBBox())
302 isContainedXy = xs >= bbox.minX
303 isContainedXy &= xs <= bbox.maxX
304 isContainedXy &= ys >= bbox.minY
305 isContainedXy &= ys <= bbox.maxY
306
307 return fakeCat[isContainedRaDec & isContainedXy]
308
309 def _getVectors(self, ras, decs):
310 """Convert ra dec to unit vectors on the sphere.
311
312 Parameters
313 ----------
314 ras : `numpy.ndarray`, (N,)
315 RA coordinates in radians.
316 decs : `numpy.ndarray`, (N,)
317 Dec coordinates in radians.
318
319 Returns
320 -------
321 vectors : `numpy.ndarray`, (N, 3)
322 Vectors on the unit sphere for the given RA/DEC values.
323 """
324 vectors = np.empty((len(ras), 3))
325
326 vectors[:, 2] = np.sin(decs)
327 vectors[:, 0] = np.cos(decs) * np.cos(ras)
328 vectors[:, 1] = np.cos(decs) * np.sin(ras)
329
330 return vectors
331
332
333class MatchVariableFakesConnections(MatchFakesConnections):
334 ccdVisitFakeMagnitudes = connTypes.Input(
335 doc="Catalog of fakes with magnitudes scattered for this ccdVisit.",
336 name="{fakesType}ccdVisitFakeMagnitudes",
337 storageClass="DataFrame",
338 dimensions=("instrument", "visit", "detector"),
339 )
340
341
342class MatchVariableFakesConfig(MatchFakesConfig,
343 pipelineConnections=MatchVariableFakesConnections):
344 """Config for MatchFakesTask.
345 """
346 pass
347
348
349class MatchVariableFakesTask(MatchFakesTask):
350 """Match injected fakes to their detected sources in the catalog and
351 compute their expected brightness in a difference image assuming perfect
352 subtraction.
353
354 This task is generally for injected sources that cannot be easily
355 identified by their footprints such as in the case of detector sources
356 post image differencing.
357 """
358 _DefaultName = "matchVariableFakes"
359 ConfigClass = MatchVariableFakesConfig
360
361 def runQuantum(self, butlerQC, inputRefs, outputRefs):
362 inputs = butlerQC.get(inputRefs)
363 inputs["band"] = butlerQC.quantum.dataId["band"]
364
365 outputs = self.run(**inputs)
366 butlerQC.put(outputs, outputRefs)
367
368 def run(self, fakeCats, ccdVisitFakeMagnitudes, skyMap, diffIm, associatedDiaSources, band):
369 """Match fakes to detected diaSources within a difference image bound.
370
371 Parameters
372 ----------
373 fakeCat : `pandas.DataFrame`
374 Catalog of fakes to match to detected diaSources.
375 diffIm : `lsst.afw.image.Exposure`
376 Difference image where ``associatedDiaSources`` were detected in.
377 associatedDiaSources : `pandas.DataFrame`
378 Catalog of difference image sources detected in ``diffIm``.
379
380 Returns
381 -------
382 result : `lsst.pipe.base.Struct`
383 Results struct with components.
384
385 - ``matchedDiaSources`` : Fakes matched to input diaSources. Has
386 length of ``fakeCat``. (`pandas.DataFrame`)
387 """
388 fakeCat = self.composeFakeCat(fakeCats, skyMap)
389 self.computeExpectedDiffMag(fakeCat, ccdVisitFakeMagnitudes, band)
390 return self._processFakes(fakeCat, diffIm, associatedDiaSources)
391
392 def computeExpectedDiffMag(self, fakeCat, ccdVisitFakeMagnitudes, band):
393 """Compute the magnitude expected in the difference image for this
394 detector/visit. Modify fakeCat in place.
395
396 Negative magnitudes indicate that the source should be detected as
397 a negative source.
398
399 Parameters
400 ----------
401 fakeCat : `pandas.DataFrame`
402 Catalog of fake sources.
403 ccdVisitFakeMagnitudes : `pandas.DataFrame`
404 Magnitudes for variable sources in this specific ccdVisit.
405 band : `str`
406 Band that this ccdVisit was observed in.
407 """
408 magName = self.config.mag_col % band
409 magnitudes = fakeCat[magName].to_numpy()
410 visitMags = ccdVisitFakeMagnitudes["variableMag"].to_numpy()
411 diffFlux = (visitMags * u.ABmag).to_value(u.nJy) - (magnitudes * u.ABmag).to_value(u.nJy)
412 diffMag = np.where(diffFlux > 0,
413 (diffFlux * u.nJy).to_value(u.ABmag),
414 -(-diffFlux * u.nJy).to_value(u.ABmag))
415
416 noVisit = ~fakeCat["isVisitSource"]
417 noTemplate = ~fakeCat["isTemplateSource"]
418 both = np.logical_and(fakeCat["isVisitSource"],
419 fakeCat["isTemplateSource"])
420
421 fakeCat.loc[noVisit, magName] = -magnitudes[noVisit]
422 fakeCat.loc[noTemplate, magName] = visitMags[noTemplate]
423 fakeCat.loc[both, magName] = diffMag[both]
A floating-point coordinate rectangle geometry.
Definition Box.h:413