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
propagateVisitFlags.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__ = ["PropagateVisitFlagsConfig", "PropagateVisitFlagsTask"]
23import numpy
24from lsst.pex.config import Config, Field, DictField
25from lsst.pipe.base import Task
26import lsst.geom as geom
27import lsst.afw.table as afwTable
28import lsst.pex.exceptions as pexExceptions
29from deprecated.sphinx import deprecated
30
31
33 """Configuration for propagating flags to coadd."""
34
35 flags = DictField(keytype=str, itemtype=float,
36 default={"calib_psf_candidate": 0.2, "calib_psf_used": 0.2, "calib_psf_reserved": 0.2,
37 "calib_astrometry_used": 0.2, "calib_photometry_used": 0.2,
38 "calib_photometry_reserved": 0.2, },
39 doc=("Source catalog flags to propagate, with the threshold of relative occurrence "
40 "(valid range: [0-1], default is 0.2). Coadd object will have flag set if the "
41 "fraction of input visits in which it is flagged is greater than the threshold."))
42 matchRadius = Field(dtype=float, default=0.2, doc="Source matching radius (arcsec)")
43 ccdName = Field(dtype=str, default='ccd', doc="Name of ccd to give to butler")
44
45
46@deprecated(reason="This task has been replaced with PropagateSourceFlagsTask",
47 version="v24.0", category=FutureWarning)
49 """Task to propagate flags from single-frame measurements to coadd measurements.
50
51 Parameters
52 ----------
53 schema : `lsst.afw.table.Schema`
54 The input schema for the reference source catalog, used to initialize
55 the output schema.
56 **kwargs
57 Additional keyword arguments.
58
59 Notes
60 -----
61 We want to be able to set a flag for sources on the coadds using flags
62 that were determined from the individual visits. A common example is sources
63 that were used for PSF determination, since we do not do any PSF determination
64 on the coadd but use the individual visits. This requires matching the coadd
65 source catalog to each of the catalogs from the inputs (see
66 PropagateVisitFlagsConfig.matchRadius), and thresholding on the number of
67 times a source is flagged on the input catalog.
68
69 An important consideration in this is that the flagging of sources in the
70 individual visits can be somewhat stochastic, e.g., the same stars may not
71 always be used for PSF determination because the field of view moves slightly
72 between visits, or the seeing changed. We there threshold on the relative
73 occurrence of the flag in the visits (see PropagateVisitFlagsConfig.flags).
74 Flagging a source that is always flagged in inputs corresponds to a threshold
75 of 1, while flagging a source that is flagged in any of the input corresponds
76 to a threshold of 0. But neither of these extrema are really useful in
77 practise.
78
79 Setting the threshold too high means that sources that are not consistently
80 flagged (e.g., due to chip gaps) will not have the flag propagated. Setting
81 that threshold too low means that random sources which are falsely flagged in
82 the inputs will start to dominate. If in doubt, we suggest making this
83 threshold relatively low, but not zero (e.g., 0.1 to 0.2 or so). The more
84 confidence in the quality of the flagging, the lower the threshold can be.
85 The relative occurrence accounts for the edge of the field-of-view of the
86 camera, but does not include chip gaps, bad or saturated pixels, etc.
87
88 Initialization
89
90 Beyond the usual Task initialization, PropagateVisitFlagsTask also requires
91 a schema for the catalog that is being constructed.
92
93 The 'run' method (described below) is the entry-point for operations. The
94 'getCcdInputs' staticmethod is provided as a convenience for retrieving the
95 'ccdInputs' (CCD inputs table) from an Exposure.
96
97 .. code-block :: none
98
99 # Requires:
100 # * butler: data butler, for retrieving the CCD catalogs
101 # * coaddCatalog: catalog of source measurements on the coadd (lsst.afw.table.SourceCatalog)
102 # * coaddExposure: coadd (lsst.afw.image.Exposure)
103 from lsst.pipe.tasks.propagateVisitFlags import PropagateVisitFlagsTask, PropagateVisitFlagsConfig
104 config = PropagateVisitFlagsConfig()
105 config.flags["calib_psf_used"] = 0.3 # Relative threshold for this flag
106 config.matchRadius = 0.5 # Matching radius in arcsec
107 task = PropagateVisitFlagsTask(coaddCatalog.schema, config=config)
108 ccdInputs = task.getCcdInputs(coaddExposure)
109 task.run(butler, coaddCatalog, ccdInputs, coaddExposure.getWcs())
110 """
111
112 ConfigClass = PropagateVisitFlagsConfig
113
114 def __init__(self, schema, **kwargs):
115 Task.__init__(self, **kwargs)
116 self.schema = schema
117 self._keys = dict((f, self.schema.addField(f, type="Flag", doc="Propagated from visits")) for
118 f in self.config.flags)
119
120 @staticmethod
121 def getCcdInputs(coaddExposure):
122 """Convenience method to retrieve the CCD inputs table from a coadd exposure.
123
124 Parameters
125 ----------
126 coaddExposure : `lsst.afw.image.Exposure`
127 The exposure we need to retrieve the CCD inputs table from.
128
129 Returns
130 -------
131 ccdInputs : ``
132 CCD inputs table from a coadd exposure.
133 """
134 return coaddExposure.getInfo().getCoaddInputs().ccds
135
136 def run(self, butler, coaddSources, ccdInputs, coaddWcs, visitCatalogs=None, wcsUpdates=None):
137 """Propagate flags from individual visit measurements to coadd.
138
139 This requires matching the coadd source catalog to each of the catalogs
140 from the inputs, and thresholding on the number of times a source is
141 flagged on the input catalog. The threshold is made on the relative
142 occurrence of the flag in each source. Flagging a source that is always
143 flagged in inputs corresponds to a threshold of 1, while flagging a
144 source that is flagged in any of the input corresponds to a threshold of
145 0. But neither of these extrema are really useful in practise.
146
147 Setting the threshold too high means that sources that are not consistently
148 flagged (e.g., due to chip gaps) will not have the flag propagated. Setting
149 that threshold too low means that random sources which are falsely flagged in
150 the inputs will start to dominate. If in doubt, we suggest making this threshold
151 relatively low, but not zero (e.g., 0.1 to 0.2 or so). The more confidence in
152 the quality of the flagging, the lower the threshold can be.
153
154 The relative occurrence accounts for the edge of the field-of-view of
155 the camera, but does not include chip gaps, bad or saturated pixels, etc.
156
157 Parameters
158 ----------
159 butler : `Unknown`
160 Data butler, for retrieving the input source catalogs.
161 coaddSources : `lsst.afw.image.SourceCatalog`
162 Source catalog from the coadd.
163 ccdInputs : `lsst.afw.table.ExposureCatalog`
164 Table of CCDs that contribute to the coadd.
165 coaddWcs : `lsst.afw.geom.SkyWcs`
166 Wcs for coadd.
167 visitCatalogs : `list` of `lsst.afw.image.SourceCatalog`, optional
168 List of loaded source catalogs for each input ccd in
169 the coadd. If provided this is used instead of this
170 method loading in the catalogs itself.
171 wcsUpdates : `list` of `lsst.afw.geom.SkyWcs`, optional
172 If visitCatalogs is a list of ccd catalogs, this
173 should be a list of updated wcs to apply.
174
175 Raises
176 ------
177 ValueError
178 Raised if any of the following occur:
179 - A list of wcs updates for each catalog is not supplied in the wcsUpdates parameter
180 and ccdInputs is a list of src catalogs.
181 - The visitCatalogs and ccdInput parameters are both `None`.
182 """
183 if len(self.config.flags) == 0:
184 return
185
186 flags = self._keys.keys()
187 counts = dict((f, numpy.zeros(len(coaddSources), dtype=int)) for f in flags)
188 indices = numpy.array([s.getId() for s in coaddSources]) # Allowing for non-contiguous data
189 radius = self.config.matchRadius*geom.arcseconds
190
191 def processCcd(ccdSources, wcsUpdate):
192 for sourceRecord in ccdSources:
193 sourceRecord.updateCoord(wcsUpdate)
194 for flag in flags:
195 # We assume that the flags will be relatively rare, so it is more efficient to match
196 # against a subset of the input catalog for each flag than it is to match once against
197 # the entire catalog. It would be best to have built a kd-tree on coaddSources and
198 # keep reusing that for the matching, but we don't have a suitable implementation.
200 mc.findOnlyClosest = False
201 matches = afwTable.matchRaDec(coaddSources, ccdSources[ccdSources.get(flag)], radius, mc)
202 for m in matches:
203 index = (numpy.where(indices == m.first.getId()))[0][0]
204 counts[flag][index] += 1
205
206 if visitCatalogs is not None:
207 if wcsUpdates is None:
208 raise pexExceptions.ValueError("If ccdInputs is a list of src catalogs, a list of wcs"
209 " updates for each catalog must be supplied in the "
210 "wcsUpdates parameter")
211 for i, ccdSource in enumerate(visitCatalogs):
212 processCcd(ccdSource, wcsUpdates[i])
213 else:
214 if ccdInputs is None:
215 raise pexExceptions.ValueError("The visitCatalogs and ccdInput parameters can't both be None")
216 visitKey = ccdInputs.schema.find("visit").key
217 ccdKey = ccdInputs.schema.find("ccd").key
218
219 self.log.info("Propagating flags %s from inputs", flags)
220
221 # Accumulate counts of flags being set
222 for ccdRecord in ccdInputs:
223 v = ccdRecord.get(visitKey)
224 c = ccdRecord.get(ccdKey)
225 dataId = {"visit": int(v), self.config.ccdName: int(c)}
226 ccdSources = butler.get("src", dataId=dataId)
227 processCcd(ccdSources, ccdRecord.getWcs())
228
229 # Apply threshold
230 for f in flags:
231 key = self._keys[f]
232 for s, num in zip(coaddSources, counts[f]):
233 numOverlaps = len(ccdInputs.subsetContaining(s.getCentroid(), coaddWcs, True))
234 s.setFlag(key, bool(num > numOverlaps*self.config.flags[f]))
235 self.log.info("Propagated %d sources with flag %s", sum(s.get(key) for s in coaddSources), f)
Pass parameters to algorithms that match list of sources.
Definition Match.h:45
run(self, butler, coaddSources, ccdInputs, coaddWcs, visitCatalogs=None, wcsUpdates=None)
std::vector< Match< typename Cat1::Record, typename Cat2::Record > > matchRaDec(Cat1 const &cat1, Cat2 const &cat2, lsst::geom::Angle radius, MatchControl const &mc=MatchControl())
Compute all tuples (s1,s2,d) where s1 belings to cat1, s2 belongs to cat2 and d, the distance between...
Definition Match.cc:158