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
measurePsf.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__ = ["MeasurePsfConfig", "MeasurePsfTask"]
23
24import lsst.afw.display as afwDisplay
25import lsst.afw.math as afwMath
26import lsst.meas.algorithms as measAlg
27import lsst.meas.algorithms.utils as maUtils
28import lsst.pex.config as pexConfig
29import lsst.pipe.base as pipeBase
31from lsst.utils.timer import timeMethod
32
33
34class MeasurePsfConfig(pexConfig.Config):
35 starSelector = measAlg.sourceSelectorRegistry.makeField(
36 "Star selection algorithm",
37 default="objectSize"
38 )
39 makePsfCandidates = pexConfig.ConfigurableField(
40 target=measAlg.MakePsfCandidatesTask,
41 doc="Task to make psf candidates from selected stars.",
42 )
43 psfDeterminer = measAlg.psfDeterminerRegistry.makeField(
44 "PSF Determination algorithm",
45 default="psfex"
46 )
47 reserve = pexConfig.ConfigurableField(
48 target=measAlg.ReserveSourcesTask,
49 doc="Reserve sources from fitting"
50 )
51
52 def validate(self):
53 super().validate()
54 if (self.psfDeterminer.name == "piff" and self.psfDeterminer["piff"].stampSize
55 and self.psfDeterminer["piff"].stampSize > self.makePsfCandidates.kernelSize):
56 msg = (f"PIFF kernelSize={self.psfDeterminer['piff'].stampSize}"
57 f" must be >= psf candidate kernelSize={self.makePsfCandidates.kernelSize}.")
58 raise pexConfig.FieldValidationError(MeasurePsfConfig.makePsfCandidates, self, msg)
59
60
61class MeasurePsfTask(pipeBase.Task):
62 """A task that selects stars from a catalog of sources and uses those to measure the PSF.
63
64 Parameters
65 ----------
66 schema : `lsst.sfw.table.Schema`
67 An `lsst.afw.table.Schema` used to create the output `lsst.afw.table.SourceCatalog`.
68 **kwargs :
69 Keyword arguments passed to lsst.pipe.base.task.Task.__init__.
70
71 Notes
72 -----
73 If schema is not None, 'calib_psf_candidate' and 'calib_psf_used' fields will be added to
74 identify which stars were employed in the PSF estimation.
75
76 This task can add fields to the schema, so any code calling this task must ensure that
77 these fields are indeed present in the input table.
78
79 The star selector is a subclass of
80 ``lsst.meas.algorithms.starSelector.BaseStarSelectorTask`` "lsst.meas.algorithms.BaseStarSelectorTask"
81 and the PSF determiner is a sublcass of
82 ``lsst.meas.algorithms.psfDeterminer.BasePsfDeterminerTask`` "lsst.meas.algorithms.BasePsfDeterminerTask"
83
84 There is no establised set of configuration parameters for these algorithms, so once you start modifying
85 parameters (as we do in @ref pipe_tasks_measurePsf_Example) your code is no longer portable.
86
87 Debugging:
88
89 .. code-block:: none
90
91 display
92 If True, display debugging plots
93 displayExposure
94 display the Exposure + spatialCells
95 displayPsfCandidates
96 show mosaic of candidates
97 showBadCandidates
98 Include bad candidates
99 displayPsfMosaic
100 show mosaic of reconstructed PSF(xy)
101 displayResiduals
102 show residuals
103 normalizeResiduals
104 Normalise residuals by object amplitude
105
106 Additionally you can enable any debug outputs that your chosen star selector and psf determiner support.
107 """
108 ConfigClass = MeasurePsfConfig
109 _DefaultName = "measurePsf"
110
111 def __init__(self, schema=None, **kwargs):
112 pipeBase.Task.__init__(self, **kwargs)
113 if schema is not None:
114 self.candidateKey = schema.addField(
115 "calib_psf_candidate", type="Flag",
116 doc=("Flag set if the source was a candidate for PSF determination, "
117 "as determined by the star selector.")
118 )
119 self.usedKey = schema.addField(
120 "calib_psf_used", type="Flag",
121 doc=("Flag set if the source was actually used for PSF determination, "
122 "as determined by the '%s' PSF determiner.") % self.config.psfDeterminer.name
123 )
124 else:
125 self.candidateKey = None
126 self.usedKey = None
127 self.makeSubtask("starSelector")
128 self.makeSubtask("makePsfCandidates")
129 self.makeSubtask("psfDeterminer", schema=schema)
130 self.makeSubtask("reserve", columnName="calib_psf", schema=schema,
131 doc="set if source was reserved from PSF determination")
132
133 @timeMethod
134 def run(self, exposure, sources, expId=0, matches=None):
135 """Measure the PSF.
136
137 Parameters
138 ----------
139 exposure : `lsst.afw.image.Exposure`
140 Exposure to process; measured PSF will be added.
141 sources : `Unknown`
142 Measured sources on exposure; flag fields will be set marking
143 stars chosen by the star selector and the PSF determiner if a schema
144 was passed to the task constructor.
145 expId : `int`, optional
146 Exposure id used for generating random seed.
147 matches : `list`, optional
148 A list of ``lsst.afw.table.ReferenceMatch`` objects
149 (i.e. of ``lsst.afw.table.Match`` with @c first being
150 of type ``lsst.afw.table.SimpleRecord`` and @c second
151 type lsst.afw.table.SourceRecord --- the reference object and detected
152 object respectively) as returned by @em e.g. the AstrometryTask.
153 Used by star selectors that choose to refer to an external catalog.
154
155 Returns
156 -------
157 measurement : `lsst.pipe.base.Struct`
158 PSF measurement as a struct with attributes:
159
160 ``psf``
161 The measured PSF (also set in the input exposure).
162 ``cellSet``
163 An `lsst.afw.math.SpatialCellSet` containing the PSF candidates
164 as returned by the psf determiner.
165 """
166 self.log.info("Measuring PSF")
167
168 import lsstDebug
169 display = lsstDebug.Info(__name__).display
170 displayExposure = lsstDebug.Info(__name__).displayExposure # display the Exposure + spatialCells
171 displayPsfMosaic = lsstDebug.Info(__name__).displayPsfMosaic # show mosaic of reconstructed PSF(x,y)
172 displayPsfCandidates = lsstDebug.Info(__name__).displayPsfCandidates # show mosaic of candidates
173 displayResiduals = lsstDebug.Info(__name__).displayResiduals # show residuals
174 showBadCandidates = lsstDebug.Info(__name__).showBadCandidates # include bad candidates
175 normalizeResiduals = lsstDebug.Info(__name__).normalizeResiduals # normalise residuals by object peak
176
177 #
178 # Run star selector
179 #
180 stars = self.starSelector.run(sourceCat=sources, matches=matches, exposure=exposure)
181 selectionResult = self.makePsfCandidates.run(stars.sourceCat, exposure=exposure)
182 self.log.info("PSF star selector found %d candidates", len(selectionResult.psfCandidates))
183 reserveResult = self.reserve.run(selectionResult.goodStarCat, expId=expId)
184 # Make list of psf candidates to send to the determiner (omitting those marked as reserved)
185 psfDeterminerList = [cand for cand, use
186 in zip(selectionResult.psfCandidates, reserveResult.use) if use]
187
188 if selectionResult.psfCandidates and self.candidateKey is not None:
189 for cand in selectionResult.psfCandidates:
190 source = cand.getSource()
191 source.set(self.candidateKey, True)
192
193 self.log.info("Sending %d candidates to PSF determiner", len(psfDeterminerList))
194
195 if display:
196 frame = 1
197 if displayExposure:
198 disp = afwDisplay.Display(frame=frame)
199 disp.mtv(exposure, title="psf determination")
200 frame += 1
201 #
202 # Determine PSF
203 #
204 psf, cellSet = self.psfDeterminer.determinePsf(exposure, psfDeterminerList, self.metadata,
205 flagKey=self.usedKey)
206 self.log.info("PSF determination using %d/%d stars.",
207 self.metadata.getScalar("numGoodStars"), self.metadata.getScalar("numAvailStars"))
208
209 exposure.setPsf(psf)
210
211 if display:
212 frame = display
213 if displayExposure:
214 disp = afwDisplay.Display(frame=frame)
215 showPsfSpatialCells(exposure, cellSet, showBadCandidates, frame=frame)
216 frame += 1
217
218 if displayPsfCandidates: # Show a mosaic of PSF candidates
219 plotPsfCandidates(cellSet, showBadCandidates=showBadCandidates, frame=frame)
220 frame += 1
221
222 if displayResiduals:
223 frame = plotResiduals(exposure, cellSet,
224 showBadCandidates=showBadCandidates,
225 normalizeResiduals=normalizeResiduals,
226 frame=frame)
227 if displayPsfMosaic:
228 disp = afwDisplay.Display(frame=frame)
229 maUtils.showPsfMosaic(exposure, psf, display=disp, showFwhm=True)
230 disp.scale("linear", 0, 1)
231 frame += 1
232
233 return pipeBase.Struct(
234 psf=psf,
235 cellSet=cellSet,
236 )
237
238 @property
239 def usesMatches(self):
240 """Return True if this task makes use of the "matches" argument to the run method"""
241 return self.starSelector.usesMatches
242
243#
244# Debug code
245#
246
247
248def showPsfSpatialCells(exposure, cellSet, showBadCandidates, frame=1):
249 disp = afwDisplay.Display(frame=frame)
250 maUtils.showPsfSpatialCells(exposure, cellSet,
251 symb="o", ctype=afwDisplay.CYAN, ctypeUnused=afwDisplay.YELLOW,
252 size=4, display=disp)
253 for cell in cellSet.getCellList():
254 for cand in cell.begin(not showBadCandidates): # maybe include bad candidates
255 status = cand.getStatus()
256 disp.dot('+', *cand.getSource().getCentroid(),
257 ctype=afwDisplay.GREEN if status == afwMath.SpatialCellCandidate.GOOD else
258 afwDisplay.YELLOW if status == afwMath.SpatialCellCandidate.UNKNOWN else afwDisplay.RED)
259
260
261def plotPsfCandidates(cellSet, showBadCandidates=False, frame=1):
262 stamps = []
263 for cell in cellSet.getCellList():
264 for cand in cell.begin(not showBadCandidates): # maybe include bad candidates
265 try:
266 im = cand.getMaskedImage()
267
268 chi2 = cand.getChi2()
269 if chi2 < 1e100:
270 chi2 = "%.1f" % chi2
271 else:
272 chi2 = float("nan")
273
274 stamps.append((im, "%d%s" %
275 (maUtils.splitId(cand.getSource().getId(), True)["objId"], chi2),
276 cand.getStatus()))
277 except Exception:
278 continue
279
280 mos = afwDisplay.utils.Mosaic()
281 disp = afwDisplay.Display(frame=frame)
282 for im, label, status in stamps:
283 im = type(im)(im, True)
284 try:
285 im /= afwMath.makeStatistics(im, afwMath.MAX).getValue()
286 except NotImplementedError:
287 pass
288
289 mos.append(im, label,
290 afwDisplay.GREEN if status == afwMath.SpatialCellCandidate.GOOD else
291 afwDisplay.YELLOW if status == afwMath.SpatialCellCandidate.UNKNOWN else afwDisplay.RED)
292
293 if mos.images:
294 disp.mtv(mos.makeMosaic(), title="Psf Candidates")
295
296
297def plotResiduals(exposure, cellSet, showBadCandidates=False, normalizeResiduals=True, frame=2):
298 psf = exposure.getPsf()
299 disp = afwDisplay.Display(frame=frame)
300 while True:
301 try:
302 maUtils.showPsfCandidates(exposure, cellSet, psf=psf, display=disp,
303 normalize=normalizeResiduals,
304 showBadCandidates=showBadCandidates)
305 frame += 1
306 maUtils.showPsfCandidates(exposure, cellSet, psf=psf, display=disp,
307 normalize=normalizeResiduals,
308 showBadCandidates=showBadCandidates,
309 variance=True)
310 frame += 1
311 except Exception:
312 if not showBadCandidates:
313 showBadCandidates = True
314 continue
315 break
316
317 return frame
__init__(self, schema=None, **kwargs)
run(self, exposure, sources, expId=0, matches=None)
Statistics makeStatistics(lsst::afw::image::Image< Pixel > const &img, lsst::afw::image::Mask< image::MaskPixel > const &msk, int const flags, StatisticsControl const &sctrl=StatisticsControl())
Handle a watered-down front-end to the constructor (no variance)
Definition Statistics.h:361
plotPsfCandidates(cellSet, showBadCandidates=False, frame=1)
plotResiduals(exposure, cellSet, showBadCandidates=False, normalizeResiduals=True, frame=2)
showPsfSpatialCells(exposure, cellSet, showBadCandidates, frame=1)