diff --git a/python/lsst/pipe/tasks/matchDiffimSourceInjected.py b/python/lsst/pipe/tasks/matchDiffimSourceInjected.py
new file mode 100644
index 000000000..d05f606e9
--- /dev/null
+++ b/python/lsst/pipe/tasks/matchDiffimSourceInjected.py
@@ -0,0 +1,430 @@
+# This file is part of ap_pipe.
+#
+# Developed for the LSST Data Management System.
+# This product includes software developed by the LSST Project
+# (https://www.lsst.org).
+# See the COPYRIGHT file at the top-level directory of this distribution
+# for details of code ownership.
+#
+# This program is free software: you can redistribute it and/or modify
+# it under the terms of the GNU General Public License as published by
+# the Free Software Foundation, either version 3 of the License, or
+# (at your option) any later version.
+#
+# This program is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+# GNU General Public License for more details.
+#
+# You should have received a copy of the GNU General Public License
+# along with this program. If not, see .
+
+__all__ = ["MatchInjectedToDiaSourceTask",
+ "MatchInjectedToDiaSourceConfig",
+ "MatchInjectedToAssocDiaSourceTask",
+ "MatchInjectedToAssocDiaSourceConfig"]
+
+import astropy.units as u
+import numpy as np
+from scipy.spatial import cKDTree
+
+from lsst.afw import table as afwTable
+from lsst import geom as lsstGeom
+import lsst.pex.config as pexConfig
+from lsst.pipe.base import PipelineTask, PipelineTaskConfig, PipelineTaskConnections, Struct
+import lsst.pipe.base.connectionTypes as connTypes
+from lsst.meas.base import ForcedMeasurementTask, ForcedMeasurementConfig
+
+
+class MatchInjectedToDiaSourceConnections(
+ PipelineTaskConnections,
+ defaultTemplates={"coaddName": "deep",
+ "fakesType": "fakes_"},
+ dimensions=("instrument",
+ "visit",
+ "detector")):
+ injectedCat = connTypes.Input(
+ doc="Catalog of sources injected in the images.",
+ name="{fakesType}_pvi_catalog",
+ storageClass="ArrowAstropy",
+ dimensions=("instrument", "visit", "detector"),
+ )
+ diffIm = connTypes.Input(
+ doc="Difference image on which the DiaSources were detected.",
+ name="{fakesType}{coaddName}Diff_differenceExp",
+ storageClass="ExposureF",
+ dimensions=("instrument", "visit", "detector"),
+ )
+ diaSources = connTypes.Input(
+ doc="A DiaSource catalog to match against fakeCat.",
+ name="{fakesType}{coaddName}Diff_diaSrc",
+ storageClass="SourceCatalog",
+ dimensions=("instrument", "visit", "detector"),
+ )
+ matchDiaSources = connTypes.Output(
+ doc="A catalog of those fakeCat sources that have a match in "
+ "diaSrc. The schema is the union of the schemas for "
+ "``fakeCat`` and ``diaSrc``.",
+ name="{fakesType}{coaddName}Diff_matchDiaSrc",
+ storageClass="DataFrame",
+ dimensions=("instrument", "visit", "detector"),
+ )
+
+
+class MatchInjectedToDiaSourceConfig(
+ PipelineTaskConfig,
+ pipelineConnections=MatchInjectedToDiaSourceConnections):
+ """Config for MatchFakesTask.
+ """
+ matchDistanceArcseconds = pexConfig.RangeField(
+ doc="Distance in arcseconds to match within.",
+ dtype=float,
+ default=0.5,
+ min=0,
+ max=10,
+ )
+ doMatchVisit = pexConfig.Field(
+ dtype=bool,
+ default=True,
+ doc="Match visit to trim the fakeCat"
+ )
+ trimBuffer = pexConfig.Field(
+ doc="Size of the pixel buffer surrounding the image."
+ "Only those fake sources with a centroid"
+ "falling within the image+buffer region will be considered matches.",
+ dtype=int,
+ default=50,
+ )
+ doForcedMeasurement = pexConfig.Field(
+ dtype=bool,
+ default=True,
+ doc="Force measurement of the fakes at the injection locations."
+ )
+ forcedMeasurement = pexConfig.ConfigurableField(
+ target=ForcedMeasurementTask,
+ doc="Task to force photometer difference image at injection locations.",
+ )
+
+
+class MatchInjectedToDiaSourceTask(PipelineTask):
+
+ _DefaultName = "matchInjectedToDiaSource"
+ ConfigClass = MatchInjectedToDiaSourceConfig
+
+ def run(self, injectedCat, diffIm, diaSources):
+ """Match injected sources to detected diaSources within a difference image bound.
+
+ Parameters
+ ----------
+ injectedCat : `astropy.table.table.Table`
+ Table of catalog of synthetic sources to match to detected diaSources.
+ diffIm : `lsst.afw.image.Exposure`
+ Difference image where ``diaSources`` were detected.
+ diaSources : `afw.table.SourceCatalog`
+ Catalog of difference image sources detected in ``diffIm``.
+ Returns
+ -------
+ result : `lsst.pipe.base.Struct`
+ Results struct with components.
+
+ - ``matchedDiaSources`` : Fakes matched to input diaSources. Has
+ length of ``injectedCalexpCat``. (`pandas.DataFrame`)
+ """
+
+ if self.config.doMatchVisit:
+ fakeCat = self._trimFakeCat(injectedCat, diffIm)
+ else:
+ fakeCat = injectedCat
+ if self.config.doForcedMeasurement:
+ self._estimateFakesSNR(fakeCat, diffIm)
+
+ return self._processFakes(fakeCat, diaSources)
+
+ def _estimateFakesSNR(self, injectedCat, diffIm):
+ """Estimate the signal-to-noise ratio of the fakes in the given catalog.
+
+ Parameters
+ ----------
+ injectedCat : `astropy.table.Table`
+ Catalog of synthetic sources to estimate the S/N of. **This table
+ will be modified in place**.
+ diffIm : `lsst.afw.image.Exposure`
+ Difference image where the sources were detected.
+ """
+ # Create a schema for the forced measurement task
+ schema = afwTable.SourceTable.makeMinimalSchema()
+ schema.addField("x", "D", "x position in image.", units="pixel")
+ schema.addField("y", "D", "y position in image.", units="pixel")
+ schema.addField("deblend_nChild", "I", "Need for minimal forced phot schema")
+
+ pluginList = [
+ "base_PixelFlags",
+ "base_SdssCentroid",
+ "base_CircularApertureFlux",
+ "base_PsfFlux",
+ "base_LocalBackground"
+ ]
+ forcedMeasConfig = ForcedMeasurementConfig(plugins=pluginList)
+ forcedMeasConfig.slots.centroid = 'base_SdssCentroid'
+ forcedMeasConfig.slots.shape = None
+
+ # Create the forced measurement task
+ forcedMeas = ForcedMeasurementTask(schema, config=forcedMeasConfig)
+
+ # Specify the columns to copy from the input catalog to the output catalog
+ forcedMeas.copyColumns = {"coord_ra": "ra", "coord_dec": "dec"}
+
+ # Create an afw table from the input catalog
+ outputCatalog = afwTable.SourceCatalog(schema)
+ outputCatalog.reserve(len(injectedCat))
+ for row in injectedCat:
+ outputRecord = outputCatalog.addNew()
+ outputRecord.setId(row['injection_id'])
+ outputRecord.setCoord(lsstGeom.SpherePoint(row["ra"], row["dec"], lsstGeom.degrees))
+ outputRecord.set("x", row["x"])
+ outputRecord.set("y", row["y"])
+
+ # Generate the forced measurement catalog
+ forcedSources = forcedMeas.generateMeasCat(diffIm, outputCatalog, diffIm.getWcs())
+ # Attach the PSF shape footprints to the forced measurement catalog
+ forcedMeas.attachPsfShapeFootprints(forcedSources, diffIm)
+
+ # Copy the x and y positions from the forced measurement catalog back
+ # to the input catalog
+ for src, tgt in zip(forcedSources, outputCatalog):
+ src.set('base_SdssCentroid_x', tgt['x'])
+ src.set('base_SdssCentroid_y', tgt['y'])
+
+ # Define the centroid for the forced measurement catalog
+ forcedSources.defineCentroid('base_SdssCentroid')
+ # Run the forced measurement task
+ forcedMeas.run(forcedSources, diffIm, outputCatalog, diffIm.getWcs())
+ # Convert the forced measurement catalog to an astropy table
+ forcedSources_table = forcedSources.asAstropy()
+
+ # Add the forced measurement columns to the input catalog
+ for column in forcedSources_table.columns:
+ if "Flux" in column or "flag" in column:
+ injectedCat["forced_"+column] = forcedSources_table[column]
+
+ # Add the SNR columns to the input catalog
+ for column in injectedCat.colnames:
+ if column.endswith("instFlux"):
+ flux = injectedCat[column]
+ fluxErr = injectedCat[column+"Err"].copy()
+ fluxErr = np.where(
+ (fluxErr <= 0) | (np.isnan(fluxErr)), np.nanmax(fluxErr), fluxErr)
+
+ injectedCat[column+"_SNR"] = flux / fluxErr
+
+ def _processFakes(self, injectedCat, diaSources):
+ """Match fakes to detected diaSources within a difference image bound.
+
+ Parameters
+ ----------
+ injectedCat : `astropy.table.table.Table`
+ Catalog of injected sources to match to detected diaSources.
+ diaSources : `afw.table.SourceCatalog`
+ Catalog of difference image sources detected in ``diffIm``.
+ associatedDiaSources : `pandas.DataFrame`
+ Catalog of associated difference image sources detected in ``diffIm``.
+
+ Returns
+ -------
+ result : `lsst.pipe.base.Struct`
+ Results struct with components.
+
+ - ``matchedDiaSources`` : Fakes matched to input diaSources. Has
+ length of ``fakeCat``. (`pandas.DataFrame`)
+ """
+ # First match the diaSrc to the injected fakes
+ injectedCat = injectedCat.to_pandas()
+ nPossibleFakes = len(injectedCat)
+
+ fakeVects = self._getVectors(
+ np.radians(injectedCat.ra),
+ np.radians(injectedCat.dec))
+ diaSrcVects = self._getVectors(
+ diaSources['coord_ra'],
+ diaSources['coord_dec'])
+
+ diaSrcTree = cKDTree(diaSrcVects)
+ dist, idxs = diaSrcTree.query(
+ fakeVects,
+ distance_upper_bound=np.radians(self.config.matchDistanceArcseconds / 3600))
+ nFakesFound = np.isfinite(dist).sum()
+
+ self.log.info("Found %d out of %d possible in diaSources.", nFakesFound, nPossibleFakes)
+
+ # assign diaSourceId to the matched fakes
+ diaSrcIds = diaSources['id'][np.where(np.isfinite(dist), idxs, 0)]
+ matchedFakes = injectedCat.assign(diaSourceId=np.where(np.isfinite(dist), diaSrcIds, 0))
+ matchedFakes['dist_diaSrc'] = np.where(np.isfinite(dist), 3600*np.rad2deg(dist), -1)
+
+ return Struct(matchDiaSources=matchedFakes)
+
+ def _getVectors(self, ras, decs):
+ """Convert ra dec to unit vectors on the sphere.
+
+ Parameters
+ ----------
+ ras : `numpy.ndarray`, (N,)
+ RA coordinates in radians.
+ decs : `numpy.ndarray`, (N,)
+ Dec coordinates in radians.
+
+ Returns
+ -------
+ vectors : `numpy.ndarray`, (N, 3)
+ Vectors on the unit sphere for the given RA/DEC values.
+ """
+ vectors = np.empty((len(ras), 3))
+
+ vectors[:, 2] = np.sin(decs)
+ vectors[:, 0] = np.cos(decs) * np.cos(ras)
+ vectors[:, 1] = np.cos(decs) * np.sin(ras)
+
+ return vectors
+
+ def _addPixCoords(self, fakeCat, image):
+ """Add pixel coordinates to the catalog of fakes.
+
+ Parameters
+ ----------
+ fakeCat : `astropy.table.table.Table`
+ The catalog of fake sources to be input
+ image : `lsst.afw.image.exposure.exposure.ExposureF`
+ The image into which the fake sources should be added
+ Returns
+ -------
+ fakeCat : `astropy.table.table.Table`
+ """
+
+ wcs = image.getWcs()
+
+ # Get x/y pixel coordinates for injected sources.
+ xs, ys = wcs.skyToPixelArray(
+ fakeCat["ra"],
+ fakeCat["dec"],
+ degrees=True
+ )
+ fakeCat["x"] = xs
+ fakeCat["y"] = ys
+
+ return fakeCat
+
+ def _trimFakeCat(self, fakeCat, image):
+ """Trim the fake cat to the exact size of the input image.
+
+ Parameters
+ ----------
+ fakeCat : `astropy.table.table.Table`
+ The catalog of fake sources that was input
+ image : `lsst.afw.image.exposure.exposure.ExposureF`
+ The image into which the fake sources were added
+ Returns
+ -------
+ fakeCat : `astropy.table.table.Table`
+ The original fakeCat trimmed to the area of the image
+ """
+
+ # fakeCat must be processed with _addPixCoords before trimming
+ fakeCat = self._addPixCoords(fakeCat, image)
+
+ # Prefilter in ra/dec to avoid cases where the wcs incorrectly maps
+ # input fakes which are really off the chip onto it.
+ ras = fakeCat["ra"] * u.deg
+ decs = fakeCat["dec"] * u.deg
+
+ isContainedRaDec = image.containsSkyCoords(ras, decs, padding=0)
+
+ # now use the exact pixel BBox to filter to only fakes that were inserted
+ xs = fakeCat["x"]
+ ys = fakeCat["y"]
+
+ bbox = lsstGeom.Box2D(image.getBBox())
+ isContainedXy = xs - self.config.trimBuffer >= bbox.minX
+ isContainedXy &= xs + self.config.trimBuffer <= bbox.maxX
+ isContainedXy &= ys - self.config.trimBuffer >= bbox.minY
+ isContainedXy &= ys + self.config.trimBuffer <= bbox.maxY
+
+ return fakeCat[isContainedRaDec & isContainedXy]
+
+
+class MatchInjectedToAssocDiaSourceConnections(
+ PipelineTaskConnections,
+ defaultTemplates={"coaddName": "deep",
+ "fakesType": "fakes_"},
+ dimensions=("instrument",
+ "visit",
+ "detector")):
+
+ assocDiaSources = connTypes.Input(
+ doc="An assocDiaSource catalog to match against fakeCat from the"
+ "diaPipe run. Assumed to be SDMified.",
+ name="{fakesType}{coaddName}Diff_assocDiaSrc",
+ storageClass="DataFrame",
+ dimensions=("instrument", "visit", "detector"),
+ )
+ matchDiaSources = connTypes.Input(
+ doc="A catalog of those fakeCat sources that have a match in "
+ "diaSrc. The schema is the union of the schemas for "
+ "``fakeCat`` and ``diaSrc``.",
+ name="{fakesType}{coaddName}Diff_matchDiaSrc",
+ storageClass="DataFrame",
+ dimensions=("instrument", "visit", "detector"),
+ )
+ matchAssocDiaSources = connTypes.Output(
+ doc="A catalog of those fakeCat sources that have a match in "
+ "associatedDiaSources. The schema is the union of the schemas for "
+ "``fakeCat`` and ``associatedDiaSources``.",
+ name="{fakesType}{coaddName}Diff_matchAssocDiaSrc",
+ storageClass="DataFrame",
+ dimensions=("instrument", "visit", "detector"),
+ )
+
+
+class MatchInjectedToAssocDiaSourceConfig(
+ PipelineTaskConfig,
+ pipelineConnections=MatchInjectedToAssocDiaSourceConnections):
+ """Config for MatchFakesTask.
+ """
+
+
+class MatchInjectedToAssocDiaSourceTask(PipelineTask):
+
+ _DefaultName = "matchInjectedToAssocDiaSource"
+ ConfigClass = MatchInjectedToAssocDiaSourceConfig
+
+ def run(self, assocDiaSources, matchDiaSources):
+ """Tag matched injected sources to associated diaSources.
+
+ Parameters
+ ----------
+ matchDiaSources : `pandas.DataFrame`
+ Catalog of matched diaSrc to injected sources
+ assocDiaSources : `pandas.DataFrame`
+ Catalog of associated difference image sources detected in ``diffIm``.
+ Returns
+ -------
+ result : `lsst.pipe.base.Struct`
+ Results struct with components.
+
+ - ``matchAssocDiaSources`` : Fakes matched to associated diaSources. Has
+ length of ``matchDiaSources``. (`pandas.DataFrame`)
+ """
+ # Match the fakes to the associated sources. For this we don't use the coordinates
+ # but instead check for the diaSources. Since they were present in the table already
+ nPossibleFakes = len(matchDiaSources)
+ matchDiaSources['isAssocDiaSource'] = matchDiaSources.diaSourceId.isin(assocDiaSources.diaSourceId)
+ assocNFakesFound = matchDiaSources.isAssocDiaSource.sum()
+ self.log.info("Found %d out of %d possible in assocDiaSources."%(assocNFakesFound, nPossibleFakes))
+
+ return Struct(
+ matchAssocDiaSources=matchDiaSources.merge(
+ assocDiaSources.reset_index(drop=True),
+ on="diaSourceId",
+ how="left",
+ suffixes=('_ssi', '_diaSrc')
+ )
+ )
diff --git a/python/lsst/pipe/tasks/matchFakes.py b/python/lsst/pipe/tasks/matchFakes.py
deleted file mode 100644
index c9cd68d12..000000000
--- a/python/lsst/pipe/tasks/matchFakes.py
+++ /dev/null
@@ -1,445 +0,0 @@
-# This file is part of pipe_tasks.
-#
-# Developed for the LSST Data Management System.
-# This product includes software developed by the LSST Project
-# (https://www.lsst.org).
-# See the COPYRIGHT file at the top-level directory of this distribution
-# for details of code ownership.
-#
-# This program is free software: you can redistribute it and/or modify
-# it under the terms of the GNU General Public License as published by
-# the Free Software Foundation, either version 3 of the License, or
-# (at your option) any later version.
-#
-# This program is distributed in the hope that it will be useful,
-# but WITHOUT ANY WARRANTY; without even the implied warranty of
-# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
-# GNU General Public License for more details.
-#
-# You should have received a copy of the GNU General Public License
-# along with this program. If not, see .
-
-__all__ = ["MatchFakesTask",
- "MatchFakesConfig",
- "MatchVariableFakesConfig",
- "MatchVariableFakesTask"]
-
-import astropy.units as u
-import numpy as np
-import pandas as pd
-from scipy.spatial import cKDTree
-
-from lsst.geom import Box2D
-import lsst.pex.config as pexConfig
-from lsst.pipe.base import PipelineTask, PipelineTaskConnections, Struct
-import lsst.pipe.base.connectionTypes as connTypes
-from lsst.skymap import BaseSkyMap
-
-from lsst.pipe.tasks.insertFakes import InsertFakesConfig
-
-from deprecated.sphinx import deprecated
-
-
-class MatchFakesConnections(PipelineTaskConnections,
- defaultTemplates={"coaddName": "deep",
- "fakesType": "fakes_"},
- dimensions=("instrument",
- "visit",
- "detector")):
- skyMap = connTypes.Input(
- doc="Input definition of geometry/bbox and projection/wcs for "
- "template exposures. Needed to test which tract to generate ",
- name=BaseSkyMap.SKYMAP_DATASET_TYPE_NAME,
- dimensions=("skymap",),
- storageClass="SkyMap",
- )
- fakeCats = connTypes.Input(
- doc="Catalog of fake sources inserted into an image.",
- name="{fakesType}fakeSourceCat",
- storageClass="DataFrame",
- dimensions=("tract", "skymap"),
- deferLoad=True,
- multiple=True
- )
- diffIm = connTypes.Input(
- doc="Difference image on which the DiaSources were detected.",
- name="{fakesType}{coaddName}Diff_differenceExp",
- storageClass="ExposureF",
- dimensions=("instrument", "visit", "detector"),
- )
- associatedDiaSources = connTypes.Input(
- doc="A DiaSource catalog to match against fakeCat. Assumed "
- "to be SDMified.",
- name="{fakesType}{coaddName}Diff_assocDiaSrc",
- storageClass="DataFrame",
- dimensions=("instrument", "visit", "detector"),
- )
- matchedDiaSources = connTypes.Output(
- doc="A catalog of those fakeCat sources that have a match in "
- "associatedDiaSources. The schema is the union of the schemas for "
- "``fakeCat`` and ``associatedDiaSources``.",
- name="{fakesType}{coaddName}Diff_matchDiaSrc",
- storageClass="DataFrame",
- dimensions=("instrument", "visit", "detector"),
- )
-
-
-@deprecated(
- reason="This task will be removed in v28.0 as it is replaced by `source_injection` tasks.",
- version="v28.0",
- category=FutureWarning,
-)
-class MatchFakesConfig(
- InsertFakesConfig,
- pipelineConnections=MatchFakesConnections):
- """Config for MatchFakesTask.
- """
- matchDistanceArcseconds = pexConfig.RangeField(
- doc="Distance in arcseconds to match within.",
- dtype=float,
- default=0.5,
- min=0,
- max=10,
- )
-
- doMatchVisit = pexConfig.Field(
- dtype=bool,
- default=False,
- doc="Match visit to trim the fakeCat"
- )
-
- trimBuffer = pexConfig.Field(
- doc="Size of the pixel buffer surrounding the image. Only those fake sources with a centroid"
- "falling within the image+buffer region will be considered matches.",
- dtype=int,
- default=100,
- )
-
-
-@deprecated(
- reason="This task will be removed in v28.0 as it is replaced by `source_injection` tasks.",
- version="v28.0",
- category=FutureWarning,
-)
-class MatchFakesTask(PipelineTask):
- """Match a pre-existing catalog of fakes to a catalog of detections on
- a difference image.
-
- This task is generally for injected sources that cannot be easily
- identified by their footprints such as in the case of detector sources
- post image differencing.
- """
-
- _DefaultName = "matchFakes"
- ConfigClass = MatchFakesConfig
-
- def run(self, fakeCats, skyMap, diffIm, associatedDiaSources):
- """Compose fakes into a single catalog and match fakes to detected
- diaSources within a difference image bound.
-
- Parameters
- ----------
- fakeCats : `pandas.DataFrame`
- List of catalog of fakes to match to detected diaSources.
- skyMap : `lsst.skymap.SkyMap`
- SkyMap defining the tracts and patches the fakes are stored over.
- diffIm : `lsst.afw.image.Exposure`
- Difference image where ``associatedDiaSources`` were detected.
- associatedDiaSources : `pandas.DataFrame`
- Catalog of difference image sources detected in ``diffIm``.
-
- Returns
- -------
- result : `lsst.pipe.base.Struct`
- Results struct with components.
-
- - ``matchedDiaSources`` : Fakes matched to input diaSources. Has
- length of ``fakeCat``. (`pandas.DataFrame`)
- """
- fakeCat = self.composeFakeCat(fakeCats, skyMap)
-
- if self.config.doMatchVisit:
- fakeCat = self.getVisitMatchedFakeCat(fakeCat, diffIm)
-
- return self._processFakes(fakeCat, diffIm, associatedDiaSources)
-
- def _processFakes(self, fakeCat, diffIm, associatedDiaSources):
- """Match fakes to detected diaSources within a difference image bound.
-
- Parameters
- ----------
- fakeCat : `pandas.DataFrame`
- Catalog of fakes to match to detected diaSources.
- diffIm : `lsst.afw.image.Exposure`
- Difference image where ``associatedDiaSources`` were detected.
- associatedDiaSources : `pandas.DataFrame`
- Catalog of difference image sources detected in ``diffIm``.
-
- Returns
- -------
- result : `lsst.pipe.base.Struct`
- Results struct with components.
-
- - ``matchedDiaSources`` : Fakes matched to input diaSources. Has
- length of ``fakeCat``. (`pandas.DataFrame`)
- """
- trimmedFakes = self._trimFakeCat(fakeCat, diffIm)
- nPossibleFakes = len(trimmedFakes)
-
- fakeVects = self._getVectors(trimmedFakes[self.config.ra_col],
- trimmedFakes[self.config.dec_col])
- diaSrcVects = self._getVectors(
- np.radians(associatedDiaSources.loc[:, "ra"]),
- np.radians(associatedDiaSources.loc[:, "dec"]))
-
- diaSrcTree = cKDTree(diaSrcVects)
- dist, idxs = diaSrcTree.query(
- fakeVects,
- distance_upper_bound=np.radians(self.config.matchDistanceArcseconds / 3600))
- nFakesFound = np.isfinite(dist).sum()
-
- self.log.info("Found %d out of %d possible.", nFakesFound, nPossibleFakes)
- diaSrcIds = associatedDiaSources.iloc[np.where(np.isfinite(dist), idxs, 0)]["diaSourceId"].to_numpy()
- matchedFakes = trimmedFakes.assign(diaSourceId=np.where(np.isfinite(dist), diaSrcIds, 0))
-
- return Struct(
- matchedDiaSources=matchedFakes.merge(
- associatedDiaSources.reset_index(drop=True), on="diaSourceId", how="left")
- )
-
- def composeFakeCat(self, fakeCats, skyMap):
- """Concatenate the fakeCats from tracts that may cover the exposure.
-
- Parameters
- ----------
- fakeCats : `list` of `lsst.daf.butler.DeferredDatasetHandle`
- Set of fake cats to concatenate.
- skyMap : `lsst.skymap.SkyMap`
- SkyMap defining the geometry of the tracts and patches.
-
- Returns
- -------
- combinedFakeCat : `pandas.DataFrame`
- All fakes that cover the inner polygon of the tracts in this
- quantum.
- """
- if len(fakeCats) == 1:
- return fakeCats[0].get()
- outputCat = []
- for fakeCatRef in fakeCats:
- cat = fakeCatRef.get()
- tractId = fakeCatRef.dataId["tract"]
- # Make sure all data is within the inner part of the tract.
- outputCat.append(cat[
- skyMap.findTractIdArray(cat[self.config.ra_col],
- cat[self.config.dec_col],
- degrees=False)
- == tractId])
-
- return pd.concat(outputCat)
-
- def getVisitMatchedFakeCat(self, fakeCat, exposure):
- """Trim the fakeCat to select particular visit
-
- Parameters
- ----------
- fakeCat : `pandas.core.frame.DataFrame`
- The catalog of fake sources to add to the exposure
- exposure : `lsst.afw.image.exposure.exposure.ExposureF`
- The exposure to add the fake sources to
-
- Returns
- -------
- movingFakeCat : `pandas.DataFrame`
- All fakes that belong to the visit
- """
- selected = exposure.getInfo().getVisitInfo().getId() == fakeCat["visit"]
-
- return fakeCat[selected]
-
- def _addPixCoords(self, fakeCat, image):
-
- """Add pixel coordinates to the catalog of fakes.
-
- Parameters
- ----------
- fakeCat : `pandas.core.frame.DataFrame`
- The catalog of fake sources to be input
- image : `lsst.afw.image.exposure.exposure.ExposureF`
- The image into which the fake sources should be added
- Returns
- -------
- fakeCat : `pandas.core.frame.DataFrame`
- """
- wcs = image.getWcs()
- ras = fakeCat[self.config.ra_col].values
- decs = fakeCat[self.config.dec_col].values
- xs, ys = wcs.skyToPixelArray(ras, decs)
- fakeCat["x"] = xs
- fakeCat["y"] = ys
-
- return fakeCat
-
- def _trimFakeCat(self, fakeCat, image):
- """Trim the fake cat to the exact size of the input image.
-
- Parameters
- ----------
- fakeCat : `pandas.core.frame.DataFrame`
- The catalog of fake sources that was input
- image : `lsst.afw.image.exposure.exposure.ExposureF`
- The image into which the fake sources were added
- Returns
- -------
- fakeCat : `pandas.core.frame.DataFrame`
- The original fakeCat trimmed to the area of the image
- """
-
- # fakeCat must be processed with _addPixCoords before trimming
- if ('x' not in fakeCat.columns) or ('y' not in fakeCat.columns):
- fakeCat = self._addPixCoords(fakeCat, image)
-
- # Prefilter in ra/dec to avoid cases where the wcs incorrectly maps
- # input fakes which are really off the chip onto it.
- ras = fakeCat[self.config.ra_col].values * u.rad
- decs = fakeCat[self.config.dec_col].values * u.rad
-
- isContainedRaDec = image.containsSkyCoords(ras, decs, padding=0)
-
- # now use the exact pixel BBox to filter to only fakes that were inserted
- xs = fakeCat["x"].values
- ys = fakeCat["y"].values
-
- bbox = Box2D(image.getBBox())
- isContainedXy = xs >= bbox.minX
- isContainedXy &= xs <= bbox.maxX
- isContainedXy &= ys >= bbox.minY
- isContainedXy &= ys <= bbox.maxY
-
- return fakeCat[isContainedRaDec & isContainedXy]
-
- def _getVectors(self, ras, decs):
- """Convert ra dec to unit vectors on the sphere.
-
- Parameters
- ----------
- ras : `numpy.ndarray`, (N,)
- RA coordinates in radians.
- decs : `numpy.ndarray`, (N,)
- Dec coordinates in radians.
-
- Returns
- -------
- vectors : `numpy.ndarray`, (N, 3)
- Vectors on the unit sphere for the given RA/DEC values.
- """
- vectors = np.empty((len(ras), 3))
-
- vectors[:, 2] = np.sin(decs)
- vectors[:, 0] = np.cos(decs) * np.cos(ras)
- vectors[:, 1] = np.cos(decs) * np.sin(ras)
-
- return vectors
-
-
-class MatchVariableFakesConnections(MatchFakesConnections):
- ccdVisitFakeMagnitudes = connTypes.Input(
- doc="Catalog of fakes with magnitudes scattered for this ccdVisit.",
- name="{fakesType}ccdVisitFakeMagnitudes",
- storageClass="DataFrame",
- dimensions=("instrument", "visit", "detector"),
- )
-
-
-@deprecated(
- reason="This task will be removed in v28.0 as it is replaced by `source_injection` tasks.",
- version="v28.0",
- category=FutureWarning,
-)
-class MatchVariableFakesConfig(MatchFakesConfig,
- pipelineConnections=MatchVariableFakesConnections):
- """Config for MatchFakesTask.
- """
- pass
-
-
-@deprecated(
- reason="This task will be removed in v28.0 as it is replaced by `source_injection` tasks.",
- version="v28.0",
- category=FutureWarning,
-)
-class MatchVariableFakesTask(MatchFakesTask):
- """Match injected fakes to their detected sources in the catalog and
- compute their expected brightness in a difference image assuming perfect
- subtraction.
-
- This task is generally for injected sources that cannot be easily
- identified by their footprints such as in the case of detector sources
- post image differencing.
- """
- _DefaultName = "matchVariableFakes"
- ConfigClass = MatchVariableFakesConfig
-
- def runQuantum(self, butlerQC, inputRefs, outputRefs):
- inputs = butlerQC.get(inputRefs)
- inputs["band"] = butlerQC.quantum.dataId["band"]
-
- outputs = self.run(**inputs)
- butlerQC.put(outputs, outputRefs)
-
- def run(self, fakeCats, ccdVisitFakeMagnitudes, skyMap, diffIm, associatedDiaSources, band):
- """Match fakes to detected diaSources within a difference image bound.
-
- Parameters
- ----------
- fakeCat : `pandas.DataFrame`
- Catalog of fakes to match to detected diaSources.
- diffIm : `lsst.afw.image.Exposure`
- Difference image where ``associatedDiaSources`` were detected in.
- associatedDiaSources : `pandas.DataFrame`
- Catalog of difference image sources detected in ``diffIm``.
-
- Returns
- -------
- result : `lsst.pipe.base.Struct`
- Results struct with components.
-
- - ``matchedDiaSources`` : Fakes matched to input diaSources. Has
- length of ``fakeCat``. (`pandas.DataFrame`)
- """
- fakeCat = self.composeFakeCat(fakeCats, skyMap)
- self.computeExpectedDiffMag(fakeCat, ccdVisitFakeMagnitudes, band)
- return self._processFakes(fakeCat, diffIm, associatedDiaSources)
-
- def computeExpectedDiffMag(self, fakeCat, ccdVisitFakeMagnitudes, band):
- """Compute the magnitude expected in the difference image for this
- detector/visit. Modify fakeCat in place.
-
- Negative magnitudes indicate that the source should be detected as
- a negative source.
-
- Parameters
- ----------
- fakeCat : `pandas.DataFrame`
- Catalog of fake sources.
- ccdVisitFakeMagnitudes : `pandas.DataFrame`
- Magnitudes for variable sources in this specific ccdVisit.
- band : `str`
- Band that this ccdVisit was observed in.
- """
- magName = self.config.mag_col % band
- magnitudes = fakeCat[magName].to_numpy()
- visitMags = ccdVisitFakeMagnitudes["variableMag"].to_numpy()
- diffFlux = (visitMags * u.ABmag).to_value(u.nJy) - (magnitudes * u.ABmag).to_value(u.nJy)
- diffMag = np.where(diffFlux > 0,
- (diffFlux * u.nJy).to_value(u.ABmag),
- -(-diffFlux * u.nJy).to_value(u.ABmag))
-
- noVisit = ~fakeCat["isVisitSource"]
- noTemplate = ~fakeCat["isTemplateSource"]
- both = np.logical_and(fakeCat["isVisitSource"],
- fakeCat["isTemplateSource"])
-
- fakeCat.loc[noVisit, magName] = -magnitudes[noVisit]
- fakeCat.loc[noTemplate, magName] = visitMags[noTemplate]
- fakeCat.loc[both, magName] = diffMag[both]
diff --git a/tests/test_matchFakes.py b/tests/test_matchFakes.py
deleted file mode 100644
index 4e878d1ba..000000000
--- a/tests/test_matchFakes.py
+++ /dev/null
@@ -1,141 +0,0 @@
-#
-# This file is part of pipe_tasks.
-#
-# Developed for the LSST Data Management System.
-# This product includes software developed by the LSST Project
-# (http://www.lsst.org).
-# See the COPYRIGHT file at the top-level directory of this distribution
-# for details of code ownership.
-#
-# This program is free software: you can redistribute it and/or modify
-# it under the terms of the GNU General Public License as published by
-# the Free Software Foundation, either version 3 of the License, or
-# (at your option) any later version.
-#
-# This program is distributed in the hope that it will be useful,
-# but WITHOUT ANY WARRANTY; without even the implied warranty of
-# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
-# GNU General Public License for more details.
-#
-# You should have received a copy of the GNU General Public License
-# along with this program. If not, see .
-#
-
-import numpy as np
-import pandas as pd
-import unittest
-import uuid
-
-import lsst.sphgeom as sphgeom
-import lsst.geom as geom
-import lsst.meas.base.tests as measTests
-import lsst.skymap as skyMap
-import lsst.utils.tests
-
-from lsst.pipe.tasks.matchFakes import MatchFakesTask, MatchFakesConfig
-
-
-class TestMatchFakes(lsst.utils.tests.TestCase):
-
- def setUp(self):
- """Create fake data to use in the tests.
- """
- self.bbox = geom.Box2I(geom.Point2I(0, 0),
- geom.Extent2I(1024, 1153))
- dataset = measTests.TestDataset(self.bbox)
- self.exposure = dataset.exposure
-
- simpleMapConfig = skyMap.discreteSkyMap.DiscreteSkyMapConfig()
- simpleMapConfig.raList = [dataset.exposure.getWcs().getSkyOrigin().getRa().asDegrees()]
- simpleMapConfig.decList = [dataset.exposure.getWcs().getSkyOrigin().getDec().asDegrees()]
- simpleMapConfig.radiusList = [0.1]
-
- self.simpleMap = skyMap.DiscreteSkyMap(simpleMapConfig)
- self.tractId = 0
- bCircle = self.simpleMap.generateTract(self.tractId).getOuterSkyPolygon().getBoundingCircle()
- bCenter = sphgeom.LonLat(bCircle.getCenter())
- bRadius = bCircle.getOpeningAngle().asRadians()
- targetSources = 10000
- self.sourceDensity = (targetSources
- / (bCircle.getArea() * (180 / np.pi) ** 2))
- self.rng = np.random.default_rng(1234)
-
- self.fakeCat = pd.DataFrame({
- "fakeId": [uuid.uuid4().int & (1 << 64) - 1 for n in range(targetSources)],
- # Quick-and-dirty values for testing
- "ra": bCenter.getLon().asRadians() + bRadius * (2.0 * self.rng.random(targetSources) - 1.0),
- "dec": bCenter.getLat().asRadians() + bRadius * (2.0 * self.rng.random(targetSources) - 1.0),
- "isVisitSource": np.concatenate([np.ones(targetSources//2, dtype="bool"),
- np.zeros(targetSources - targetSources//2, dtype="bool")]),
- "isTemplateSource": np.concatenate([np.zeros(targetSources//2, dtype="bool"),
- np.ones(targetSources - targetSources//2, dtype="bool")]),
- **{band: self.rng.uniform(20, 30, size=targetSources)
- for band in {"u", "g", "r", "i", "z", "y"}},
- "DiskHalfLightRadius": np.ones(targetSources, dtype="float"),
- "BulgeHalfLightRadius": np.ones(targetSources, dtype="float"),
- "disk_n": np.ones(targetSources, dtype="float"),
- "bulge_n": np.ones(targetSources, dtype="float"),
- "a_d": np.ones(targetSources, dtype="float"),
- "a_b": np.ones(targetSources, dtype="float"),
- "b_d": np.ones(targetSources, dtype="float"),
- "b_b": np.ones(targetSources, dtype="float"),
- "pa_disk": np.ones(targetSources, dtype="float"),
- "pa_bulge": np.ones(targetSources, dtype="float"),
- "sourceType": targetSources * ["star"],
- })
-
- self.inExp = np.zeros(len(self.fakeCat), dtype=bool)
- bbox = geom.Box2D(self.exposure.getBBox())
- for idx, row in self.fakeCat.iterrows():
- coord = geom.SpherePoint(row["ra"],
- row["dec"],
- geom.radians)
- cent = self.exposure.getWcs().skyToPixel(coord)
- self.inExp[idx] = bbox.contains(cent)
-
- tmpCat = self.fakeCat[self.inExp].iloc[:int(self.inExp.sum() / 2)]
- extraColumnData = self.rng.integers(0, 100, size=len(tmpCat))
- self.sourceCat = pd.DataFrame(
- data={"ra": np.degrees(tmpCat["ra"]),
- "dec": np.degrees(tmpCat["dec"]),
- "diaObjectId": np.arange(1, len(tmpCat) + 1, dtype=int),
- "band": "g",
- "diaSourceId": np.arange(1, len(tmpCat) + 1, dtype=int),
- "extraColumn": extraColumnData})
- self.sourceCat.set_index(["diaObjectId", "band", "extraColumn"],
- drop=False,
- inplace=True)
-
- def testProcessFakes(self):
- """Test the run method.
- """
- matchFakesConfig = MatchFakesConfig()
- matchFakesConfig.matchDistanceArcseconds = 0.1
- matchFakes = MatchFakesTask(config=matchFakesConfig)
- result = matchFakes._processFakes(self.fakeCat,
- self.exposure,
- self.sourceCat)
- self.assertEqual(self.inExp.sum(), len(result.matchedDiaSources))
- self.assertEqual(
- len(self.sourceCat),
- np.sum(np.isfinite(result.matchedDiaSources["extraColumn"])))
-
- def testTrimCat(self):
- """Test that the correct number of sources are in the ccd area.
- """
- matchTask = MatchFakesTask()
- result = matchTask._trimFakeCat(self.fakeCat, self.exposure)
- self.assertEqual(len(result), self.inExp.sum())
-
-
-class MemoryTester(lsst.utils.tests.MemoryTestCase):
- pass
-
-
-def setup_module(module):
- lsst.utils.tests.init()
-
-
-if __name__ == "__main__":
- lsst.utils.tests.init()
- unittest.main()
diff --git a/tests/test_matchSourceInjected.py b/tests/test_matchSourceInjected.py
new file mode 100644
index 000000000..9d2ee58f6
--- /dev/null
+++ b/tests/test_matchSourceInjected.py
@@ -0,0 +1,188 @@
+#
+# This file is part of ap_pipe.
+#
+# Developed for the LSST Data Management System.
+# This product includes software developed by the LSST Project
+# (http://www.lsst.org).
+# See the COPYRIGHT file at the top-level directory of this distribution
+# for details of code ownership.
+#
+# This program is free software: you can redistribute it and/or modify
+# it under the terms of the GNU General Public License as published by
+# the Free Software Foundation, either version 3 of the License, or
+# (at your option) any later version.
+#
+# This program is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+# GNU General Public License for more details.
+#
+# You should have received a copy of the GNU General Public License
+# along with this program. If not, see .
+#
+
+import numpy as np
+import pandas as pd
+import unittest
+
+from astropy.table import Table
+from lsst.afw.table import SourceCatalog, SourceTable
+from lsst.afw.image import ExposureF
+from lsst.afw.geom import makeSkyWcs
+from lsst.geom import SpherePoint, degrees, Point2D, Extent2I
+import lsst.utils.tests
+
+from lsst.pipe.tasks.matchDiffimSourceInjected import (
+ MatchInjectedToDiaSourceTask,
+ MatchInjectedToAssocDiaSourceTask,
+ MatchInjectedToDiaSourceConfig,
+ MatchInjectedToAssocDiaSourceConfig,
+)
+
+
+class BaseTestMatchInjected(lsst.utils.tests.TestCase):
+ def setUp(self):
+ rng = np.random.RandomState(6)
+ # 0.35 arcsec = 0.5/np.sqrt(2) arcsec (half of the diagonal of a 0.5 arcsec square)
+ offsetFactor = 1./3600 * 0.35
+
+ # Create a mock injected catalog
+ self.injectedCat = Table(
+ {
+ "injection_id": [1, 2, 3, 5, 6, 7, 12, 21, 31, 49],
+ "injection_flag": np.repeat(0, 10),
+ # random positions with 10 arcmin size
+ # ra centered at 0, dec centered at -30
+ "ra": (1/6.) * rng.random(size=10),
+ "dec": -30 + (1/6.) * rng.random(size=10),
+ "mag": 20.25 - 0.5 * rng.random(size=10), # random magnitudes
+ "source_type": np.repeat("DeltaFunction", 10)
+ }
+ )
+
+ # Create a mock diaSources catalog
+ schema = SourceTable.makeMinimalSchema()
+ self.diaSources = SourceCatalog(schema)
+ for i in range(5):
+ record = self.diaSources.addNew()
+ record.setId(100 + i)
+ record.setCoord(
+ # Random posisions of diaSources
+ SpherePoint((1/6.) * rng.random(), -30 + (1/6.) * rng.random(), degrees)
+ )
+
+ for i, (ra, dec) in enumerate(self.injectedCat[['ra', 'dec']][:8]):
+ record = self.diaSources.addNew()
+ sign = rng.choice([-1, 1], size=2)
+ record.setId(i + 200)
+ record.setCoord(
+ SpherePoint(
+ ra+sign[0]*rng.random()*offsetFactor,
+ dec+sign[1]*rng.random()*offsetFactor,
+ degrees
+ )
+ )
+
+ # Create a mock difference image
+ self.diffIm = ExposureF(Extent2I(4096, 4096))
+ crpix = Point2D(0, 0)
+ crval = SpherePoint(0, -30, degrees)
+ cdMatrix = np.array([[5.19513851e-05, 2.81124812e-07],
+ [3.25186974e-07, 5.19112119e-05]])
+ wcs = makeSkyWcs(crpix, crval, cdMatrix)
+ self.diffIm.setWcs(wcs)
+
+ # add a fake source outside of the image
+ self.injectedCatForTrimming = self.injectedCat.copy()
+ self.injectedCatForTrimming.add_row(
+ {
+ 'injection_id': 50,
+ 'injection_flag': 0,
+ 'ra': 360. - 0.1/6.,
+ 'dec': -30 - 0.1/6.,
+ 'mag': 20.0,
+ 'source_type': 'DeltaFunction'
+ }
+ )
+
+ # only 4 injected sources are associated
+ self.assocDiaSources = pd.DataFrame(
+ {
+ "diaSourceId": [101, 102, 103, 201, 202, 205, 207],
+ "band": np.repeat("r", 7),
+ "visit": np.repeat(410001, 7),
+ "detector": np.repeat(0, 7),
+ "diaObjectId": np.arange(7),
+ }
+ )
+
+
+class TestMatchInjectedToDiaSourceTask(BaseTestMatchInjected):
+
+ def test_run(self):
+ config = MatchInjectedToDiaSourceConfig()
+ config.matchDistanceArcseconds = 0.5
+ config.doMatchVisit = False
+ config.doForcedMeasurement = False
+
+ task = MatchInjectedToDiaSourceTask(config=config)
+
+ result = task.run(self.injectedCat, self.diffIm, self.diaSources)
+ self.assertEqual(len(result.matchDiaSources), len(self.injectedCat))
+ self.assertEqual(np.sum(result.matchDiaSources['diaSourceId'] > 0), 8)
+ self.assertEqual(np.sum(result.matchDiaSources['dist_diaSrc'] > 0), 8)
+ self.assertEqual(
+ np.sum(np.abs(result.matchDiaSources['dist_diaSrc']) < config.matchDistanceArcseconds), 8
+ )
+
+ def test_run_trimming(self):
+ config = MatchInjectedToDiaSourceConfig()
+ config.matchDistanceArcseconds = 0.5
+ config.doMatchVisit = True
+ config.doForcedMeasurement = False
+
+ task = MatchInjectedToDiaSourceTask(config=config)
+ result = task.run(self.injectedCatForTrimming, self.diffIm, self.diaSources)
+
+ self.assertEqual(len(result.matchDiaSources), len(self.injectedCatForTrimming) - 1)
+ self.assertEqual(np.sum(result.matchDiaSources['diaSourceId'] > 0), 8)
+ self.assertEqual(np.sum(result.matchDiaSources['dist_diaSrc'] > 0), 8)
+ self.assertEqual(
+ np.sum(np.abs(result.matchDiaSources['dist_diaSrc']) < config.matchDistanceArcseconds), 8
+ )
+
+ def test_getVectors(self):
+ config = MatchInjectedToDiaSourceConfig()
+ config.matchDistanceArcseconds = 0.5
+ config.doMatchVisit = False
+ config.doForcedMeasurement = False
+
+ task = MatchInjectedToDiaSourceTask(config=config)
+
+ ras = np.radians(self.injectedCat['ra'])
+ decs = np.radians(self.injectedCat['dec'])
+ vectors = task._getVectors(ras, decs)
+ self.assertEqual(vectors.shape, (10, 3))
+
+
+class TestMatchInjectedToAssocDiaSourceTask(BaseTestMatchInjected):
+
+ def test_run(self):
+ config = MatchInjectedToDiaSourceConfig()
+ config.matchDistanceArcseconds = 0.5
+ config.doMatchVisit = False
+ config.doForcedMeasurement = False
+
+ task = MatchInjectedToDiaSourceTask(config=config)
+
+ result = task.run(self.injectedCat, self.diffIm, self.diaSources)
+
+ configAssoc = MatchInjectedToAssocDiaSourceConfig()
+ taskAssoc = MatchInjectedToAssocDiaSourceTask(config=configAssoc)
+ resultAssoc = taskAssoc.run(self.assocDiaSources, result.matchDiaSources)
+ self.assertEqual(len(resultAssoc.matchAssocDiaSources), len(self.injectedCat))
+ self.assertEqual(np.sum(resultAssoc.matchAssocDiaSources['isAssocDiaSource']), 4)
+
+
+if __name__ == "__main__":
+ unittest.main()