Coverage for python/lsst/ip/isr/crosstalk.py : 10%

Hot-keys on this page
r m x p toggle line displays
j k next/prev highlighted chunk
0 (zero) top of page
1 (one) first highlighted chunk
1#
2# LSST Data Management System
3# Copyright 2008-2017 AURA/LSST.
4#
5# This product includes software developed by the
6# LSST Project (http://www.lsst.org/).
7#
8# This program is free software: you can redistribute it and/or modify
9# it under the terms of the GNU General Public License as published by
10# the Free Software Foundation, either version 3 of the License, or
11# (at your option) any later version.
12#
13# This program is distributed in the hope that it will be useful,
14# but WITHOUT ANY WARRANTY; without even the implied warranty of
15# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16# GNU General Public License for more details.
17#
18# You should have received a copy of the LSST License Statement and
19# the GNU General Public License along with this program. If not,
20# see <https://www.lsstcorp.org/LegalNotices/>.
21#
22"""
23Apply intra-detector crosstalk corrections
24"""
25import numpy as np
26from astropy.table import Table
28import lsst.afw.math
29import lsst.afw.detection
30from lsst.pex.config import Config, Field, ChoiceField, ListField
31from lsst.pipe.base import Task
33from lsst.ip.isr import IsrCalib
36__all__ = ["CrosstalkCalib", "CrosstalkConfig", "CrosstalkTask",
37 "NullCrosstalkTask"]
40class CrosstalkCalib(IsrCalib):
41 """Calibration of amp-to-amp crosstalk coefficients.
43 Parameters
44 ----------
45 detector : `lsst.afw.cameraGeom.Detector`, optional
46 Detector to use to pull coefficients from.
47 nAmp : `int`, optional
48 Number of amplifiers to initialize.
49 log : `lsst.log.Log`, optional
50 Log to write messages to.
51 **kwargs :
52 Parameters to pass to parent constructor.
54 Notes
55 -----
56 The crosstalk attributes stored are:
58 hasCrosstalk : `bool`
59 Whether there is crosstalk defined for this detector.
60 nAmp : `int`
61 Number of amplifiers in this detector.
62 crosstalkShape : `tuple` [`int`, `int`]
63 A tuple containing the shape of the ``coeffs`` matrix. This
64 should be equivalent to (``nAmp``, ``nAmp``).
65 coeffs : `np.ndarray`
66 A matrix containing the crosstalk coefficients. coeff[i][j]
67 contains the coefficients to calculate the contribution
68 amplifier_j has on amplifier_i (each row[i] contains the
69 corrections for detector_i).
70 coeffErr : `np.ndarray`, optional
71 A matrix (as defined by ``coeffs``) containing the standard
72 distribution of the crosstalk measurements.
73 coeffNum : `np.ndarray`, optional
74 A matrix containing the number of pixel pairs used to measure
75 the ``coeffs`` and ``coeffErr``.
76 coeffValid : `np.ndarray`, optional
77 A matrix of Boolean values indicating if the coefficient is
78 valid, defined as abs(coeff) > coeffErr / sqrt(coeffNum).
79 interChip : `dict` [`np.ndarray`]
80 A dictionary keyed by detectorName containing ``coeffs``
81 matrices used to correct for inter-chip crosstalk with a
82 source on the detector indicated.
84 """
85 _OBSTYPE = 'CROSSTALK'
86 _SCHEMA = 'Gen3 Crosstalk'
87 _VERSION = 1.0
89 def __init__(self, detector=None, nAmp=0, **kwargs):
90 self.hasCrosstalk = False
91 self.nAmp = nAmp if nAmp else 0
92 self.crosstalkShape = (self.nAmp, self.nAmp)
94 self.coeffs = np.zeros(self.crosstalkShape) if self.nAmp else None
95 self.coeffErr = np.zeros(self.crosstalkShape) if self.nAmp else None
96 self.coeffNum = np.zeros(self.crosstalkShape,
97 dtype=int) if self.nAmp else None
98 self.coeffValid = np.zeros(self.crosstalkShape,
99 dtype=bool) if self.nAmp else None
100 self.interChip = {}
102 super().__init__(**kwargs)
103 self.requiredAttributes.update(['hasCrosstalk', 'nAmp', 'coeffs',
104 'coeffErr', 'coeffNum', 'coeffValid',
105 'interChip'])
106 if detector:
107 self.fromDetector(detector)
109 def updateMetadata(self, setDate=False, **kwargs):
110 """Update calibration metadata.
112 This calls the base class's method after ensuring the required
113 calibration keywords will be saved.
115 Parameters
116 ----------
117 setDate : `bool`, optional
118 Update the CALIBDATE fields in the metadata to the current
119 time. Defaults to False.
120 kwargs :
121 Other keyword parameters to set in the metadata.
122 """
123 kwargs['DETECTOR'] = self._detectorId
124 kwargs['DETECTOR_NAME'] = self._detectorName
125 kwargs['DETECTOR_SERIAL'] = self._detectorSerial
126 kwargs['HAS_CROSSTALK'] = self.hasCrosstalk
127 kwargs['NAMP'] = self.nAmp
128 self.crosstalkShape = (self.nAmp, self.nAmp)
129 kwargs['CROSSTALK_SHAPE'] = self.crosstalkShape
131 super().updateMetadata(setDate=setDate, **kwargs)
133 def fromDetector(self, detector, coeffVector=None):
134 """Set calibration parameters from the detector.
136 Parameters
137 ----------
138 detector : `lsst.afw.cameraGeom.Detector`
139 Detector to use to set parameters from.
140 coeffVector : `numpy.array`, optional
141 Use the detector geometry (bounding boxes and flip
142 information), but use ``coeffVector`` instead of the
143 output of ``detector.getCrosstalk()``.
145 Returns
146 -------
147 calib : `lsst.ip.isr.CrosstalkCalib`
148 The calibration constructed from the detector.
150 """
151 if detector.hasCrosstalk() or coeffVector:
152 self._detectorName = detector.getName()
153 self._detectorSerial = detector.getSerial()
155 self.nAmp = len(detector)
156 self.crosstalkShape = (self.nAmp, self.nAmp)
158 if coeffVector is not None:
159 crosstalkCoeffs = coeffVector
160 else:
161 crosstalkCoeffs = detector.getCrosstalk()
162 if len(crosstalkCoeffs) == 1 and crosstalkCoeffs[0] == 0.0:
163 return self
164 self.coeffs = np.array(crosstalkCoeffs).reshape(self.crosstalkShape)
166 if self.coeffs.shape != self.crosstalkShape:
167 raise RuntimeError("Crosstalk coefficients do not match detector shape. "
168 f"{self.crosstalkShape} {self.nAmp}")
170 self.interChip = {}
171 self.hasCrosstalk = True
172 self.updateMetadata()
173 return self
175 @classmethod
176 def fromDict(cls, dictionary):
177 """Construct a calibration from a dictionary of properties.
179 Must be implemented by the specific calibration subclasses.
181 Parameters
182 ----------
183 dictionary : `dict`
184 Dictionary of properties.
186 Returns
187 -------
188 calib : `lsst.ip.isr.CalibType`
189 Constructed calibration.
191 Raises
192 ------
193 RuntimeError :
194 Raised if the supplied dictionary is for a different
195 calibration.
196 """
197 calib = cls()
199 if calib._OBSTYPE != dictionary['metadata']['OBSTYPE']:
200 raise RuntimeError(f"Incorrect crosstalk supplied. Expected {calib._OBSTYPE}, "
201 f"found {dictionary['metadata']['OBSTYPE']}")
203 calib.setMetadata(dictionary['metadata'])
205 if 'detectorName' in dictionary:
206 calib._detectorName = dictionary.get('detectorName')
207 elif 'DETECTOR_NAME' in dictionary:
208 calib._detectorName = dictionary.get('DETECTOR_NAME')
209 elif 'DET_NAME' in dictionary['metadata']:
210 calib._detectorName = dictionary['metadata']['DET_NAME']
211 else:
212 calib._detectorName = None
214 if 'detectorSerial' in dictionary:
215 calib._detectorSerial = dictionary.get('detectorSerial')
216 elif 'DETECTOR_SERIAL' in dictionary:
217 calib._detectorSerial = dictionary.get('DETECTOR_SERIAL')
218 elif 'DET_SER' in dictionary['metadata']:
219 calib._detectorSerial = dictionary['metadata']['DET_SER']
220 else:
221 calib._detectorSerial = None
223 if 'detectorId' in dictionary:
224 calib._detectorId = dictionary.get('detectorId')
225 elif 'DETECTOR' in dictionary:
226 calib._detectorId = dictionary.get('DETECTOR')
227 elif 'DETECTOR' in dictionary['metadata']:
228 calib._detectorId = dictionary['metadata']['DETECTOR']
229 elif calib._detectorSerial:
230 calib._detectorId = calib._detectorSerial
231 else:
232 calib._detectorId = None
234 if 'instrument' in dictionary:
235 calib._instrument = dictionary.get('instrument')
236 elif 'INSTRUME' in dictionary['metadata']:
237 calib._instrument = dictionary['metadata']['INSTRUME']
238 else:
239 calib._instrument = None
241 calib.hasCrosstalk = dictionary.get('hasCrosstalk',
242 dictionary['metadata'].get('HAS_CROSSTALK', False))
243 if calib.hasCrosstalk:
244 calib.nAmp = dictionary.get('nAmp', dictionary['metadata'].get('NAMP', 0))
245 calib.crosstalkShape = (calib.nAmp, calib.nAmp)
246 calib.coeffs = np.array(dictionary['coeffs']).reshape(calib.crosstalkShape)
247 if 'coeffErr' in dictionary:
248 calib.coeffErr = np.array(dictionary['coeffErr']).reshape(calib.crosstalkShape)
249 else:
250 calib.coeffErr = np.zeros_like(calib.coeffs)
251 if 'coeffNum' in dictionary:
252 calib.coeffNum = np.array(dictionary['coeffNum']).reshape(calib.crosstalkShape)
253 else:
254 calib.coeffNum = np.zeros_like(calib.coeffs, dtype=int)
255 if 'coeffValid' in dictionary:
256 calib.coeffValid = np.array(dictionary['coeffValid']).reshape(calib.crosstalkShape)
257 else:
258 calib.coeffValid = np.ones_like(calib.coeffs, dtype=bool)
260 calib.interChip = dictionary.get('interChip', None)
261 if calib.interChip:
262 for detector in calib.interChip:
263 coeffVector = calib.interChip[detector]
264 calib.interChip[detector] = np.array(coeffVector).reshape(calib.crosstalkShape)
266 calib.updateMetadata()
267 return calib
269 def toDict(self):
270 """Return a dictionary containing the calibration properties.
272 The dictionary should be able to be round-tripped through
273 `fromDict`.
275 Returns
276 -------
277 dictionary : `dict`
278 Dictionary of properties.
279 """
280 self.updateMetadata()
282 outDict = {}
283 metadata = self.getMetadata()
284 outDict['metadata'] = metadata
286 outDict['hasCrosstalk'] = self.hasCrosstalk
287 outDict['nAmp'] = self.nAmp
288 outDict['crosstalkShape'] = self.crosstalkShape
290 ctLength = self.nAmp*self.nAmp
291 outDict['coeffs'] = self.coeffs.reshape(ctLength).tolist()
293 if self.coeffErr is not None:
294 outDict['coeffErr'] = self.coeffErr.reshape(ctLength).tolist()
295 if self.coeffNum is not None:
296 outDict['coeffNum'] = self.coeffNum.reshape(ctLength).tolist()
297 if self.coeffValid is not None:
298 outDict['coeffValid'] = self.coeffValid.reshape(ctLength).tolist()
300 if self.interChip:
301 outDict['interChip'] = dict()
302 for detector in self.interChip:
303 outDict['interChip'][detector] = self.interChip[detector].reshape(ctLength).tolist()
305 return outDict
307 @classmethod
308 def fromTable(cls, tableList):
309 """Construct calibration from a list of tables.
311 This method uses the `fromDict` method to create the
312 calibration, after constructing an appropriate dictionary from
313 the input tables.
315 Parameters
316 ----------
317 tableList : `list` [`lsst.afw.table.Table`]
318 List of tables to use to construct the crosstalk
319 calibration.
321 Returns
322 -------
323 calib : `lsst.ip.isr.CrosstalkCalib`
324 The calibration defined in the tables.
326 """
327 coeffTable = tableList[0]
329 metadata = coeffTable.meta
330 inDict = dict()
331 inDict['metadata'] = metadata
332 inDict['hasCrosstalk'] = metadata['HAS_CROSSTALK']
333 inDict['nAmp'] = metadata['NAMP']
335 inDict['coeffs'] = coeffTable['CT_COEFFS']
336 if 'CT_ERRORS' in coeffTable:
337 inDict['coeffErr'] = coeffTable['CT_ERRORS']
338 if 'CT_COUNTS' in coeffTable:
339 inDict['coeffNum'] = coeffTable['CT_COUNTS']
340 if 'CT_VALID' in coeffTable:
341 inDict['coeffValid'] = coeffTable['CT_VALID']
343 if len(tableList) > 1:
344 inDict['interChip'] = dict()
345 interChipTable = tableList[1]
346 for record in interChipTable:
347 inDict['interChip'][record['IC_SOURCE_DET']] = record['IC_COEFFS']
349 return cls().fromDict(inDict)
351 def toTable(self):
352 """Construct a list of tables containing the information in this calibration.
354 The list of tables should create an identical calibration
355 after being passed to this class's fromTable method.
357 Returns
358 -------
359 tableList : `list` [`lsst.afw.table.Table`]
360 List of tables containing the crosstalk calibration
361 information.
363 """
364 tableList = []
365 self.updateMetadata()
366 catalog = Table([{'CT_COEFFS': self.coeffs.reshape(self.nAmp*self.nAmp),
367 'CT_ERRORS': self.coeffErr.reshape(self.nAmp*self.nAmp),
368 'CT_COUNTS': self.coeffNum.reshape(self.nAmp*self.nAmp),
369 'CT_VALID': self.coeffValid.reshape(self.nAmp*self.nAmp),
370 }])
372 catalog.meta = self.getMetadata().toDict()
373 tableList.append(catalog)
375 if self.interChip:
376 interChipTable = Table([{'IC_SOURCE_DET': sourceDet,
377 'IC_COEFFS': self.interChip[sourceDet].reshape(self.nAmp*self.nAmp)}
378 for sourceDet in self.interChip.keys()])
379 tableList.append(interChipTable)
380 return tableList
382 # Implementation methods.
383 @staticmethod
384 def extractAmp(image, amp, ampTarget, isTrimmed=False):
385 """Extract the image data from an amp, flipped to match ampTarget.
387 Parameters
388 ----------
389 image : `lsst.afw.image.Image` or `lsst.afw.image.MaskedImage`
390 Image containing the amplifier of interest.
391 amp : `lsst.afw.cameraGeom.Amplifier`
392 Amplifier on image to extract.
393 ampTarget : `lsst.afw.cameraGeom.Amplifier`
394 Target amplifier that the extracted image will be flipped
395 to match.
396 isTrimmed : `bool`
397 The image is already trimmed.
398 TODO : DM-15409 will resolve this.
400 Returns
401 -------
402 output : `lsst.afw.image.Image`
403 Image of the amplifier in the desired configuration.
404 """
405 X_FLIP = {lsst.afw.cameraGeom.ReadoutCorner.LL: False,
406 lsst.afw.cameraGeom.ReadoutCorner.LR: True,
407 lsst.afw.cameraGeom.ReadoutCorner.UL: False,
408 lsst.afw.cameraGeom.ReadoutCorner.UR: True}
409 Y_FLIP = {lsst.afw.cameraGeom.ReadoutCorner.LL: False,
410 lsst.afw.cameraGeom.ReadoutCorner.LR: False,
411 lsst.afw.cameraGeom.ReadoutCorner.UL: True,
412 lsst.afw.cameraGeom.ReadoutCorner.UR: True}
414 output = image[amp.getBBox() if isTrimmed else amp.getRawDataBBox()]
415 thisAmpCorner = amp.getReadoutCorner()
416 targetAmpCorner = ampTarget.getReadoutCorner()
418 # Flipping is necessary only if the desired configuration doesn't match what we currently have
419 xFlip = X_FLIP[targetAmpCorner] ^ X_FLIP[thisAmpCorner]
420 yFlip = Y_FLIP[targetAmpCorner] ^ Y_FLIP[thisAmpCorner]
421 return lsst.afw.math.flipImage(output, xFlip, yFlip)
423 @staticmethod
424 def calculateBackground(mi, badPixels=["BAD"]):
425 """Estimate median background in image.
427 Getting a great background model isn't important for crosstalk correction,
428 since the crosstalk is at a low level. The median should be sufficient.
430 Parameters
431 ----------
432 mi : `lsst.afw.image.MaskedImage`
433 MaskedImage for which to measure background.
434 badPixels : `list` of `str`
435 Mask planes to ignore.
436 Returns
437 -------
438 bg : `float`
439 Median background level.
440 """
441 mask = mi.getMask()
442 stats = lsst.afw.math.StatisticsControl()
443 stats.setAndMask(mask.getPlaneBitMask(badPixels))
444 return lsst.afw.math.makeStatistics(mi, lsst.afw.math.MEDIAN, stats).getValue()
446 def subtractCrosstalk(self, thisExposure, sourceExposure=None, crosstalkCoeffs=None,
447 badPixels=["BAD"], minPixelToMask=45000,
448 crosstalkStr="CROSSTALK", isTrimmed=False,
449 backgroundMethod="None"):
450 """Subtract the crosstalk from thisExposure, optionally using a different source.
452 We set the mask plane indicated by ``crosstalkStr`` in a target amplifier
453 for pixels in a source amplifier that exceed ``minPixelToMask``. Note that
454 the correction is applied to all pixels in the amplifier, but only those
455 that have a substantial crosstalk are masked with ``crosstalkStr``.
457 The uncorrected image is used as a template for correction. This is good
458 enough if the crosstalk is small (e.g., coefficients < ~ 1e-3), but if it's
459 larger you may want to iterate.
461 Parameters
462 ----------
463 thisExposure : `lsst.afw.image.Exposure`
464 Exposure for which to subtract crosstalk.
465 sourceExposure : `lsst.afw.image.Exposure`, optional
466 Exposure to use as the source of the crosstalk. If not set,
467 thisExposure is used as the source (intra-detector crosstalk).
468 crosstalkCoeffs : `numpy.ndarray`, optional.
469 Coefficients to use to correct crosstalk.
470 badPixels : `list` of `str`
471 Mask planes to ignore.
472 minPixelToMask : `float`
473 Minimum pixel value (relative to the background level) in
474 source amplifier for which to set ``crosstalkStr`` mask plane
475 in target amplifier.
476 crosstalkStr : `str`
477 Mask plane name for pixels greatly modified by crosstalk
478 (above minPixelToMask).
479 isTrimmed : `bool`
480 The image is already trimmed.
481 This should no longer be needed once DM-15409 is resolved.
482 backgroundMethod : `str`
483 Method used to subtract the background. "AMP" uses
484 amplifier-by-amplifier background levels, "DETECTOR" uses full
485 exposure/maskedImage levels. Any other value results in no
486 background subtraction.
487 """
488 mi = thisExposure.getMaskedImage()
489 mask = mi.getMask()
490 detector = thisExposure.getDetector()
491 if self.hasCrosstalk is False:
492 self.fromDetector(detector, coeffVector=crosstalkCoeffs)
494 numAmps = len(detector)
495 if numAmps != self.nAmp:
496 raise RuntimeError(f"Crosstalk built for {self.nAmp} in {self._detectorName}, received "
497 f"{numAmps} in {detector.getName()}")
499 if sourceExposure:
500 source = sourceExposure.getMaskedImage()
501 sourceDetector = sourceExposure.getDetector()
502 else:
503 source = mi
504 sourceDetector = detector
506 if crosstalkCoeffs is not None:
507 coeffs = crosstalkCoeffs
508 else:
509 coeffs = self.coeffs
510 self.log.debug("CT COEFF: %s", coeffs)
511 # Set background level based on the requested method. The
512 # thresholdBackground holds the offset needed so that we only mask
513 # pixels high relative to the background, not in an absolute
514 # sense.
515 thresholdBackground = self.calculateBackground(source, badPixels)
517 backgrounds = [0.0 for amp in sourceDetector]
518 if backgroundMethod is None:
519 pass
520 elif backgroundMethod == "AMP":
521 backgrounds = [self.calculateBackground(source[amp.getBBox()], badPixels)
522 for amp in sourceDetector]
523 elif backgroundMethod == "DETECTOR":
524 backgrounds = [self.calculateBackground(source, badPixels) for amp in sourceDetector]
526 # Set the crosstalkStr bit for the bright pixels (those which will have
527 # significant crosstalk correction)
528 crosstalkPlane = mask.addMaskPlane(crosstalkStr)
529 footprints = lsst.afw.detection.FootprintSet(source,
530 lsst.afw.detection.Threshold(minPixelToMask
531 + thresholdBackground))
532 footprints.setMask(mask, crosstalkStr)
533 crosstalk = mask.getPlaneBitMask(crosstalkStr)
535 # Define a subtrahend image to contain all the scaled crosstalk signals
536 subtrahend = source.Factory(source.getBBox())
537 subtrahend.set((0, 0, 0))
539 coeffs = coeffs.transpose()
540 for ii, iAmp in enumerate(sourceDetector):
541 iImage = subtrahend[iAmp.getBBox() if isTrimmed else iAmp.getRawDataBBox()]
542 for jj, jAmp in enumerate(detector):
543 if coeffs[ii, jj] == 0.0:
544 continue
545 jImage = self.extractAmp(mi, jAmp, iAmp, isTrimmed)
546 jImage.getMask().getArray()[:] &= crosstalk # Remove all other masks
547 jImage -= backgrounds[jj]
548 iImage.scaledPlus(coeffs[ii, jj], jImage)
550 # Set crosstalkStr bit only for those pixels that have been significantly modified (i.e., those
551 # masked as such in 'subtrahend'), not necessarily those that are bright originally.
552 mask.clearMaskPlane(crosstalkPlane)
553 mi -= subtrahend # also sets crosstalkStr bit for bright pixels
556class CrosstalkConfig(Config):
557 """Configuration for intra-detector crosstalk removal."""
558 minPixelToMask = Field(
559 dtype=float,
560 doc="Set crosstalk mask plane for pixels over this value.",
561 default=45000
562 )
563 crosstalkMaskPlane = Field(
564 dtype=str,
565 doc="Name for crosstalk mask plane.",
566 default="CROSSTALK"
567 )
568 crosstalkBackgroundMethod = ChoiceField(
569 dtype=str,
570 doc="Type of background subtraction to use when applying correction.",
571 default="None",
572 allowed={
573 "None": "Do no background subtraction.",
574 "AMP": "Subtract amplifier-by-amplifier background levels.",
575 "DETECTOR": "Subtract detector level background."
576 },
577 )
578 useConfigCoefficients = Field(
579 dtype=bool,
580 doc="Ignore the detector crosstalk information in favor of CrosstalkConfig values?",
581 default=False,
582 )
583 crosstalkValues = ListField(
584 dtype=float,
585 doc=("Amplifier-indexed crosstalk coefficients to use. This should be arranged as a 1 x nAmp**2 "
586 "list of coefficients, such that when reshaped by crosstalkShape, the result is nAmp x nAmp. "
587 "This matrix should be structured so CT * [amp0 amp1 amp2 ...]^T returns the column "
588 "vector [corr0 corr1 corr2 ...]^T."),
589 default=[0.0],
590 )
591 crosstalkShape = ListField(
592 dtype=int,
593 doc="Shape of the coefficient array. This should be equal to [nAmp, nAmp].",
594 default=[1],
595 )
597 def getCrosstalk(self, detector=None):
598 """Return a 2-D numpy array of crosstalk coefficients in the proper shape.
600 Parameters
601 ----------
602 detector : `lsst.afw.cameraGeom.detector`
603 Detector that is to be crosstalk corrected.
605 Returns
606 -------
607 coeffs : `numpy.ndarray`
608 Crosstalk coefficients that can be used to correct the detector.
610 Raises
611 ------
612 RuntimeError
613 Raised if no coefficients could be generated from this detector/configuration.
614 """
615 if self.useConfigCoefficients is True:
616 coeffs = np.array(self.crosstalkValues).reshape(self.crosstalkShape)
617 if detector is not None:
618 nAmp = len(detector)
619 if coeffs.shape != (nAmp, nAmp):
620 raise RuntimeError("Constructed crosstalk coeffients do not match detector shape. "
621 f"{coeffs.shape} {nAmp}")
622 return coeffs
623 elif detector is not None and detector.hasCrosstalk() is True:
624 # Assume the detector defines itself consistently.
625 return detector.getCrosstalk()
626 else:
627 raise RuntimeError("Attempted to correct crosstalk without crosstalk coefficients")
629 def hasCrosstalk(self, detector=None):
630 """Return a boolean indicating if crosstalk coefficients exist.
632 Parameters
633 ----------
634 detector : `lsst.afw.cameraGeom.detector`
635 Detector that is to be crosstalk corrected.
637 Returns
638 -------
639 hasCrosstalk : `bool`
640 True if this detector/configuration has crosstalk coefficients defined.
641 """
642 if self.useConfigCoefficients is True and self.crosstalkValues is not None:
643 return True
644 elif detector is not None and detector.hasCrosstalk() is True:
645 return True
646 else:
647 return False
650class CrosstalkTask(Task):
651 """Apply intra-detector crosstalk correction."""
652 ConfigClass = CrosstalkConfig
653 _DefaultName = 'isrCrosstalk'
655 def prepCrosstalk(self, dataRef, crosstalk=None):
656 """Placeholder for crosstalk preparation method, e.g., for inter-detector crosstalk.
658 Parameters
659 ----------
660 dataRef : `daf.persistence.butlerSubset.ButlerDataRef`
661 Butler reference of the detector data to be processed.
662 crosstalk : `~lsst.ip.isr.CrosstalkConfig`
663 Crosstalk calibration that will be used.
665 See also
666 --------
667 lsst.obs.decam.crosstalk.DecamCrosstalkTask.prepCrosstalk
668 """
669 return
671 def run(self, exposure, crosstalk=None,
672 crosstalkSources=None, isTrimmed=False):
673 """Apply intra-detector crosstalk correction
675 Parameters
676 ----------
677 exposure : `lsst.afw.image.Exposure`
678 Exposure for which to remove crosstalk.
679 crosstalkCalib : `lsst.ip.isr.CrosstalkCalib`, optional
680 External crosstalk calibration to apply. Constructed from
681 detector if not found.
682 crosstalkSources : `defaultdict`, optional
683 Image data for other detectors that are sources of
684 crosstalk in exposure. The keys are expected to be names
685 of the other detectors, with the values containing
686 `lsst.afw.image.Exposure` at the same level of processing
687 as ``exposure``.
688 The default for intra-detector crosstalk here is None.
689 isTrimmed : `bool`
690 The image is already trimmed.
691 This should no longer be needed once DM-15409 is resolved.
693 Raises
694 ------
695 RuntimeError
696 Raised if called for a detector that does not have a
697 crosstalk correction.
698 """
699 if not crosstalk:
700 crosstalk = CrosstalkCalib(log=self.log)
701 crosstalk = crosstalk.fromDetector(exposure.getDetector(),
702 coeffVector=self.config.crosstalkValues)
703 if not crosstalk.log:
704 crosstalk.log = self.log
705 if not crosstalk.hasCrosstalk:
706 raise RuntimeError("Attempted to correct crosstalk without crosstalk coefficients.")
708 else:
709 self.log.info("Applying crosstalk correction.")
710 crosstalk.subtractCrosstalk(exposure, crosstalkCoeffs=crosstalk.coeffs,
711 minPixelToMask=self.config.minPixelToMask,
712 crosstalkStr=self.config.crosstalkMaskPlane, isTrimmed=isTrimmed,
713 backgroundMethod=self.config.crosstalkBackgroundMethod)
715 if crosstalk.interChip:
716 if crosstalkSources:
717 for detName in crosstalk.interChip:
718 if isinstance(crosstalkSources[0], 'lsst.afw.image.Exposure'):
719 # Received afwImage.Exposure
720 sourceNames = [exp.getDetector().getName() for exp in crosstalkSources]
721 else:
722 # Received dafButler.DeferredDatasetHandle
723 sourceNames = [expRef.get(datasetType='isrOscanCorr').getDetector().getName()
724 for expRef in crosstalkSources]
725 if detName not in sourceNames:
726 self.log.warn("Crosstalk lists %s, not found in sources: %s",
727 detName, sourceNames)
728 continue
729 interChipCoeffs = crosstalk.interChip[detName]
730 sourceExposure = crosstalkSources[sourceNames.index(detName)]
731 crosstalk.subtractCrosstalk(exposure, sourceExposure=sourceExposure,
732 crosstalkCoeffs=interChipCoeffs,
733 minPixelToMask=self.config.minPixelToMask,
734 crosstalkStr=self.config.crosstalkMaskPlane,
735 isTrimmed=isTrimmed,
736 backgroundMethod=self.config.crosstalkBackgroundMethod)
737 else:
738 self.log.warn("Crosstalk contains interChip coefficients, but no sources found!")
741class NullCrosstalkTask(CrosstalkTask):
742 def run(self, exposure, crosstalkSources=None):
743 self.log.info("Not performing any crosstalk correction")