LSST Applications g0603fd7c41+022847dfd1,g0aad566f14+f45185db35,g180d380827+40e913b07a,g2079a07aa2+86d27d4dc4,g2305ad1205+696e5f3872,g2bbee38e9b+047b288a59,g337abbeb29+047b288a59,g33d1c0ed96+047b288a59,g3a166c0a6a+047b288a59,g3d1719c13e+f45185db35,g3de15ee5c7+5201731f0d,g487adcacf7+19f9b77d7d,g50ff169b8f+96c6868917,g52b1c1532d+585e252eca,g591dd9f2cf+248b16177b,g63cd9335cc+585e252eca,g858d7b2824+f45185db35,g88963caddf+0cb8e002cc,g991b906543+f45185db35,g99cad8db69+1747e75aa3,g9b9dfce982+78139cbddb,g9ddcbc5298+9a081db1e4,ga1e77700b3+a912195c07,gae0086650b+585e252eca,gb0e22166c9+60f28cb32d,gb3a676b8dc+b4feba26a1,gb4b16eec92+f82f04eb54,gba4ed39666+c2a2e4ac27,gbb8dafda3b+215b19b0ab,gc120e1dc64+b0284b5341,gc28159a63d+047b288a59,gc3e9b769f7+dcad4ace9a,gcf0d15dbbd+78139cbddb,gdaeeff99f8+f9a426f77a,ge79ae78c31+047b288a59,w.2024.19
LSST Data Management Base Package
Loading...
Searching...
No Matches
isrTaskLSST.py
Go to the documentation of this file.
1__all__ = ["IsrTaskLSST", "IsrTaskLSSTConfig"]
2
3import numpy
4import math
5
6from . import isrFunctions
7from . import isrQa
8from . import linearize
9from .defects import Defects
10
11from contextlib import contextmanager
12from lsst.afw.cameraGeom import NullLinearityType
13import lsst.pex.config as pexConfig
14import lsst.afw.math as afwMath
15import lsst.pipe.base as pipeBase
16import lsst.afw.image as afwImage
17import lsst.pipe.base.connectionTypes as cT
18from lsst.meas.algorithms.detection import SourceDetectionTask
19
20from .overscan import SerialOverscanCorrectionTask, ParallelOverscanCorrectionTask
21from .overscanAmpConfig import OverscanCameraConfig
22from .assembleCcdTask import AssembleCcdTask
23from .deferredCharge import DeferredChargeTask
24from .crosstalk import CrosstalkTask
25from .masking import MaskingTask
26from .isrStatistics import IsrStatisticsTask
27from .isr import maskNans
28
29
30class IsrTaskLSSTConnections(pipeBase.PipelineTaskConnections,
31 dimensions={"instrument", "exposure", "detector"},
32 defaultTemplates={}):
33 ccdExposure = cT.Input(
34 name="raw",
35 doc="Input exposure to process.",
36 storageClass="Exposure",
37 dimensions=["instrument", "exposure", "detector"],
38 )
39 camera = cT.PrerequisiteInput(
40 name="camera",
41 storageClass="Camera",
42 doc="Input camera to construct complete exposures.",
43 dimensions=["instrument"],
44 isCalibration=True,
45 )
46 dnlLUT = cT.PrerequisiteInput(
47 name="dnlLUT",
48 doc="Look-up table for differential non-linearity.",
49 storageClass="IsrCalib",
50 dimensions=["instrument", "exposure", "detector"],
51 isCalibration=True,
52 # TODO DM 36636
53 )
54 bias = cT.PrerequisiteInput(
55 name="bias",
56 doc="Input bias calibration.",
57 storageClass="ExposureF",
58 dimensions=["instrument", "detector"],
59 isCalibration=True,
60 )
61 deferredChargeCalib = cT.PrerequisiteInput(
62 name="cpCtiCalib",
63 doc="Deferred charge/CTI correction dataset.",
64 storageClass="IsrCalib",
65 dimensions=["instrument", "detector"],
66 isCalibration=True,
67 )
68 linearizer = cT.PrerequisiteInput(
69 name='linearizer',
70 storageClass="Linearizer",
71 doc="Linearity correction calibration.",
72 dimensions=["instrument", "detector"],
73 isCalibration=True,
74 )
75 ptc = cT.PrerequisiteInput(
76 name="ptc",
77 doc="Input Photon Transfer Curve dataset",
78 storageClass="PhotonTransferCurveDataset",
79 dimensions=["instrument", "detector"],
80 isCalibration=True,
81 )
82 crosstalk = cT.PrerequisiteInput(
83 name="crosstalk",
84 doc="Input crosstalk object",
85 storageClass="CrosstalkCalib",
86 dimensions=["instrument", "detector"],
87 isCalibration=True,
88 )
89 defects = cT.PrerequisiteInput(
90 name='defects',
91 doc="Input defect tables.",
92 storageClass="Defects",
93 dimensions=["instrument", "detector"],
94 isCalibration=True,
95 )
96 bfKernel = cT.PrerequisiteInput(
97 name='brighterFatterKernel',
98 doc="Complete kernel + gain solutions.",
99 storageClass="BrighterFatterKernel",
100 dimensions=["instrument", "detector"],
101 isCalibration=True,
102 )
103 dark = cT.PrerequisiteInput(
104 name='dark',
105 doc="Input dark calibration.",
106 storageClass="ExposureF",
107 dimensions=["instrument", "detector"],
108 isCalibration=True,
109 )
110 outputExposure = cT.Output(
111 name='postISRCCD',
112 doc="Output ISR processed exposure.",
113 storageClass="Exposure",
114 dimensions=["instrument", "exposure", "detector"],
115 )
116 preInterpExposure = cT.Output(
117 name='preInterpISRCCD',
118 doc="Output ISR processed exposure, with pixels left uninterpolated.",
119 storageClass="ExposureF",
120 dimensions=["instrument", "exposure", "detector"],
121 )
122 outputBin1Exposure = cT.Output(
123 name="postIsrBin1",
124 doc="First binned image.",
125 storageClass="ExposureF",
126 dimensions=["instrument", "exposure", "detector"],
127 )
128 outputBin2Exposure = cT.Output(
129 name="postIsrBin2",
130 doc="Second binned image.",
131 storageClass="ExposureF",
132 dimensions=["instrument", "exposure", "detector"],
133 )
134
135 outputStatistics = cT.Output(
136 name="isrStatistics",
137 doc="Output of additional statistics table.",
138 storageClass="StructuredDataDict",
139 dimensions=["instrument", "exposure", "detector"],
140 )
141
142 def __init__(self, *, config=None):
143 super().__init__(config=config)
144
145 if config.doDiffNonLinearCorrection is not True:
146 del self.dnlLUT
147 if config.doBias is not True:
148 del self.bias
149 if config.doDeferredCharge is not True:
150 del self.deferredChargeCalib
151 if config.doLinearize is not True:
152 del self.linearizer
153 if not config.doCrosstalk and not config.overscanCamera.doAnyParallelOverscanCrosstalk:
154 del self.crosstalk
155 if config.doDefect is not True:
156 del self.defects
157 if config.doBrighterFatter is not True:
158 del self.bfKernel
159 if config.doDark is not True:
160 del self.dark
161
162 if config.doBinnedExposures is not True:
163 del self.outputBin1Exposure
164 del self.outputBin2Exposure
165 if config.doSaveInterpPixels is not True:
166 del self.preInterpExposure
167
168 if config.doCalculateStatistics is not True:
169 del self.outputStatistics
170
171
172class IsrTaskLSSTConfig(pipeBase.PipelineTaskConfig,
173 pipelineConnections=IsrTaskLSSTConnections):
174 """Configuration parameters for IsrTaskLSST.
175
176 Items are grouped in the order in which they are executed by the task.
177 """
178 expectWcs = pexConfig.Field(
179 dtype=bool,
180 default=True,
181 doc="Expect input science images to have a WCS (set False for e.g. spectrographs)."
182 )
183 qa = pexConfig.ConfigField(
184 dtype=isrQa.IsrQaConfig,
185 doc="QA related configuration options.",
186 )
187 doHeaderProvenance = pexConfig.Field(
188 dtype=bool,
189 default=True,
190 doc="Write calibration identifiers into output exposure header.",
191 )
192
193 # Differential non-linearity correction.
194 doDiffNonLinearCorrection = pexConfig.Field(
195 dtype=bool,
196 doc="Do differential non-linearity correction?",
197 default=False,
198 )
199
200 overscanCamera = pexConfig.ConfigField(
201 dtype=OverscanCameraConfig,
202 doc="Per-detector and per-amplifier overscan configurations.",
203 )
204
205 # Amplifier to CCD assembly configuration.
206 doAssembleCcd = pexConfig.Field(
207 dtype=bool,
208 default=True,
209 doc="Assemble amp-level exposures into a ccd-level exposure?"
210 )
211 assembleCcd = pexConfig.ConfigurableField(
212 target=AssembleCcdTask,
213 doc="CCD assembly task.",
214 )
215
216 # Bias subtraction.
217 doBias = pexConfig.Field(
218 dtype=bool,
219 doc="Apply bias frame correction?",
220 default=True,
221 )
222
223 # Deferred charge correction.
224 doDeferredCharge = pexConfig.Field(
225 dtype=bool,
226 doc="Apply deferred charge correction?",
227 default=True,
228 )
229 deferredChargeCorrection = pexConfig.ConfigurableField(
230 target=DeferredChargeTask,
231 doc="Deferred charge correction task.",
232 )
233
234 # Linearization.
235 doLinearize = pexConfig.Field(
236 dtype=bool,
237 doc="Correct for nonlinearity of the detector's response?",
238 default=True,
239 )
240
241 # Gains.
242 doGainsCorrection = pexConfig.Field(
243 dtype=bool,
244 doc="Apply temperature correction to the gains?",
245 default=False,
246 )
247 doApplyGains = pexConfig.Field(
248 dtype=bool,
249 doc="Apply gains to the image?",
250 default=True,
251 )
252
253 # Variance construction.
254 doVariance = pexConfig.Field(
255 dtype=bool,
256 doc="Calculate variance?",
257 default=True
258 )
259 gain = pexConfig.Field(
260 dtype=float,
261 doc="The gain to use if no Detector is present in the Exposure (ignored if NaN).",
262 default=float("NaN"),
263 )
264 maskNegativeVariance = pexConfig.Field(
265 dtype=bool,
266 doc="Mask pixels that claim a negative variance. This likely indicates a failure "
267 "in the measurement of the overscan at an edge due to the data falling off faster "
268 "than the overscan model can account for it.",
269 default=True,
270 )
271 negativeVarianceMaskName = pexConfig.Field(
272 dtype=str,
273 doc="Mask plane to use to mark pixels with negative variance, if `maskNegativeVariance` is True.",
274 default="BAD",
275 )
276 doSaturation = pexConfig.Field(
277 dtype=bool,
278 doc="Mask saturated pixels? NB: this is totally independent of the"
279 " interpolation option - this is ONLY setting the bits in the mask."
280 " To have them interpolated make sure doSaturationInterpolation=True",
281 default=True,
282 )
283 saturation = pexConfig.Field(
284 dtype=float,
285 doc="The saturation level to use if no Detector is present in the Exposure (ignored if NaN)",
286 default=float("NaN"),
287 )
288 saturatedMaskName = pexConfig.Field(
289 dtype=str,
290 doc="Name of mask plane to use in saturation detection and interpolation.",
291 default="SAT",
292 )
293 doSuspect = pexConfig.Field(
294 dtype=bool,
295 doc="Mask suspect pixels?",
296 default=False,
297 )
298 suspectMaskName = pexConfig.Field(
299 dtype=str,
300 doc="Name of mask plane to use for suspect pixels.",
301 default="SUSPECT",
302 )
303
304 # Crosstalk.
305 doCrosstalk = pexConfig.Field(
306 dtype=bool,
307 doc="Apply intra-CCD crosstalk correction?",
308 default=True,
309 )
310 crosstalk = pexConfig.ConfigurableField(
311 target=CrosstalkTask,
312 doc="Intra-CCD crosstalk correction.",
313 )
314
315 # Masking options.
316 doDefect = pexConfig.Field(
317 dtype=bool,
318 doc="Apply correction for CCD defects, e.g. hot pixels?",
319 default=True,
320 )
321 doNanMasking = pexConfig.Field(
322 dtype=bool,
323 doc="Mask non-finite (NAN, inf) pixels.",
324 default=True,
325 )
326 doWidenSaturationTrails = pexConfig.Field(
327 dtype=bool,
328 doc="Widen bleed trails based on their width.",
329 default=True,
330 )
331 masking = pexConfig.ConfigurableField(
332 target=MaskingTask,
333 doc="Masking task."
334 )
335
336 # Interpolation options.
337 doInterpolate = pexConfig.Field(
338 dtype=bool,
339 doc="Interpolate masked pixels?",
340 default=True,
341 )
342 maskListToInterpolate = pexConfig.ListField(
343 dtype=str,
344 doc="List of mask planes that should be interpolated.",
345 default=['SAT', 'BAD'],
346 )
347 doSaveInterpPixels = pexConfig.Field(
348 dtype=bool,
349 doc="Save a copy of the pre-interpolated pixel values?",
350 default=False,
351 )
352
353 # Initial masking options.
354 doSetBadRegions = pexConfig.Field(
355 dtype=bool,
356 doc="Should we set the level of all BAD patches of the chip to the chip's average value?",
357 default=True,
358 )
359
360 # Brighter-Fatter correction.
361 doBrighterFatter = pexConfig.Field(
362 dtype=bool,
363 doc="Apply the brighter-fatter correction?",
364 default=True,
365 )
366 brighterFatterLevel = pexConfig.ChoiceField(
367 dtype=str,
368 doc="The level at which to correct for brighter-fatter.",
369 allowed={
370 "AMP": "Every amplifier treated separately.",
371 "DETECTOR": "One kernel per detector.",
372 },
373 default="DETECTOR",
374 )
375 brighterFatterMaxIter = pexConfig.Field(
376 dtype=int,
377 doc="Maximum number of iterations for the brighter-fatter correction.",
378 default=10,
379 )
380 brighterFatterThreshold = pexConfig.Field(
381 dtype=float,
382 doc="Threshold used to stop iterating the brighter-fatter correction. It is the "
383 "absolute value of the difference between the current corrected image and the one "
384 "from the previous iteration summed over all the pixels.",
385 default=1000,
386 )
387 brighterFatterApplyGain = pexConfig.Field(
388 dtype=bool,
389 doc="Should the gain be applied when applying the brighter-fatter correction?",
390 default=True,
391 )
392 brighterFatterMaskListToInterpolate = pexConfig.ListField(
393 dtype=str,
394 doc="List of mask planes that should be interpolated over when applying the brighter-fatter "
395 "correction.",
396 default=["SAT", "BAD", "NO_DATA", "UNMASKEDNAN"],
397 )
398 brighterFatterMaskGrowSize = pexConfig.Field(
399 dtype=int,
400 doc="Number of pixels to grow the masks listed in config.brighterFatterMaskListToInterpolate "
401 "when brighter-fatter correction is applied.",
402 default=0,
403 )
404 brighterFatterFwhmForInterpolation = pexConfig.Field(
405 dtype=float,
406 doc="FWHM of PSF in arcseconds used for interpolation in brighter-fatter correction "
407 "(currently unused).",
408 default=1.0,
409 )
410 growSaturationFootprintSize = pexConfig.Field(
411 dtype=int,
412 doc="Number of pixels by which to grow the saturation footprints.",
413 default=1,
414 )
415 brighterFatterMaskListToInterpolate = pexConfig.ListField(
416 dtype=str,
417 doc="List of mask planes that should be interpolated over when applying the brighter-fatter."
418 "correction.",
419 default=["SAT", "BAD", "NO_DATA", "UNMASKEDNAN"],
420 )
421
422 # Dark subtraction.
423 doDark = pexConfig.Field(
424 dtype=bool,
425 doc="Apply dark frame correction.",
426 default=True,
427 )
428
429 # Flat correction.
430 doFlat = pexConfig.Field(
431 dtype=bool,
432 doc="Apply flat field correction.",
433 default=True,
434 )
435 flatScalingType = pexConfig.ChoiceField(
436 dtype=str,
437 doc="The method for scaling the flat on the fly.",
438 default='USER',
439 allowed={
440 "USER": "Scale by flatUserScale",
441 "MEAN": "Scale by the inverse of the mean",
442 "MEDIAN": "Scale by the inverse of the median",
443 },
444 )
445 flatUserScale = pexConfig.Field(
446 dtype=float,
447 doc="If flatScalingType is 'USER' then scale flat by this amount; ignored otherwise.",
448 default=1.0,
449 )
450
451 # Calculate image quality statistics?
452 doStandardStatistics = pexConfig.Field(
453 dtype=bool,
454 doc="Should standard image quality statistics be calculated?",
455 default=True,
456 )
457 # Calculate additional statistics?
458 doCalculateStatistics = pexConfig.Field(
459 dtype=bool,
460 doc="Should additional ISR statistics be calculated?",
461 default=True,
462 )
463 isrStats = pexConfig.ConfigurableField(
464 target=IsrStatisticsTask,
465 doc="Task to calculate additional statistics.",
466 )
467
468 # Make binned images?
469 doBinnedExposures = pexConfig.Field(
470 dtype=bool,
471 doc="Should binned exposures be calculated?",
472 default=False,
473 )
474 binFactor1 = pexConfig.Field(
475 dtype=int,
476 doc="Binning factor for first binned exposure. This is intended for a finely binned output.",
477 default=8,
478 check=lambda x: x > 1,
479 )
480 binFactor2 = pexConfig.Field(
481 dtype=int,
482 doc="Binning factor for second binned exposure. This is intended for a coarsely binned output.",
483 default=64,
484 check=lambda x: x > 1,
485 )
486
487 def validate(self):
488 super().validate()
489
490 if self.doCalculateStatistics and self.isrStats.doCtiStatistics:
491 # DM-41912: Implement doApplyGains in LSST IsrTask
492 # if self.doApplyGains !=
493 # self.isrStats.doApplyGainsForCtiStatistics:
494 raise ValueError("doApplyGains must match isrStats.applyGainForCtiStatistics.")
495
496 def setDefaults(self):
497 super().setDefaults()
498
499
500class IsrTaskLSST(pipeBase.PipelineTask):
501 ConfigClass = IsrTaskLSSTConfig
502 _DefaultName = "isr"
503
504 def __init__(self, **kwargs):
505 super().__init__(**kwargs)
506 self.makeSubtask("assembleCcd")
507 self.makeSubtask("deferredChargeCorrection")
508 self.makeSubtask("crosstalk")
509 self.makeSubtask("masking")
510 self.makeSubtask("isrStats")
511
512 def runQuantum(self, butlerQC, inputRefs, outputRefs):
513
514 inputs = butlerQC.get(inputRefs)
515 self.validateInput(inputs)
516 super().runQuantum(butlerQC, inputRefs, outputRefs)
517
518 def validateInput(self, inputs):
519 """
520 This is a check that all the inputs required by the config
521 are available.
522 """
523
524 doCrosstalk = self.config.doCrosstalk or self.config.overscanCamera.doAnyParallelOverscanCrosstalk
525
526 inputMap = {'dnlLUT': self.config.doDiffNonLinearCorrection,
527 'bias': self.config.doBias,
528 'deferredChargeCalib': self.config.doDeferredCharge,
529 'linearizer': self.config.doLinearize,
530 'ptc': self.config.doApplyGains,
531 'crosstalk': doCrosstalk,
532 'defects': self.config.doDefect,
533 'bfKernel': self.config.doBrighterFatter,
534 'dark': self.config.doDark,
535 }
536
537 for calibrationFile, configValue in inputMap.items():
538 if configValue and inputs[calibrationFile] is None:
539 raise RuntimeError("Must supply ", calibrationFile)
540
541 def diffNonLinearCorrection(self, ccdExposure, dnlLUT, **kwargs):
542 # TODO DM 36636
543 # isrFunctions.diffNonLinearCorrection
544 pass
545
546 def maskFullDefectAmplifiers(self, ccdExposure, detector, defects):
547 """
548 Check for fully masked bad amplifiers and mask them.
549
550 Full defect masking happens later to allow for defects which
551 cross amplifier boundaries.
552
553 Parameters
554 ----------
555 ccdExposure : `lsst.afw.image.Exposure`
556 Input exposure to be masked.
557 detector : `lsst.afw.cameraGeom.Detector`
558 Detector object.
559 defects : `lsst.ip.isr.Defects`
560 List of defects. Used to determine if an entire
561 amplifier is bad.
562
563 Returns
564 -------
565 badAmpDict : `str`[`bool`]
566 Dictionary of amplifiers, keyed by name, value is True if
567 amplifier is fully masked.
568 """
569 badAmpDict = {}
570
571 maskedImage = ccdExposure.getMaskedImage()
572
573 for amp in detector:
574 ampName = amp.getName()
575 badAmpDict[ampName] = False
576
577 # Check if entire amp region is defined as a defect
578 # NB: need to use amp.getBBox() for correct comparison with current
579 # defects definition.
580 if defects is not None:
581 badAmpDict[ampName] = bool(sum([v.getBBox().contains(amp.getBBox()) for v in defects]))
582
583 # In the case of a bad amp, we will set mask to "BAD"
584 # (here use amp.getRawBBox() for correct association with pixels in
585 # current ccdExposure).
586 if badAmpDict[ampName]:
587 dataView = afwImage.MaskedImageF(maskedImage, amp.getRawBBox(),
588 afwImage.PARENT)
589 maskView = dataView.getMask()
590 maskView |= maskView.getPlaneBitMask("BAD")
591 del maskView
592
593 self.log.warning("Amplifier %s is bad (completely covered with defects)", ampName)
594
595 return badAmpDict
596
597 def maskSaturatedPixels(self, badAmpDict, ccdExposure, detector):
598 """
599 Mask SATURATED and SUSPECT pixels and check if any amplifiers
600 are fully masked.
601
602 Parameters
603 ----------
604 badAmpDict : `str` [`bool`]
605 Dictionary of amplifiers, keyed by name, value is True if
606 amplifier is fully masked.
607 ccdExposure : `lsst.afw.image.Exposure`
608 Input exposure to be masked.
609 detector : `lsst.afw.cameraGeom.Detector`
610 Detector object.
611 defects : `lsst.ip.isr.Defects`
612 List of defects. Used to determine if an entire
613 amplifier is bad.
614
615 Returns
616 -------
617 badAmpDict : `str`[`bool`]
618 Dictionary of amplifiers, keyed by name.
619 """
620 maskedImage = ccdExposure.getMaskedImage()
621
622 for amp in detector:
623 ampName = amp.getName()
624
625 if badAmpDict[ampName]:
626 # No need to check fully bad amplifiers.
627 continue
628
629 # Mask saturated and suspect pixels.
630 limits = {}
631 if self.config.doSaturation:
632 # Set to the default from the camera model.
633 limits.update({self.config.saturatedMaskName: amp.getSaturation()})
634 # And update if it is set in the config.
635 if math.isfinite(self.config.saturation):
636 limits.update({self.config.saturatedMaskName: self.config.saturation})
637 if self.config.doSuspect:
638 limits.update({self.config.suspectMaskName: amp.getSuspectLevel()})
639
640 for maskName, maskThreshold in limits.items():
641 if not math.isnan(maskThreshold):
642 dataView = maskedImage.Factory(maskedImage, amp.getRawBBox())
643 isrFunctions.makeThresholdMask(
644 maskedImage=dataView,
645 threshold=maskThreshold,
646 growFootprints=0,
647 maskName=maskName
648 )
649
650 # Determine if we've fully masked this amplifier with SUSPECT and
651 # SAT pixels.
652 maskView = afwImage.Mask(maskedImage.getMask(), amp.getRawDataBBox(),
653 afwImage.PARENT)
654 maskVal = maskView.getPlaneBitMask([self.config.saturatedMaskName,
655 self.config.suspectMaskName])
656 if numpy.all(maskView.getArray() & maskVal > 0):
657 self.log.warning("Amplifier %s is bad (completely SATURATED or SUSPECT)", ampName)
658 badAmpDict[ampName] = True
659 maskView |= maskView.getPlaneBitMask("BAD")
660
661 return badAmpDict
662
663 def overscanCorrection(self, mode, detectorConfig, detector, badAmpDict, ccdExposure):
664 """Apply serial overscan correction in place to all amps.
665
666 The actual overscan subtraction is performed by the
667 `lsst.ip.isr.overscan.OverscanTask`, which is called here.
668
669 Parameters
670 ----------
671 mode : `str`
672 Must be `SERIAL` or `PARALLEL`.
673 detectorConfig : `lsst.ip.isr.OverscanDetectorConfig`
674 Per-amplifier configurations.
675 detector : `lsst.afw.cameraGeom.Detector`
676 Detector object.
677 badAmpDict : `dict`
678 Dictionary of amp name to whether it is a bad amp.
679 ccdExposure : `lsst.afw.image.Exposure`
680 Exposure to have overscan correction performed.
681
682 Returns
683 -------
684 overscans : `list` [`lsst.pipe.base.Struct` or None]
685 Each result struct has components:
686
687 ``imageFit``
688 Value or fit subtracted from the amplifier image data.
689 (scalar or `lsst.afw.image.Image`)
690 ``overscanFit``
691 Value or fit subtracted from the overscan image data.
692 (scalar or `lsst.afw.image.Image`)
693 ``overscanImage``
694 Image of the overscan region with the overscan
695 correction applied. This quantity is used to estimate
696 the amplifier read noise empirically.
697 (`lsst.afw.image.Image`)
698 ``overscanMean``
699 Mean overscan fit value. (`float`)
700 ``overscanMedian``
701 Median overscan fit value. (`float`)
702 ``overscanSigma``
703 Clipped standard deviation of the overscan fit. (`float`)
704 ``residualMean``
705 Mean of the overscan after fit subtraction. (`float`)
706 ``residualMedian``
707 Median of the overscan after fit subtraction. (`float`)
708 ``residualSigma``
709 Clipped standard deviation of the overscan after fit
710 subtraction. (`float`)
711
712 See Also
713 --------
714 lsst.ip.isr.overscan.OverscanTask
715 """
716 if mode not in ["SERIAL", "PARALLEL"]:
717 raise ValueError("Mode must be SERIAL or PARALLEL")
718
719 # This returns a list in amp order, with None for uncorrected amps.
720 overscans = []
721
722 for i, amp in enumerate(detector):
723 ampName = amp.getName()
724
725 ampConfig = detectorConfig.getOverscanAmpConfig(amp)
726
727 if mode == "SERIAL" and not ampConfig.doSerialOverscan:
728 self.log.debug(
729 "ISR_OSCAN: Amplifier %s/%s configured to skip serial overscan.",
730 detector.getName(),
731 ampName,
732 )
733 results = None
734 elif mode == "PARALLEL" and not ampConfig.doParallelOverscan:
735 self.log.debug(
736 "ISR_OSCAN: Amplifier %s configured to skip parallel overscan.",
737 detector.getName(),
738 ampName,
739 )
740 results = None
741 elif badAmpDict[ampName] or not ccdExposure.getBBox().contains(amp.getBBox()):
742 results = None
743 else:
744 # This check is to confirm that we are not trying to run
745 # overscan on an already trimmed image. Therefore, always
746 # checking just the horizontal overscan bounding box is
747 # sufficient.
748 if amp.getRawHorizontalOverscanBBox().isEmpty():
749 self.log.warning(
750 "ISR_OSCAN: No overscan region for amp %s. Not performing overscan correction.",
751 ampName,
752 )
753 results = None
754 else:
755 if mode == "SERIAL":
756 # We need to set up the subtask here with a custom
757 # configuration.
758 serialOverscan = SerialOverscanCorrectionTask(config=ampConfig.serialOverscanConfig)
759 results = serialOverscan.run(ccdExposure, amp)
760 else:
761 parallelOverscan = ParallelOverscanCorrectionTask(
762 config=ampConfig.parallelOverscanConfig,
763 )
764 results = parallelOverscan.run(ccdExposure, amp)
765
766 metadata = ccdExposure.getMetadata()
767 keyBase = "LSST ISR OVERSCAN"
768 metadata[f"{keyBase} {mode} MEAN {ampName}"] = results.overscanMean
769 metadata[f"{keyBase} {mode} MEDIAN {ampName}"] = results.overscanMedian
770 metadata[f"{keyBase} {mode} STDEV {ampName}"] = results.overscanSigma
771
772 metadata[f"{keyBase} RESIDUAL {mode} MEAN {ampName}"] = results.residualMean
773 metadata[f"{keyBase} RESIDUAL {mode} MEDIAN {ampName}"] = results.residualMedian
774 metadata[f"{keyBase} RESIDUAL {mode} STDEV {ampName}"] = results.residualSigma
775
776 overscans.append(results)
777
778 # Question: should this be finer grained?
779 ccdExposure.getMetadata().set("OVERSCAN", "Overscan corrected")
780
781 return overscans
782
783 def getLinearizer(self, detector):
784 # Here we assume linearizer as dict or LUT are not supported
785 # TODO DM 28741
786
787 # TODO construct isrcalib input
788 linearizer = linearize.Linearizer(detector=detector, log=self.log)
789 self.log.warning("Constructing linearizer from cameraGeom information.")
790
791 return linearizer
792
793 def gainsCorrection(self, **kwargs):
794 # TODO DM 36639
795 gains = []
796 readNoise = []
797
798 return gains, readNoise
799
800 def updateVariance(self, ampExposure, amp, ptcDataset=None):
801 """Set the variance plane using the gain and read noise.
802
803 Parameters
804 ----------
805 ampExposure : `lsst.afw.image.Exposure`
806 Exposure to process.
807 amp : `lsst.afw.cameraGeom.Amplifier` or `FakeAmp`
808 Amplifier detector data.
809 ptcDataset : `lsst.ip.isr.PhotonTransferCurveDataset`, optional
810 PTC dataset containing the gains and read noise.
811
812 Raises
813 ------
814 RuntimeError
815 Raised if ptcDataset is not provided.
816
817 See also
818 --------
819 lsst.ip.isr.isrFunctions.updateVariance
820 """
821 # Get gains from PTC
822 if ptcDataset is None:
823 raise RuntimeError("No ptcDataset provided to use PTC gains.")
824 else:
825 gain = ptcDataset.gain[amp.getName()]
826 self.log.debug("Getting gain from Photon Transfer Curve.")
827
828 if math.isnan(gain):
829 gain = 1.0
830 self.log.warning("Gain set to NAN! Updating to 1.0 to generate Poisson variance.")
831 elif gain <= 0:
832 patchedGain = 1.0
833 self.log.warning("Gain for amp %s == %g <= 0; setting to %f.",
834 amp.getName(), gain, patchedGain)
835 gain = patchedGain
836
837 # Get read noise from PTC
838 if ptcDataset is None:
839 raise RuntimeError("No ptcDataset provided to use PTC readnoise.")
840 else:
841 readNoise = ptcDataset.noise[amp.getName()]
842 self.log.debug("Getting read noise from Photon Transfer Curve.")
843
844 metadata = ampExposure.getMetadata()
845 metadata[f'LSST GAIN {amp.getName()}'] = gain
846 metadata[f'LSST READNOISE {amp.getName()}'] = readNoise
847
848 isrFunctions.updateVariance(
849 maskedImage=ampExposure.getMaskedImage(),
850 gain=gain,
851 readNoise=readNoise,
852 )
853
854 def maskNegativeVariance(self, exposure):
855 """Identify and mask pixels with negative variance values.
856
857 Parameters
858 ----------
859 exposure : `lsst.afw.image.Exposure`
860 Exposure to process.
861
862 See Also
863 --------
864 lsst.ip.isr.isrFunctions.updateVariance
865 """
866 maskPlane = exposure.getMask().getPlaneBitMask(self.config.negativeVarianceMaskName)
867 bad = numpy.where(exposure.getVariance().getArray() <= 0.0)
868 exposure.mask.array[bad] |= maskPlane
869
870 def variancePlane(self, ccdExposure, ccd, ptc):
871 for amp in ccd:
872 if ccdExposure.getBBox().contains(amp.getBBox()):
873 self.log.debug("Constructing variance map for amplifer %s.", amp.getName())
874 ampExposure = ccdExposure.Factory(ccdExposure, amp.getBBox())
875
876 self.updateVariance(ampExposure, amp, ptcDataset=ptc)
877
878 if self.config.qa is not None and self.config.qa.saveStats is True:
879 qaStats = afwMath.makeStatistics(ampExposure.getVariance(),
880 afwMath.MEDIAN | afwMath.STDEVCLIP)
881 self.log.debug(" Variance stats for amplifer %s: %f +/- %f.",
882 amp.getName(), qaStats.getValue(afwMath.MEDIAN),
883 qaStats.getValue(afwMath.STDEVCLIP))
884 if self.config.maskNegativeVariance:
885 self.maskNegativeVariance(ccdExposure)
886
887 def maskDefect(self, exposure, defectBaseList):
888 """Mask defects using mask plane "BAD", in place.
889
890 Parameters
891 ----------
892 exposure : `lsst.afw.image.Exposure`
893 Exposure to process.
894
895 defectBaseList : defect-type
896 List of defects to mask. Can be of type `lsst.ip.isr.Defects`
897 or `list` of `lsst.afw.image.DefectBase`.
898 """
899 maskedImage = exposure.getMaskedImage()
900 if not isinstance(defectBaseList, Defects):
901 # Promotes DefectBase to Defect
902 defectList = Defects(defectBaseList)
903 else:
904 defectList = defectBaseList
905 defectList.maskPixels(maskedImage, maskName="BAD")
906
907 def maskEdges(self, exposure, numEdgePixels=0, maskPlane="SUSPECT", level='DETECTOR'):
908 """Mask edge pixels with applicable mask plane.
909
910 Parameters
911 ----------
912 exposure : `lsst.afw.image.Exposure`
913 Exposure to process.
914 numEdgePixels : `int`, optional
915 Number of edge pixels to mask.
916 maskPlane : `str`, optional
917 Mask plane name to use.
918 level : `str`, optional
919 Level at which to mask edges.
920 """
921 maskedImage = exposure.getMaskedImage()
922 maskBitMask = maskedImage.getMask().getPlaneBitMask(maskPlane)
923
924 if numEdgePixels > 0:
925 if level == 'DETECTOR':
926 boxes = [maskedImage.getBBox()]
927 elif level == 'AMP':
928 boxes = [amp.getBBox() for amp in exposure.getDetector()]
929
930 for box in boxes:
931 # This makes a bbox numEdgeSuspect pixels smaller than the
932 # image on each side
933 subImage = maskedImage[box]
934 box.grow(-numEdgePixels)
935 # Mask pixels outside box
936 SourceDetectionTask.setEdgeBits(
937 subImage,
938 box,
939 maskBitMask)
940
941 def maskNan(self, exposure):
942 """Mask NaNs using mask plane "UNMASKEDNAN", in place.
943
944 Parameters
945 ----------
946 exposure : `lsst.afw.image.Exposure`
947 Exposure to process.
948
949 Notes
950 -----
951 We mask over all non-finite values (NaN, inf), including those
952 that are masked with other bits (because those may or may not be
953 interpolated over later, and we want to remove all NaN/infs).
954 Despite this behaviour, the "UNMASKEDNAN" mask plane is used to
955 preserve the historical name.
956 """
957 maskedImage = exposure.getMaskedImage()
958
959 # Find and mask NaNs
960 maskedImage.getMask().addMaskPlane("UNMASKEDNAN")
961 maskVal = maskedImage.getMask().getPlaneBitMask("UNMASKEDNAN")
962 numNans = maskNans(maskedImage, maskVal)
963 self.metadata["NUMNANS"] = numNans
964 if numNans > 0:
965 self.log.warning("There were %d unmasked NaNs.", numNans)
966
967 def countBadPixels(self, exposure):
968 """
969 Notes
970 -----
971 Reset and interpolate bad pixels.
972
973 Large contiguous bad regions (which should have the BAD mask
974 bit set) should have their values set to the image median.
975 This group should include defects and bad amplifiers. As the
976 area covered by these defects are large, there's little
977 reason to expect that interpolation would provide a more
978 useful value.
979
980 Smaller defects can be safely interpolated after the larger
981 regions have had their pixel values reset. This ensures
982 that the remaining defects adjacent to bad amplifiers (as an
983 example) do not attempt to interpolate extreme values.
984 """
985 badPixelCount, badPixelValue = isrFunctions.setBadRegions(exposure)
986 if badPixelCount > 0:
987 self.log.info("Set %d BAD pixels to %f.", badPixelCount, badPixelValue)
988
989 @contextmanager
990 def flatContext(self, exp, flat, dark=None):
991 """Context manager that applies and removes flats and darks,
992 if the task is configured to apply them.
993
994 Parameters
995 ----------
996 exp : `lsst.afw.image.Exposure`
997 Exposure to process.
998 flat : `lsst.afw.image.Exposure`
999 Flat exposure the same size as ``exp``.
1000 dark : `lsst.afw.image.Exposure`, optional
1001 Dark exposure the same size as ``exp``.
1002
1003 Yields
1004 ------
1005 exp : `lsst.afw.image.Exposure`
1006 The flat and dark corrected exposure.
1007 """
1008 if self.config.doDark and dark is not None:
1009 self.darkCorrection(exp, dark)
1010 if self.config.doFlat and flat is not None:
1011 self.flatCorrection(exp, flat)
1012 try:
1013 yield exp
1014 finally:
1015 if self.config.doFlat and flat is not None:
1016 self.flatCorrection(exp, flat, invert=True)
1017 if self.config.doDark and dark is not None:
1018 self.darkCorrection(exp, dark, invert=True)
1019
1020 def getBrighterFatterKernel(self, detector, bfKernel):
1021 detName = detector.getName()
1022
1023 # This is expected to be a dictionary of amp-wise gains.
1024 bfGains = bfKernel.gain
1025 if bfKernel.level == 'DETECTOR':
1026 if detName in bfKernel.detKernels:
1027 bfKernelOut = bfKernel.detKernels[detName]
1028 return bfKernelOut, bfGains
1029 else:
1030 raise RuntimeError("Failed to extract kernel from new-style BF kernel.")
1031 elif bfKernel.level == 'AMP':
1032 self.log.warning("Making DETECTOR level kernel from AMP based brighter "
1033 "fatter kernels.")
1034 bfKernel.makeDetectorKernelFromAmpwiseKernels(detName)
1035 bfKernelOut = bfKernel.detKernels[detName]
1036 return bfKernelOut, bfGains
1037
1038 def applyBrighterFatterCorrection(self, ccdExposure, flat, dark, bfKernel, bfGains):
1039 # We need to apply flats and darks before we can interpolate, and
1040 # we need to interpolate before we do B-F, but we do B-F without
1041 # the flats and darks applied so we can work in units of electrons
1042 # or holes. This context manager applies and then removes the darks
1043 # and flats.
1044 #
1045 # We also do not want to interpolate values here, so operate on
1046 # temporary images so we can apply only the BF-correction and roll
1047 # back the interpolation.
1048 # This won't be necessary once the gain normalization
1049 # is done appropriately.
1050 interpExp = ccdExposure.clone()
1051 with self.flatContext(interpExp, flat, dark):
1052 isrFunctions.interpolateFromMask(
1053 maskedImage=interpExp.getMaskedImage(),
1054 fwhm=self.config.brighterFatterFwhmForInterpolation,
1055 growSaturatedFootprints=self.config.growSaturationFootprintSize,
1056 maskNameList=list(self.config.brighterFatterMaskListToInterpolate)
1057 )
1058 bfExp = interpExp.clone()
1059 self.log.info("Applying brighter-fatter correction using kernel type %s / gains %s.",
1060 type(bfKernel), type(bfGains))
1061 bfResults = isrFunctions.brighterFatterCorrection(bfExp, bfKernel,
1062 self.config.brighterFatterMaxIter,
1063 self.config.brighterFatterThreshold,
1064 self.config.brighterFatterApplyGain,
1065 bfGains)
1066 if bfResults[1] == self.config.brighterFatterMaxIter:
1067 self.log.warning("Brighter-fatter correction did not converge, final difference %f.",
1068 bfResults[0])
1069 else:
1070 self.log.info("Finished brighter-fatter correction in %d iterations.",
1071 bfResults[1])
1072
1073 image = ccdExposure.getMaskedImage().getImage()
1074 bfCorr = bfExp.getMaskedImage().getImage()
1075 bfCorr -= interpExp.getMaskedImage().getImage()
1076 image += bfCorr
1077
1078 # Applying the brighter-fatter correction applies a
1079 # convolution to the science image. At the edges this
1080 # convolution may not have sufficient valid pixels to
1081 # produce a valid correction. Mark pixels within the size
1082 # of the brighter-fatter kernel as EDGE to warn of this
1083 # fact.
1084 self.log.info("Ensuring image edges are masked as EDGE to the brighter-fatter kernel size.")
1085 self.maskEdges(ccdExposure, numEdgePixels=numpy.max(bfKernel.shape) // 2,
1086 maskPlane="EDGE")
1087
1088 if self.config.brighterFatterMaskGrowSize > 0:
1089 self.log.info("Growing masks to account for brighter-fatter kernel convolution.")
1090 for maskPlane in self.config.brighterFatterMaskListToInterpolate:
1091 isrFunctions.growMasks(ccdExposure.getMask(),
1092 radius=self.config.brighterFatterMaskGrowSize,
1093 maskNameList=maskPlane,
1094 maskValue=maskPlane)
1095
1096 return ccdExposure
1097
1098 def darkCorrection(self, exposure, darkExposure, invert=False):
1099 """Apply dark correction in place.
1100
1101 Parameters
1102 ----------
1103 exposure : `lsst.afw.image.Exposure`
1104 Exposure to process.
1105 darkExposure : `lsst.afw.image.Exposure`
1106 Dark exposure of the same size as ``exposure``.
1107 invert : `Bool`, optional
1108 If True, re-add the dark to an already corrected image.
1109
1110 Raises
1111 ------
1112 RuntimeError
1113 Raised if either ``exposure`` or ``darkExposure`` do not
1114 have their dark time defined.
1115
1116 See Also
1117 --------
1118 lsst.ip.isr.isrFunctions.darkCorrection
1119 """
1120 expScale = exposure.getInfo().getVisitInfo().getDarkTime()
1121 if math.isnan(expScale):
1122 raise RuntimeError("Exposure darktime is NAN.")
1123 if darkExposure.getInfo().getVisitInfo() is not None \
1124 and not math.isnan(darkExposure.getInfo().getVisitInfo().getDarkTime()):
1125 darkScale = darkExposure.getInfo().getVisitInfo().getDarkTime()
1126 else:
1127 # DM-17444: darkExposure.getInfo.getVisitInfo() is None
1128 # so getDarkTime() does not exist.
1129 self.log.warning("darkExposure.getInfo().getVisitInfo() does not exist. Using darkScale = 1.0.")
1130 darkScale = 1.0
1131
1132 isrFunctions.darkCorrection(
1133 maskedImage=exposure.getMaskedImage(),
1134 darkMaskedImage=darkExposure.getMaskedImage(),
1135 expScale=expScale,
1136 darkScale=darkScale,
1137 invert=invert,
1138 )
1139
1140 @staticmethod
1142 """Extract common calibration metadata values that will be written to
1143 output header.
1144
1145 Parameters
1146 ----------
1147 calib : `lsst.afw.image.Exposure` or `lsst.ip.isr.IsrCalib`
1148 Calibration to pull date information from.
1149
1150 Returns
1151 -------
1152 dateString : `str`
1153 Calibration creation date string to add to header.
1154 """
1155 if hasattr(calib, "getMetadata"):
1156 if 'CALIB_CREATION_DATE' in calib.getMetadata():
1157 return " ".join((calib.getMetadata().get("CALIB_CREATION_DATE", "Unknown"),
1158 calib.getMetadata().get("CALIB_CREATION_TIME", "Unknown")))
1159 else:
1160 return " ".join((calib.getMetadata().get("CALIB_CREATE_DATE", "Unknown"),
1161 calib.getMetadata().get("CALIB_CREATE_TIME", "Unknown")))
1162 else:
1163 return "Unknown Unknown"
1164
1165 def doLinearize(self, detector):
1166 """Check if linearization is needed for the detector cameraGeom.
1167
1168 Checks config.doLinearize and the linearity type of the first
1169 amplifier.
1170
1171 Parameters
1172 ----------
1173 detector : `lsst.afw.cameraGeom.Detector`
1174 Detector to get linearity type from.
1175
1176 Returns
1177 -------
1178 doLinearize : `Bool`
1179 If True, linearization should be performed.
1180 """
1181 return self.config.doLinearize and \
1182 detector.getAmplifiers()[0].getLinearityType() != NullLinearityType
1183
1184 def flatCorrection(self, exposure, flatExposure, invert=False):
1185 """Apply flat correction in place.
1186
1187 Parameters
1188 ----------
1189 exposure : `lsst.afw.image.Exposure`
1190 Exposure to process.
1191 flatExposure : `lsst.afw.image.Exposure`
1192 Flat exposure of the same size as ``exposure``.
1193 invert : `Bool`, optional
1194 If True, unflatten an already flattened image.
1195
1196 See Also
1197 --------
1198 lsst.ip.isr.isrFunctions.flatCorrection
1199 """
1200 isrFunctions.flatCorrection(
1201 maskedImage=exposure.getMaskedImage(),
1202 flatMaskedImage=flatExposure.getMaskedImage(),
1203 scalingType=self.config.flatScalingType,
1204 userScale=self.config.flatUserScale,
1205 invert=invert
1206 )
1207
1208 def makeBinnedImages(self, exposure):
1209 """Make visualizeVisit style binned exposures.
1210
1211 Parameters
1212 ----------
1213 exposure : `lsst.afw.image.Exposure`
1214 Exposure to bin.
1215
1216 Returns
1217 -------
1218 bin1 : `lsst.afw.image.Exposure`
1219 Binned exposure using binFactor1.
1220 bin2 : `lsst.afw.image.Exposure`
1221 Binned exposure using binFactor2.
1222 """
1223 mi = exposure.getMaskedImage()
1224
1225 bin1 = afwMath.binImage(mi, self.config.binFactor1)
1226 bin2 = afwMath.binImage(mi, self.config.binFactor2)
1227
1228 return bin1, bin2
1229
1230 def run(self, ccdExposure, *, dnlLUT=None, bias=None, deferredChargeCalib=None, linearizer=None,
1231 ptc=None, crosstalk=None, defects=None, bfKernel=None, bfGains=None, dark=None,
1232 flat=None, camera=None, **kwargs
1233 ):
1234
1235 detector = ccdExposure.getDetector()
1236
1237 overscanDetectorConfig = self.config.overscanCamera.getOverscanDetectorConfig(detector)
1238
1239 gains = ptc.gain
1240
1241 if self.config.doHeaderProvenance:
1242 # Inputs have been validated, so we can add their date
1243 # information to the output header.
1244 exposureMetadata = ccdExposure.getMetadata()
1245 exposureMetadata["LSST CALIB OVERSCAN HASH"] = overscanDetectorConfig.md5
1246 exposureMetadata["LSST CALIB DATE PTC"] = self.extractCalibDate(ptc)
1247 if self.config.doDiffNonLinearCorrection:
1248 exposureMetadata["LSST CALIB DATE DNL"] = self.extractCalibDate(dnlLUT)
1249 if self.config.doBias:
1250 exposureMetadata["LSST CALIB DATE BIAS"] = self.extractCalibDate(bias)
1251 if self.config.doDeferredCharge:
1252 exposureMetadata["LSST CALIB DATE CTI"] = self.extractCalibDate(deferredChargeCalib)
1253 if self.doLinearize(detector):
1254 exposureMetadata["LSST CALIB DATE LINEARIZER"] = self.extractCalibDate(linearizer)
1255 if self.config.doCrosstalk or overscanDetectorConfig.doAnyParallelOverscanCrosstalk:
1256 exposureMetadata["LSST CALIB DATE CROSSTALK"] = self.extractCalibDate(crosstalk)
1257 if self.config.doDefect:
1258 exposureMetadata["LSST CALIB DATE DEFECTS"] = self.extractCalibDate(defects)
1259 if self.config.doBrighterFatter:
1260 exposureMetadata["LSST CALIB DATE BFK"] = self.extractCalibDate(bfKernel)
1261 if self.config.doDark:
1262 exposureMetadata["LSST CALIB DATE DARK"] = self.extractCalibDate(dark)
1263 if self.config.doFlat:
1264 exposureMetadata["LSST CALIB DATE FLAT"] = self.extractCalibDate(flat)
1265
1266 # First we mark which amplifiers are completely bad from defects.
1267 badAmpDict = self.maskFullDefectAmplifiers(ccdExposure, detector, defects)
1268
1269 if self.config.doDiffNonLinearCorrection:
1270 self.diffNonLinearCorrection(ccdExposure, dnlLUT)
1271
1272 if overscanDetectorConfig.doAnySerialOverscan:
1273 # Input units: ADU
1274 serialOverscans = self.overscanCorrection(
1275 "SERIAL",
1276 overscanDetectorConfig,
1277 detector,
1278 badAmpDict,
1279 ccdExposure,
1280 )
1281 else:
1282 serialOverscans = [None]*len(detector)
1283
1284 if overscanDetectorConfig.doAnyParallelOverscanCrosstalk:
1285 # Input units: ADU
1286 # Make sure that the units here are consistent with later
1287 # application.
1288 self.crosstalk.run(
1289 ccdExposure,
1290 crosstalk=crosstalk,
1291 camera=camera,
1292 parallelOverscanRegion=True,
1293 detectorConfig=overscanDetectorConfig,
1294 )
1295
1296 # After serial overscan correction, we can mask SATURATED and
1297 # SUSPECT pixels. This updates badAmpDict if any amplifier
1298 # is fully saturated after serial overscan correction.
1299 badAmpDict = self.maskSaturatedPixels(badAmpDict, ccdExposure, detector)
1300
1301 if overscanDetectorConfig.doAnyParallelOverscan:
1302 # Input units: ADU
1303 # At the moment we do not use the parallelOverscans return value.
1304 _ = self.overscanCorrection(
1305 "PARALLEL",
1306 overscanDetectorConfig,
1307 detector,
1308 badAmpDict,
1309 ccdExposure,
1310 )
1311
1312 if self.config.doAssembleCcd:
1313 # Input units: ADU
1314 self.log.info("Assembling CCD from amplifiers.")
1315 ccdExposure = self.assembleCcd.assembleCcd(ccdExposure)
1316
1317 if self.config.expectWcs and not ccdExposure.getWcs():
1318 self.log.warning("No WCS found in input exposure.")
1319
1320 if self.config.doLinearize:
1321 # Input units: ADU
1322 self.log.info("Applying linearizer.")
1323 linearizer = self.getLinearizer(detector=detector)
1324 linearizer.applyLinearity(image=ccdExposure.getMaskedImage().getImage(),
1325 detector=detector, log=self.log)
1326
1327 if self.config.doCrosstalk:
1328 # Input units: ADU
1329 self.log.info("Applying crosstalk correction.")
1330 self.crosstalk.run(ccdExposure, crosstalk=crosstalk, isTrimmed=True)
1331
1332 if self.config.doBias:
1333 # Input units: ADU
1334 self.log.info("Applying bias correction.")
1335 isrFunctions.biasCorrection(ccdExposure.getMaskedImage(), bias.getMaskedImage())
1336
1337 if self.config.doGainsCorrection:
1338 # TODO DM 36639
1339 self.log.info("Apply temperature dependence to the gains.")
1340 gains, readNoise = self.gainsCorrection(**kwargs)
1341
1342 if self.config.doApplyGains:
1343 # Input units: ADU
1344 # Output units: electrons
1345 self.log.info("Apply PTC gains (temperature corrected or not) to the image.")
1346 isrFunctions.applyGains(ccdExposure, normalizeGains=False, ptcGains=gains)
1347
1348 if self.config.doDeferredCharge:
1349 # Input units: electrons
1350 self.log.info("Applying deferred charge/CTI correction.")
1351 self.deferredChargeCorrection.run(ccdExposure, deferredChargeCalib)
1352
1353 if self.config.doVariance:
1354 # Input units: electrons
1355 self.variancePlane(ccdExposure, detector, ptc)
1356
1357 # Masking block (defects, NAN pixels and trails).
1358 # Saturated and suspect pixels have already been masked.
1359 if self.config.doDefect:
1360 # Input units: electrons
1361 self.log.info("Applying defects masking.")
1362 self.maskDefect(ccdExposure, defects)
1363
1364 if self.config.doNanMasking:
1365 self.log.info("Masking non-finite (NAN, inf) value pixels.")
1366 self.maskNan(ccdExposure)
1367
1368 if self.config.doWidenSaturationTrails:
1369 self.log.info("Widening saturation trails.")
1370 isrFunctions.widenSaturationTrails(ccdExposure.getMaskedImage().getMask())
1371
1372 if self.config.doDark:
1373 # Input units: electrons
1374 self.log.info("Applying dark subtraction.")
1375 self.darkCorrection(ccdExposure, dark)
1376
1377 if self.config.doBrighterFatter:
1378 # Input units: electrons
1379 self.log.info("Applying Bright-Fatter kernels.")
1380 bfKernelOut, bfGains = self.getBrighterFatterKernel(detector, bfKernel)
1381 ccdExposure = self.applyBrighterFatterCorrection(ccdExposure, flat, dark, bfKernelOut, bfGains)
1382
1383 if self.config.doFlat:
1384 # Input units: electrons
1385 self.log.info("Applying flat correction.")
1386 # Placeholder while the LSST flat procedure is done.
1387 # The flat here would be a background flat.
1388 self.flatCorrection(ccdExposure, flat)
1389
1390 # Pixel values for masked regions are set here
1391 preInterpExp = None
1392 if self.config.doSaveInterpPixels:
1393 preInterpExp = ccdExposure.clone()
1394
1395 if self.config.doSetBadRegions:
1396 self.log.info('Counting pixels in BAD regions.')
1397 self.countBadPixels(ccdExposure)
1398
1399 if self.config.doInterpolate:
1400 self.log.info("Interpolating masked pixels.")
1401 isrFunctions.interpolateFromMask(
1402 maskedImage=ccdExposure.getMaskedImage(),
1403 fwhm=self.config.brighterFatterFwhmForInterpolation,
1404 growSaturatedFootprints=self.config.growSaturationFootprintSize,
1405 maskNameList=list(self.config.maskListToInterpolate)
1406 )
1407
1408 # Calculate standard image quality statistics
1409 if self.config.doStandardStatistics:
1410 metadata = ccdExposure.getMetadata()
1411 for amp in detector:
1412 ampExposure = ccdExposure.Factory(ccdExposure, amp.getBBox())
1413 ampName = amp.getName()
1414 metadata[f"LSST ISR MASK SAT {ampName}"] = isrFunctions.countMaskedPixels(
1415 ampExposure.getMaskedImage(),
1416 [self.config.saturatedMaskName]
1417 )
1418 metadata[f"LSST ISR MASK BAD {ampName}"] = isrFunctions.countMaskedPixels(
1419 ampExposure.getMaskedImage(),
1420 ["BAD"]
1421 )
1422 qaStats = afwMath.makeStatistics(ampExposure.getImage(),
1423 afwMath.MEAN | afwMath.MEDIAN | afwMath.STDEVCLIP)
1424
1425 metadata[f"LSST ISR FINAL MEAN {ampName}"] = qaStats.getValue(afwMath.MEAN)
1426 metadata[f"LSST ISR FINAL MEDIAN {ampName}"] = qaStats.getValue(afwMath.MEDIAN)
1427 metadata[f"LSST ISR FINAL STDEV {ampName}"] = qaStats.getValue(afwMath.STDEVCLIP)
1428
1429 k1 = f"LSST ISR FINAL MEDIAN {ampName}"
1430 k2 = f"LSST ISR OVERSCAN SERIAL MEDIAN {ampName}"
1431 if overscanDetectorConfig.doAnySerialOverscan and k1 in metadata and k2 in metadata:
1432 metadata[f"LSST ISR LEVEL {ampName}"] = metadata[k1] - metadata[k2]
1433 else:
1434 metadata[f"LSST ISR LEVEL {ampName}"] = numpy.nan
1435
1436 # calculate additional statistics.
1437 outputStatistics = None
1438 if self.config.doCalculateStatistics:
1439 outputStatistics = self.isrStats.run(ccdExposure, overscanResults=serialOverscans,
1440 bias=bias, dark=dark, flat=flat, ptc=ptc,
1441 defects=defects).results
1442
1443 # do image binning.
1444 outputBin1Exposure = None
1445 outputBin2Exposure = None
1446 if self.config.doBinnedExposures:
1447 outputBin1Exposure, outputBin2Exposure = self.makeBinnedImages(ccdExposure)
1448
1449 return pipeBase.Struct(
1450 exposure=ccdExposure,
1451
1452 outputBin1Exposure=outputBin1Exposure,
1453 outputBin2Exposure=outputBin2Exposure,
1454
1455 preInterpExposure=preInterpExp,
1456 outputExposure=ccdExposure,
1457 outputStatistics=outputStatistics,
1458 )
Represent a 2-dimensional array of bitmask pixels.
Definition Mask.h:77
runQuantum(self, butlerQC, inputRefs, outputRefs)
variancePlane(self, ccdExposure, ccd, ptc)
flatContext(self, exp, flat, dark=None)
flatCorrection(self, exposure, flatExposure, invert=False)
updateVariance(self, ampExposure, amp, ptcDataset=None)
maskSaturatedPixels(self, badAmpDict, ccdExposure, detector)
getBrighterFatterKernel(self, detector, bfKernel)
darkCorrection(self, exposure, darkExposure, invert=False)
overscanCorrection(self, mode, detectorConfig, detector, badAmpDict, ccdExposure)
diffNonLinearCorrection(self, ccdExposure, dnlLUT, **kwargs)
maskFullDefectAmplifiers(self, ccdExposure, detector, defects)
maskEdges(self, exposure, numEdgePixels=0, maskPlane="SUSPECT", level='DETECTOR')
maskDefect(self, exposure, defectBaseList)
applyBrighterFatterCorrection(self, ccdExposure, flat, dark, bfKernel, bfGains)
run(self, ccdExposure, *dnlLUT=None, bias=None, deferredChargeCalib=None, linearizer=None, ptc=None, crosstalk=None, defects=None, bfKernel=None, bfGains=None, dark=None, flat=None, camera=None, **kwargs)
daf::base::PropertySet * set
Definition fits.cc:931
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
std::shared_ptr< ImageT > binImage(ImageT const &inImage, int const binX, int const binY, lsst::afw::math::Property const flags=lsst::afw::math::MEAN)
Definition binImage.cc:44