LSST Applications g2079a07aa2+86d27d4dc4,g2305ad1205+a659bff248,g2bbee38e9b+3c60f8fe34,g337abbeb29+3c60f8fe34,g33d1c0ed96+3c60f8fe34,g3502564af9+d77d6d1350,g3a166c0a6a+3c60f8fe34,g487adcacf7+25d9892218,g4be5004598+d77d6d1350,g50ff169b8f+96c6868917,g52b1c1532d+585e252eca,g591dd9f2cf+4d81263f9a,g5cd07815a0+980d2b1c3b,g607f77f49a+d77d6d1350,g858d7b2824+d77d6d1350,g88963caddf+83e433e629,g99cad8db69+a4d3c48eeb,g9ddcbc5298+9a081db1e4,ga1e77700b3+bcf1af89ad,ga57fefb910+9a39d7b2d7,gae0086650b+585e252eca,gb065fddaf9+4f9fd82a2c,gb0e22166c9+60f28cb32d,gb363559e06+d84b1d3d07,gb3b7280ab2+4563d032e1,gb4b16eec92+babe958938,gba4ed39666+c2a2e4ac27,gbb8dafda3b+ed6854b564,gc120e1dc64+b72d212f87,gc28159a63d+3c60f8fe34,gc3e9b769f7+921dbcd359,gcf0d15dbbd+9a39d7b2d7,gdaeeff99f8+f9a426f77a,gddc38dedce+585e252eca,ge79ae78c31+3c60f8fe34,w.2024.21
LSST Data Management Base Package
Loading...
Searching...
No Matches
computeExposureSummaryStats.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__ = ["ComputeExposureSummaryStatsTask", "ComputeExposureSummaryStatsConfig"]
23
24import warnings
25import numpy as np
26from scipy.stats import median_abs_deviation as sigmaMad
27import astropy.units as units
28from astropy.time import Time
29from astropy.coordinates import AltAz, SkyCoord, EarthLocation
30from lsst.daf.base import DateTime
31
32import lsst.pipe.base as pipeBase
33import lsst.pex.config as pexConfig
34import lsst.afw.math as afwMath
35import lsst.afw.image as afwImage
36import lsst.geom as geom
37from lsst.meas.algorithms import ScienceSourceSelectorTask
38from lsst.utils.timer import timeMethod
39
40
41class ComputeExposureSummaryStatsConfig(pexConfig.Config):
42 """Config for ComputeExposureSummaryTask"""
43 sigmaClip = pexConfig.Field(
44 dtype=float,
45 doc="Sigma for outlier rejection for sky noise.",
46 default=3.0,
47 )
48 clipIter = pexConfig.Field(
49 dtype=int,
50 doc="Number of iterations of outlier rejection for sky noise.",
51 default=2,
52 )
53 badMaskPlanes = pexConfig.ListField(
54 dtype=str,
55 doc="Mask planes that, if set, the associated pixel should not be included sky noise calculation.",
56 default=("NO_DATA", "SUSPECT"),
57 )
58 starSelection = pexConfig.Field(
59 doc="Field to select full list of sources used for PSF modeling.",
60 dtype=str,
61 default="calib_psf_used",
62 )
63 starSelector = pexConfig.ConfigurableField(
64 target=ScienceSourceSelectorTask,
65 doc="Selection of sources to compute PSF star statistics.",
66 )
67 starShape = pexConfig.Field(
68 doc="Base name of columns to use for the source shape in the PSF statistics computation.",
69 dtype=str,
70 default="slot_Shape"
71 )
72 psfShape = pexConfig.Field(
73 doc="Base name of columns to use for the PSF shape in the PSF statistics computation.",
74 dtype=str,
75 default="slot_PsfShape"
76 )
77 psfSampling = pexConfig.Field(
78 dtype=int,
79 doc="Sampling rate in pixels in each dimension for the maxDistToNearestPsf metric "
80 "caclulation grid (the tradeoff is between adequate sampling versus speed).",
81 default=8,
82 )
83 psfGridSampling = pexConfig.Field(
84 dtype=int,
85 doc="Sampling rate in pixels in each dimension for PSF model robustness metric "
86 "caclulations grid (the tradeoff is between adequate sampling versus speed).",
87 default=96,
88 )
89 psfBadMaskPlanes = pexConfig.ListField(
90 dtype=str,
91 doc="Mask planes that, if set, the associated pixel should not be included in the PSF model "
92 "robutsness metric calculations (namely, maxDistToNearestPsf and psfTraceRadiusDelta).",
93 default=("BAD", "CR", "EDGE", "INTRP", "NO_DATA", "SAT", "SUSPECT"),
94 )
95 fiducialSkyBackground = pexConfig.DictField(
96 keytype=str,
97 itemtype=float,
98 doc="Fiducial sky background level (ADU/s) assumed when calculating effective exposure time. "
99 "Keyed by band.",
100 default={'u': 1.0, 'g': 1.0, 'r': 1.0, 'i': 1.0, 'z': 1.0, 'y': 1.0},
101 )
102 fiducialPsfSigma = pexConfig.DictField(
103 keytype=str,
104 itemtype=float,
105 doc="Fiducial PSF sigma (pixels) assumed when calculating effective exposure time. "
106 "Keyed by band.",
107 default={'u': 1.0, 'g': 1.0, 'r': 1.0, 'i': 1.0, 'z': 1.0, 'y': 1.0},
108 )
109 fiducialZeroPoint = pexConfig.DictField(
110 keytype=str,
111 itemtype=float,
112 doc="Fiducial zero point assumed when calculating effective exposure time. "
113 "Keyed by band.",
114 default={'u': 25.0, 'g': 25.0, 'r': 25.0, 'i': 25.0, 'z': 25.0, 'y': 25.0},
115 )
116 maxEffectiveTransparency = pexConfig.Field(
117 dtype=float,
118 doc="Maximum value allowed for effective transparency scale factor (often inf or 1.0).",
119 default=float('inf')
120 )
121
122 def setDefaults(self):
123 super().setDefaults()
124
126 self.starSelector.doFlags = True
127 self.starSelector.doSignalToNoise = True
128 self.starSelector.doUnresolved = False
129 self.starSelector.doIsolated = False
130 self.starSelector.doRequireFiniteRaDec = False
131 self.starSelector.doRequirePrimary = False
132
133 self.starSelector.signalToNoise.minimum = 50.0
134 self.starSelector.signalToNoise.maximum = 1000.0
135
136 self.starSelector.flags.bad = ["slot_Shape_flag", "slot_PsfFlux_flag"]
137 # Select stars used for PSF modeling.
138 self.starSelector.flags.good = ["calib_psf_used"]
139
140 self.starSelector.signalToNoise.fluxField = "slot_PsfFlux_instFlux"
141 self.starSelector.signalToNoise.errField = "slot_PsfFlux_instFluxErr"
142
143
145 """Task to compute exposure summary statistics.
146
147 This task computes various quantities suitable for DPDD and other
148 downstream processing at the detector centers, including:
149 - psfSigma
150 - psfArea
151 - psfIxx
152 - psfIyy
153 - psfIxy
154 - ra
155 - dec
156 - zenithDistance
157 - zeroPoint
158 - skyBg
159 - skyNoise
160 - meanVar
161 - raCorners
162 - decCorners
163 - astromOffsetMean
164 - astromOffsetStd
165
166 These additional quantities are computed from the stars in the detector:
167 - psfStarDeltaE1Median
168 - psfStarDeltaE2Median
169 - psfStarDeltaE1Scatter
170 - psfStarDeltaE2Scatter
171 - psfStarDeltaSizeMedian
172 - psfStarDeltaSizeScatter
173 - psfStarScaledDeltaSizeScatter
174
175 These quantities are computed based on the PSF model and image mask
176 to assess the robustness of the PSF model across a given detector
177 (against, e.g., extrapolation instability):
178 - maxDistToNearestPsf
179 - psfTraceRadiusDelta
180
181 These quantities are computed to assess depth:
182 - effTime
183 - effTimePsfSigmaScale
184 - effTimeSkyBgScale
185 - effTimeZeroPointScale
186 """
187 ConfigClass = ComputeExposureSummaryStatsConfig
188 _DefaultName = "computeExposureSummaryStats"
189
190 def __init__(self, **kwargs):
191 super().__init__(**kwargs)
192
193 self.makeSubtask("starSelector")
194
195 @timeMethod
196 def run(self, exposure, sources, background):
197 """Measure exposure statistics from the exposure, sources, and
198 background.
199
200 Parameters
201 ----------
202 exposure : `lsst.afw.image.ExposureF`
203 sources : `lsst.afw.table.SourceCatalog`
204 background : `lsst.afw.math.BackgroundList`
205
206 Returns
207 -------
208 summary : `lsst.afw.image.ExposureSummary`
209 """
210 self.log.info("Measuring exposure statistics")
211
213
214 bbox = exposure.getBBox()
215
216 psf = exposure.getPsf()
217 self.update_psf_stats(summary, psf, bbox, sources, image_mask=exposure.mask)
218
219 wcs = exposure.getWcs()
220 visitInfo = exposure.getInfo().getVisitInfo()
221 self.update_wcs_stats(summary, wcs, bbox, visitInfo)
222
223 photoCalib = exposure.getPhotoCalib()
224 self.update_photo_calib_stats(summary, photoCalib)
225
226 self.update_background_stats(summary, background)
227
228 self.update_masked_image_stats(summary, exposure.getMaskedImage())
229
230 self.update_effective_time_stats(summary, exposure)
231
232 md = exposure.getMetadata()
233 if 'SFM_ASTROM_OFFSET_MEAN' in md:
234 summary.astromOffsetMean = md['SFM_ASTROM_OFFSET_MEAN']
235 summary.astromOffsetStd = md['SFM_ASTROM_OFFSET_STD']
236
237 return summary
238
239 def update_psf_stats(self, summary, psf, bbox, sources=None, image_mask=None, sources_is_astropy=False):
240 """Compute all summary-statistic fields that depend on the PSF model.
241
242 Parameters
243 ----------
244 summary : `lsst.afw.image.ExposureSummaryStats`
245 Summary object to update in-place.
246 psf : `lsst.afw.detection.Psf` or `None`
247 Point spread function model. If `None`, all fields that depend on
248 the PSF will be reset (generally to NaN).
249 bbox : `lsst.geom.Box2I`
250 Bounding box of the image for which summary stats are being
251 computed.
252 sources : `lsst.afw.table.SourceCatalog` or `astropy.table.Table`
253 Catalog for quantities that are computed from source table columns.
254 If `None`, these quantities will be reset (generally to NaN).
255 The type of this table must correspond to the
256 ``sources_is_astropy`` argument.
257 image_mask : `lsst.afw.image.Mask`, optional
258 Mask image that may be used to compute distance-to-nearest-star
259 metrics.
260 sources_is_astropy : `bool`, optional
261 Whether ``sources`` is an `astropy.table.Table` instance instead
262 of an `lsst.afw.table.Catalog` instance. Default is `False` (the
263 latter).
264 """
265 nan = float("nan")
266 summary.psfSigma = nan
267 summary.psfIxx = nan
268 summary.psfIyy = nan
269 summary.psfIxy = nan
270 summary.psfArea = nan
271 summary.nPsfStar = 0
272 summary.psfStarDeltaE1Median = nan
273 summary.psfStarDeltaE2Median = nan
274 summary.psfStarDeltaE1Scatter = nan
275 summary.psfStarDeltaE2Scatter = nan
276 summary.psfStarDeltaSizeMedian = nan
277 summary.psfStarDeltaSizeScatter = nan
278 summary.psfStarScaledDeltaSizeScatter = nan
279 summary.maxDistToNearestPsf = nan
280 summary.psfTraceRadiusDelta = nan
281
282 if psf is None:
283 return
284 shape = psf.computeShape(bbox.getCenter())
285 summary.psfSigma = shape.getDeterminantRadius()
286 summary.psfIxx = shape.getIxx()
287 summary.psfIyy = shape.getIyy()
288 summary.psfIxy = shape.getIxy()
289 im = psf.computeKernelImage(bbox.getCenter())
290 # The calculation of effective psf area is taken from
291 # meas_base/src/PsfFlux.cc#L112. See
292 # https://github.com/lsst/meas_base/blob/
293 # 750bffe6620e565bda731add1509507f5c40c8bb/src/PsfFlux.cc#L112
294 summary.psfArea = float(np.sum(im.array)/np.sum(im.array**2.))
295
296 if image_mask is not None:
297 psfTraceRadiusDelta = psf_trace_radius_delta(
298 image_mask,
299 psf,
300 sampling=self.config.psfGridSampling,
301 bad_mask_bits=self.config.psfBadMaskPlanes
302 )
303 summary.psfTraceRadiusDelta = float(psfTraceRadiusDelta)
304
305 if sources is None:
306 # No sources are available (as in some tests and rare cases where
307 # the selection criteria in finalizeCharacterization lead to no
308 # good sources).
309 return
310
311 # Count the total number of psf stars used (prior to stats selection).
312 nPsfStar = sources[self.config.starSelection].sum()
313 summary.nPsfStar = int(nPsfStar)
314
315 psf_mask = self.starSelector.run(sources).selected
316 nPsfStarsUsedInStats = psf_mask.sum()
317
318 if nPsfStarsUsedInStats == 0:
319 # No stars to measure statistics, so we must return the defaults
320 # of 0 stars and NaN values.
321 return
322
323 if sources_is_astropy:
324 psf_cat = sources[psf_mask]
325 else:
326 psf_cat = sources[psf_mask].copy(deep=True)
327
328 starXX = psf_cat[self.config.starShape + '_xx']
329 starYY = psf_cat[self.config.starShape + '_yy']
330 starXY = psf_cat[self.config.starShape + '_xy']
331 psfXX = psf_cat[self.config.psfShape + '_xx']
332 psfYY = psf_cat[self.config.psfShape + '_yy']
333 psfXY = psf_cat[self.config.psfShape + '_xy']
334
335 # Use the trace radius for the star size.
336 starSize = np.sqrt(starXX/2. + starYY/2.)
337
338 starE1 = (starXX - starYY)/(starXX + starYY)
339 starE2 = 2*starXY/(starXX + starYY)
340 starSizeMedian = np.median(starSize)
341
342 # Use the trace radius for the psf size.
343 psfSize = np.sqrt(psfXX/2. + psfYY/2.)
344 psfE1 = (psfXX - psfYY)/(psfXX + psfYY)
345 psfE2 = 2*psfXY/(psfXX + psfYY)
346
347 psfStarDeltaE1Median = np.median(starE1 - psfE1)
348 psfStarDeltaE1Scatter = sigmaMad(starE1 - psfE1, scale='normal')
349 psfStarDeltaE2Median = np.median(starE2 - psfE2)
350 psfStarDeltaE2Scatter = sigmaMad(starE2 - psfE2, scale='normal')
351
352 psfStarDeltaSizeMedian = np.median(starSize - psfSize)
353 psfStarDeltaSizeScatter = sigmaMad(starSize - psfSize, scale='normal')
354 psfStarScaledDeltaSizeScatter = psfStarDeltaSizeScatter/starSizeMedian
355
356 summary.psfStarDeltaE1Median = float(psfStarDeltaE1Median)
357 summary.psfStarDeltaE2Median = float(psfStarDeltaE2Median)
358 summary.psfStarDeltaE1Scatter = float(psfStarDeltaE1Scatter)
359 summary.psfStarDeltaE2Scatter = float(psfStarDeltaE2Scatter)
360 summary.psfStarDeltaSizeMedian = float(psfStarDeltaSizeMedian)
361 summary.psfStarDeltaSizeScatter = float(psfStarDeltaSizeScatter)
362 summary.psfStarScaledDeltaSizeScatter = float(psfStarScaledDeltaSizeScatter)
363
364 if image_mask is not None:
365 maxDistToNearestPsf = maximum_nearest_psf_distance(
366 image_mask,
367 psf_cat,
368 sampling=self.config.psfSampling,
369 bad_mask_bits=self.config.psfBadMaskPlanes
370 )
371 summary.maxDistToNearestPsf = float(maxDistToNearestPsf)
372
373 def update_wcs_stats(self, summary, wcs, bbox, visitInfo):
374 """Compute all summary-statistic fields that depend on the WCS model.
375
376 Parameters
377 ----------
378 summary : `lsst.afw.image.ExposureSummaryStats`
379 Summary object to update in-place.
380 wcs : `lsst.afw.geom.SkyWcs` or `None`
381 Astrometric calibration model. If `None`, all fields that depend
382 on the WCS will be reset (generally to NaN).
383 bbox : `lsst.geom.Box2I`
384 Bounding box of the image for which summary stats are being
385 computed.
386 visitInfo : `lsst.afw.image.VisitInfo`
387 Observation information used in together with ``wcs`` to compute
388 the zenith distance.
389 """
390 nan = float("nan")
391 summary.raCorners = [nan]*4
392 summary.decCorners = [nan]*4
393 summary.ra = nan
394 summary.dec = nan
395 summary.zenithDistance = nan
396
397 if wcs is None:
398 return
399
400 sph_pts = wcs.pixelToSky(geom.Box2D(bbox).getCorners())
401 summary.raCorners = [float(sph.getRa().asDegrees()) for sph in sph_pts]
402 summary.decCorners = [float(sph.getDec().asDegrees()) for sph in sph_pts]
403
404 sph_pt = wcs.pixelToSky(bbox.getCenter())
405 summary.ra = sph_pt.getRa().asDegrees()
406 summary.dec = sph_pt.getDec().asDegrees()
407
408 date = visitInfo.getDate()
409
410 if date.isValid():
411 # We compute the zenithDistance at the center of the detector
412 # rather than use the boresight value available via the visitInfo,
413 # because the zenithDistance may vary significantly over a large
414 # field of view.
415 observatory = visitInfo.getObservatory()
416 loc = EarthLocation(lat=observatory.getLatitude().asDegrees()*units.deg,
417 lon=observatory.getLongitude().asDegrees()*units.deg,
418 height=observatory.getElevation()*units.m)
419 obstime = Time(visitInfo.getDate().get(system=DateTime.MJD),
420 location=loc, format='mjd')
421 coord = SkyCoord(
422 summary.ra*units.degree,
423 summary.dec*units.degree,
424 obstime=obstime,
425 location=loc,
426 )
427 with warnings.catch_warnings():
428 warnings.simplefilter('ignore')
429 altaz = coord.transform_to(AltAz)
430
431 summary.zenithDistance = float(90.0 - altaz.alt.degree)
432
433 def update_photo_calib_stats(self, summary, photo_calib):
434 """Compute all summary-statistic fields that depend on the photometric
435 calibration model.
436
437 Parameters
438 ----------
439 summary : `lsst.afw.image.ExposureSummaryStats`
440 Summary object to update in-place.
441 photo_calib : `lsst.afw.image.PhotoCalib` or `None`
442 Photometric calibration model. If `None`, all fields that depend
443 on the photometric calibration will be reset (generally to NaN).
444 """
445 if photo_calib is not None:
446 summary.zeroPoint = float(2.5*np.log10(photo_calib.getInstFluxAtZeroMagnitude()))
447 else:
448 summary.zeroPoint = float("nan")
449
450 def update_background_stats(self, summary, background):
451 """Compute summary-statistic fields that depend only on the
452 background model.
453
454 Parameters
455 ----------
456 summary : `lsst.afw.image.ExposureSummaryStats`
457 Summary object to update in-place.
458 background : `lsst.afw.math.BackgroundList` or `None`
459 Background model. If `None`, all fields that depend on the
460 background will be reset (generally to NaN).
461
462 Notes
463 -----
464 This does not include fields that depend on the background-subtracted
465 masked image; when the background changes, it should generally be
466 applied to the image and `update_masked_image_stats` should be called
467 as well.
468 """
469 if background is not None:
470 bgStats = (bg[0].getStatsImage().getImage().array
471 for bg in background)
472 summary.skyBg = float(sum(np.median(bg[np.isfinite(bg)]) for bg in bgStats))
473 else:
474 summary.skyBg = float("nan")
475
476 def update_masked_image_stats(self, summary, masked_image):
477 """Compute summary-statistic fields that depend on the masked image
478 itself.
479
480 Parameters
481 ----------
482 summary : `lsst.afw.image.ExposureSummaryStats`
483 Summary object to update in-place.
484 masked_image : `lsst.afw.image.MaskedImage` or `None`
485 Masked image. If `None`, all fields that depend
486 on the masked image will be reset (generally to NaN).
487 """
488 nan = float("nan")
489 if masked_image is None:
490 summary.skyNoise = nan
491 summary.meanVar = nan
492 return
493 statsCtrl = afwMath.StatisticsControl()
494 statsCtrl.setNumSigmaClip(self.config.sigmaClip)
495 statsCtrl.setNumIter(self.config.clipIter)
496 statsCtrl.setAndMask(afwImage.Mask.getPlaneBitMask(self.config.badMaskPlanes))
497 statsCtrl.setNanSafe(True)
498
499 statObj = afwMath.makeStatistics(masked_image, afwMath.STDEVCLIP, statsCtrl)
500 skyNoise, _ = statObj.getResult(afwMath.STDEVCLIP)
501 summary.skyNoise = skyNoise
502
503 statObj = afwMath.makeStatistics(masked_image.variance, masked_image.mask, afwMath.MEANCLIP,
504 statsCtrl)
505 meanVar, _ = statObj.getResult(afwMath.MEANCLIP)
506 summary.meanVar = meanVar
507
508 def update_effective_time_stats(self, summary, exposure):
509 """Compute effective exposure time statistics to estimate depth.
510
511 The effective exposure time is the equivalent shutter open
512 time that would be needed under nominal conditions to give the
513 same signal-to-noise for a point source as what is achieved by
514 the observation of interest. This metric combines measurements
515 of the point-spread function, the sky brightness, and the
516 transparency.
517
518 .. _teff_definitions:
519
520 The effective exposure time and its subcomponents are defined in [1]_
521
522 References
523 ----------
524
525 .. [1] Neilsen, E.H., Bernstein, G., Gruendl, R., and Kent, S. (2016).
526 Limiting Magnitude, \tau, teff, and Image Quality in DES Year 1
527 https://www.osti.gov/biblio/1250877/
528
529
530 Parameters
531 ----------
532 summary : `lsst.afw.image.ExposureSummaryStats`
533 Summary object to update in-place.
534 exposure : `lsst.afw.image.ExposureF`
535 Exposure to grab band and exposure time metadata
536
537 """
538 self.log.info("Updating effective exposure time")
539
540 nan = float("nan")
541 summary.effTime = nan
542 summary.effTimePsfSigmaScale = nan
543 summary.effTimeSkyBgScale = nan
544 summary.effTimeZeroPointScale = nan
545
546 exposureTime = exposure.getInfo().getVisitInfo().getExposureTime()
547 filterLabel = exposure.getFilter()
548 if (filterLabel is None) or (not filterLabel.hasBandLabel):
549 band = None
550 else:
551 band = filterLabel.bandLabel
552
553 if band is None:
554 self.log.warn("No band associated with exposure; effTime not calculated.")
555 return
556
557 # PSF component
558 if np.isnan(summary.psfSigma):
559 self.log.debug("PSF sigma is NaN")
560 f_eff = nan
561 elif band not in self.config.fiducialPsfSigma:
562 self.log.debug(f"Fiducial PSF value not found for {band}")
563 f_eff = nan
564 else:
565 fiducialPsfSigma = self.config.fiducialPsfSigma[band]
566 f_eff = (summary.psfSigma / fiducialPsfSigma)**-2
567
568 # Transparency component (note that exposure time may be removed from zeropoint)
569 if np.isnan(summary.zeroPoint):
570 self.log.debug("Zero point is NaN")
571 c_eff = nan
572 elif band not in self.config.fiducialZeroPoint:
573 self.log.debug(f"Fiducial zero point value not found for {band}")
574 c_eff = nan
575 else:
576 fiducialZeroPoint = self.config.fiducialZeroPoint[band]
577 zeroPointDiff = fiducialZeroPoint - (summary.zeroPoint - 2.5*np.log10(exposureTime))
578 c_eff = min(10**(-2.0*(zeroPointDiff)/2.5), self.config.maxEffectiveTransparency)
579
580 # Sky brightness component (convert to cts/s)
581 if np.isnan(summary.skyBg):
582 self.log.debug("Sky background is NaN")
583 b_eff = nan
584 elif band not in self.config.fiducialSkyBackground:
585 self.log.debug(f"Fiducial sky background value not found for {band}")
586 b_eff = nan
587 else:
588 fiducialSkyBackground = self.config.fiducialSkyBackground[band]
589 b_eff = fiducialSkyBackground/(summary.skyBg/exposureTime)
590
591 # Effective exposure time scale factor
592 t_eff = f_eff * c_eff * b_eff
593
594 # Effective exposure time (seconds)
595 effectiveTime = t_eff * exposureTime
596
597 # Output quantities
598 summary.effTime = float(effectiveTime)
599 summary.effTimePsfSigmaScale = float(f_eff)
600 summary.effTimeSkyBgScale = float(b_eff)
601 summary.effTimeZeroPointScale = float(c_eff)
602
603
605 image_mask,
606 psf_cat,
607 sampling=8,
608 bad_mask_bits=["BAD", "CR", "INTRP", "SAT", "SUSPECT", "NO_DATA", "EDGE"],
609):
610 """Compute the maximum distance of an unmasked pixel to its nearest PSF.
611
612 Parameters
613 ----------
614 image_mask : `lsst.afw.image.Mask`
615 The mask plane associated with the exposure.
616 psf_cat : `lsst.afw.table.SourceCatalog` or `astropy.table.Table`
617 Catalog containing only the stars used in the PSF modeling.
618 sampling : `int`
619 Sampling rate in each dimension to create the grid of points on which
620 to evaluate the distance to the nearest PSF star. The tradeoff is
621 between adequate sampling versus speed.
622 bad_mask_bits : `list` [`str`]
623 Mask bits required to be absent for a pixel to be considered
624 "unmasked".
625
626 Returns
627 -------
628 max_dist_to_nearest_psf : `float`
629 The maximum distance (in pixels) of an unmasked pixel to its nearest
630 PSF model star.
631 """
632 mask_arr = image_mask.array[::sampling, ::sampling]
633 bitmask = image_mask.getPlaneBitMask(bad_mask_bits)
634 good = ((mask_arr & bitmask) == 0)
635
636 x = np.arange(good.shape[1]) * sampling
637 y = np.arange(good.shape[0]) * sampling
638 xx, yy = np.meshgrid(x, y)
639
640 dist_to_nearest_psf = np.full(good.shape, np.inf)
641 for psf in psf_cat:
642 x_psf = psf["slot_Centroid_x"]
643 y_psf = psf["slot_Centroid_y"]
644 dist_to_nearest_psf = np.minimum(dist_to_nearest_psf, np.hypot(xx - x_psf, yy - y_psf))
645 unmasked_dists = dist_to_nearest_psf * good
646 max_dist_to_nearest_psf = np.max(unmasked_dists)
647
648 return max_dist_to_nearest_psf
649
650
652 image_mask,
653 image_psf,
654 sampling=96,
655 bad_mask_bits=["BAD", "CR", "INTRP", "SAT", "SUSPECT", "NO_DATA", "EDGE"],
656):
657 """Compute the delta between the maximum and minimum model PSF trace radius
658 values evaluated on a grid of points lying in the unmasked region of the
659 image.
660
661 Parameters
662 ----------
663 image_mask : `lsst.afw.image.Mask`
664 The mask plane associated with the exposure.
665 image_psf : `lsst.afw.detection.Psf`
666 The PSF model associated with the exposure.
667 sampling : `int`
668 Sampling rate in each dimension to create the grid of points at which
669 to evaluate ``image_psf``s trace radius value. The tradeoff is between
670 adequate sampling versus speed.
671 bad_mask_bits : `list` [`str`]
672 Mask bits required to be absent for a pixel to be considered
673 "unmasked".
674
675 Returns
676 -------
677 psf_trace_radius_delta : `float`
678 The delta (in pixels) between the maximum and minimum model PSF trace
679 radius values evaluated on the x,y-grid subsampled on the unmasked
680 detector pixels by a factor of ``sampling``. If any model PSF trace
681 radius value on the grid evaluates to NaN, then NaN is returned
682 immediately.
683 """
684 psf_trace_radius_list = []
685 mask_arr = image_mask.array[::sampling, ::sampling]
686 bitmask = image_mask.getPlaneBitMask(bad_mask_bits)
687 good = ((mask_arr & bitmask) == 0)
688
689 x = np.arange(good.shape[1]) * sampling
690 y = np.arange(good.shape[0]) * sampling
691 xx, yy = np.meshgrid(x, y)
692
693 for x_mesh, y_mesh, good_mesh in zip(xx, yy, good):
694 for x_point, y_point, is_good in zip(x_mesh, y_mesh, good_mesh):
695 if is_good:
696 psf_trace_radius = image_psf.computeShape(geom.Point2D(x_point, y_point)).getTraceRadius()
697 if ~np.isfinite(psf_trace_radius):
698 return float("nan")
699 psf_trace_radius_list.append(psf_trace_radius)
700
701 psf_trace_radius_delta = np.max(psf_trace_radius_list) - np.min(psf_trace_radius_list)
702
703 return psf_trace_radius_delta
int min
Pass parameters to a Statistics object.
Definition Statistics.h:83
A floating-point coordinate rectangle geometry.
Definition Box.h:413
update_psf_stats(self, summary, psf, bbox, sources=None, image_mask=None, sources_is_astropy=False)
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
psf_trace_radius_delta(image_mask, image_psf, sampling=96, bad_mask_bits=["BAD", "CR", "INTRP", "SAT", "SUSPECT", "NO_DATA", "EDGE"])
maximum_nearest_psf_distance(image_mask, psf_cat, sampling=8, bad_mask_bits=["BAD", "CR", "INTRP", "SAT", "SUSPECT", "NO_DATA", "EDGE"])