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 }])
371 # filter None, because astropy can't deal.
372 inMeta = self.getMetadata().toDict()
373 outMeta = {k: v for k, v in inMeta.items() if v is not None}
374 outMeta.update({k: "" for k, v in inMeta.items() if v is None})
375 catalog.meta = outMeta
376 tableList.append(catalog)
378 if self.interChip:
379 interChipTable = Table([{'IC_SOURCE_DET': sourceDet,
380 'IC_COEFFS': self.interChip[sourceDet].reshape(self.nAmp*self.nAmp)}
381 for sourceDet in self.interChip.keys()])
382 tableList.append(interChipTable)
383 return tableList
385 # Implementation methods.
386 @staticmethod
387 def extractAmp(image, amp, ampTarget, isTrimmed=False):
388 """Extract the image data from an amp, flipped to match ampTarget.
390 Parameters
391 ----------
392 image : `lsst.afw.image.Image` or `lsst.afw.image.MaskedImage`
393 Image containing the amplifier of interest.
394 amp : `lsst.afw.cameraGeom.Amplifier`
395 Amplifier on image to extract.
396 ampTarget : `lsst.afw.cameraGeom.Amplifier`
397 Target amplifier that the extracted image will be flipped
398 to match.
399 isTrimmed : `bool`
400 The image is already trimmed.
401 TODO : DM-15409 will resolve this.
403 Returns
404 -------
405 output : `lsst.afw.image.Image`
406 Image of the amplifier in the desired configuration.
407 """
408 X_FLIP = {lsst.afw.cameraGeom.ReadoutCorner.LL: False,
409 lsst.afw.cameraGeom.ReadoutCorner.LR: True,
410 lsst.afw.cameraGeom.ReadoutCorner.UL: False,
411 lsst.afw.cameraGeom.ReadoutCorner.UR: True}
412 Y_FLIP = {lsst.afw.cameraGeom.ReadoutCorner.LL: False,
413 lsst.afw.cameraGeom.ReadoutCorner.LR: False,
414 lsst.afw.cameraGeom.ReadoutCorner.UL: True,
415 lsst.afw.cameraGeom.ReadoutCorner.UR: True}
417 output = image[amp.getBBox() if isTrimmed else amp.getRawDataBBox()]
418 thisAmpCorner = amp.getReadoutCorner()
419 targetAmpCorner = ampTarget.getReadoutCorner()
421 # Flipping is necessary only if the desired configuration doesn't match what we currently have
422 xFlip = X_FLIP[targetAmpCorner] ^ X_FLIP[thisAmpCorner]
423 yFlip = Y_FLIP[targetAmpCorner] ^ Y_FLIP[thisAmpCorner]
424 return lsst.afw.math.flipImage(output, xFlip, yFlip)
426 @staticmethod
427 def calculateBackground(mi, badPixels=["BAD"]):
428 """Estimate median background in image.
430 Getting a great background model isn't important for crosstalk correction,
431 since the crosstalk is at a low level. The median should be sufficient.
433 Parameters
434 ----------
435 mi : `lsst.afw.image.MaskedImage`
436 MaskedImage for which to measure background.
437 badPixels : `list` of `str`
438 Mask planes to ignore.
439 Returns
440 -------
441 bg : `float`
442 Median background level.
443 """
444 mask = mi.getMask()
445 stats = lsst.afw.math.StatisticsControl()
446 stats.setAndMask(mask.getPlaneBitMask(badPixels))
447 return lsst.afw.math.makeStatistics(mi, lsst.afw.math.MEDIAN, stats).getValue()
449 def subtractCrosstalk(self, thisExposure, sourceExposure=None, crosstalkCoeffs=None,
450 badPixels=["BAD"], minPixelToMask=45000,
451 crosstalkStr="CROSSTALK", isTrimmed=False,
452 backgroundMethod="None"):
453 """Subtract the crosstalk from thisExposure, optionally using a different source.
455 We set the mask plane indicated by ``crosstalkStr`` in a target amplifier
456 for pixels in a source amplifier that exceed ``minPixelToMask``. Note that
457 the correction is applied to all pixels in the amplifier, but only those
458 that have a substantial crosstalk are masked with ``crosstalkStr``.
460 The uncorrected image is used as a template for correction. This is good
461 enough if the crosstalk is small (e.g., coefficients < ~ 1e-3), but if it's
462 larger you may want to iterate.
464 Parameters
465 ----------
466 thisExposure : `lsst.afw.image.Exposure`
467 Exposure for which to subtract crosstalk.
468 sourceExposure : `lsst.afw.image.Exposure`, optional
469 Exposure to use as the source of the crosstalk. If not set,
470 thisExposure is used as the source (intra-detector crosstalk).
471 crosstalkCoeffs : `numpy.ndarray`, optional.
472 Coefficients to use to correct crosstalk.
473 badPixels : `list` of `str`
474 Mask planes to ignore.
475 minPixelToMask : `float`
476 Minimum pixel value (relative to the background level) in
477 source amplifier for which to set ``crosstalkStr`` mask plane
478 in target amplifier.
479 crosstalkStr : `str`
480 Mask plane name for pixels greatly modified by crosstalk
481 (above minPixelToMask).
482 isTrimmed : `bool`
483 The image is already trimmed.
484 This should no longer be needed once DM-15409 is resolved.
485 backgroundMethod : `str`
486 Method used to subtract the background. "AMP" uses
487 amplifier-by-amplifier background levels, "DETECTOR" uses full
488 exposure/maskedImage levels. Any other value results in no
489 background subtraction.
490 """
491 mi = thisExposure.getMaskedImage()
492 mask = mi.getMask()
493 detector = thisExposure.getDetector()
494 if self.hasCrosstalk is False:
495 self.fromDetector(detector, coeffVector=crosstalkCoeffs)
497 numAmps = len(detector)
498 if numAmps != self.nAmp:
499 raise RuntimeError(f"Crosstalk built for {self.nAmp} in {self._detectorName}, received "
500 f"{numAmps} in {detector.getName()}")
502 if sourceExposure:
503 source = sourceExposure.getMaskedImage()
504 sourceDetector = sourceExposure.getDetector()
505 else:
506 source = mi
507 sourceDetector = detector
509 if crosstalkCoeffs is not None:
510 coeffs = crosstalkCoeffs
511 else:
512 coeffs = self.coeffs
513 self.log.debug("CT COEFF: %s", coeffs)
514 # Set background level based on the requested method. The
515 # thresholdBackground holds the offset needed so that we only mask
516 # pixels high relative to the background, not in an absolute
517 # sense.
518 thresholdBackground = self.calculateBackground(source, badPixels)
520 backgrounds = [0.0 for amp in sourceDetector]
521 if backgroundMethod is None:
522 pass
523 elif backgroundMethod == "AMP":
524 backgrounds = [self.calculateBackground(source[amp.getBBox()], badPixels)
525 for amp in sourceDetector]
526 elif backgroundMethod == "DETECTOR":
527 backgrounds = [self.calculateBackground(source, badPixels) for amp in sourceDetector]
529 # Set the crosstalkStr bit for the bright pixels (those which will have
530 # significant crosstalk correction)
531 crosstalkPlane = mask.addMaskPlane(crosstalkStr)
532 footprints = lsst.afw.detection.FootprintSet(source,
533 lsst.afw.detection.Threshold(minPixelToMask
534 + thresholdBackground))
535 footprints.setMask(mask, crosstalkStr)
536 crosstalk = mask.getPlaneBitMask(crosstalkStr)
538 # Define a subtrahend image to contain all the scaled crosstalk signals
539 subtrahend = source.Factory(source.getBBox())
540 subtrahend.set((0, 0, 0))
542 coeffs = coeffs.transpose()
543 for ii, iAmp in enumerate(sourceDetector):
544 iImage = subtrahend[iAmp.getBBox() if isTrimmed else iAmp.getRawDataBBox()]
545 for jj, jAmp in enumerate(detector):
546 if coeffs[ii, jj] == 0.0:
547 continue
548 jImage = self.extractAmp(mi, jAmp, iAmp, isTrimmed)
549 jImage.getMask().getArray()[:] &= crosstalk # Remove all other masks
550 jImage -= backgrounds[jj]
551 iImage.scaledPlus(coeffs[ii, jj], jImage)
553 # Set crosstalkStr bit only for those pixels that have been significantly modified (i.e., those
554 # masked as such in 'subtrahend'), not necessarily those that are bright originally.
555 mask.clearMaskPlane(crosstalkPlane)
556 mi -= subtrahend # also sets crosstalkStr bit for bright pixels
559class CrosstalkConfig(Config):
560 """Configuration for intra-detector crosstalk removal."""
561 minPixelToMask = Field(
562 dtype=float,
563 doc="Set crosstalk mask plane for pixels over this value.",
564 default=45000
565 )
566 crosstalkMaskPlane = Field(
567 dtype=str,
568 doc="Name for crosstalk mask plane.",
569 default="CROSSTALK"
570 )
571 crosstalkBackgroundMethod = ChoiceField(
572 dtype=str,
573 doc="Type of background subtraction to use when applying correction.",
574 default="None",
575 allowed={
576 "None": "Do no background subtraction.",
577 "AMP": "Subtract amplifier-by-amplifier background levels.",
578 "DETECTOR": "Subtract detector level background."
579 },
580 )
581 useConfigCoefficients = Field(
582 dtype=bool,
583 doc="Ignore the detector crosstalk information in favor of CrosstalkConfig values?",
584 default=False,
585 )
586 crosstalkValues = ListField(
587 dtype=float,
588 doc=("Amplifier-indexed crosstalk coefficients to use. This should be arranged as a 1 x nAmp**2 "
589 "list of coefficients, such that when reshaped by crosstalkShape, the result is nAmp x nAmp. "
590 "This matrix should be structured so CT * [amp0 amp1 amp2 ...]^T returns the column "
591 "vector [corr0 corr1 corr2 ...]^T."),
592 default=[0.0],
593 )
594 crosstalkShape = ListField(
595 dtype=int,
596 doc="Shape of the coefficient array. This should be equal to [nAmp, nAmp].",
597 default=[1],
598 )
600 def getCrosstalk(self, detector=None):
601 """Return a 2-D numpy array of crosstalk coefficients in the proper shape.
603 Parameters
604 ----------
605 detector : `lsst.afw.cameraGeom.detector`
606 Detector that is to be crosstalk corrected.
608 Returns
609 -------
610 coeffs : `numpy.ndarray`
611 Crosstalk coefficients that can be used to correct the detector.
613 Raises
614 ------
615 RuntimeError
616 Raised if no coefficients could be generated from this detector/configuration.
617 """
618 if self.useConfigCoefficients is True:
619 coeffs = np.array(self.crosstalkValues).reshape(self.crosstalkShape)
620 if detector is not None:
621 nAmp = len(detector)
622 if coeffs.shape != (nAmp, nAmp):
623 raise RuntimeError("Constructed crosstalk coeffients do not match detector shape. "
624 f"{coeffs.shape} {nAmp}")
625 return coeffs
626 elif detector is not None and detector.hasCrosstalk() is True:
627 # Assume the detector defines itself consistently.
628 return detector.getCrosstalk()
629 else:
630 raise RuntimeError("Attempted to correct crosstalk without crosstalk coefficients")
632 def hasCrosstalk(self, detector=None):
633 """Return a boolean indicating if crosstalk coefficients exist.
635 Parameters
636 ----------
637 detector : `lsst.afw.cameraGeom.detector`
638 Detector that is to be crosstalk corrected.
640 Returns
641 -------
642 hasCrosstalk : `bool`
643 True if this detector/configuration has crosstalk coefficients defined.
644 """
645 if self.useConfigCoefficients is True and self.crosstalkValues is not None:
646 return True
647 elif detector is not None and detector.hasCrosstalk() is True:
648 return True
649 else:
650 return False
653class CrosstalkTask(Task):
654 """Apply intra-detector crosstalk correction."""
655 ConfigClass = CrosstalkConfig
656 _DefaultName = 'isrCrosstalk'
658 def prepCrosstalk(self, dataRef, crosstalk=None):
659 """Placeholder for crosstalk preparation method, e.g., for inter-detector crosstalk.
661 Parameters
662 ----------
663 dataRef : `daf.persistence.butlerSubset.ButlerDataRef`
664 Butler reference of the detector data to be processed.
665 crosstalk : `~lsst.ip.isr.CrosstalkConfig`
666 Crosstalk calibration that will be used.
668 See also
669 --------
670 lsst.obs.decam.crosstalk.DecamCrosstalkTask.prepCrosstalk
671 """
672 return
674 def run(self, exposure, crosstalk=None,
675 crosstalkSources=None, isTrimmed=False):
676 """Apply intra-detector crosstalk correction
678 Parameters
679 ----------
680 exposure : `lsst.afw.image.Exposure`
681 Exposure for which to remove crosstalk.
682 crosstalkCalib : `lsst.ip.isr.CrosstalkCalib`, optional
683 External crosstalk calibration to apply. Constructed from
684 detector if not found.
685 crosstalkSources : `defaultdict`, optional
686 Image data for other detectors that are sources of
687 crosstalk in exposure. The keys are expected to be names
688 of the other detectors, with the values containing
689 `lsst.afw.image.Exposure` at the same level of processing
690 as ``exposure``.
691 The default for intra-detector crosstalk here is None.
692 isTrimmed : `bool`
693 The image is already trimmed.
694 This should no longer be needed once DM-15409 is resolved.
696 Raises
697 ------
698 RuntimeError
699 Raised if called for a detector that does not have a
700 crosstalk correction.
701 """
702 if not crosstalk:
703 crosstalk = CrosstalkCalib(log=self.log)
704 crosstalk = crosstalk.fromDetector(exposure.getDetector(),
705 coeffVector=self.config.crosstalkValues)
706 if not crosstalk.log:
707 crosstalk.log = self.log
708 if not crosstalk.hasCrosstalk:
709 raise RuntimeError("Attempted to correct crosstalk without crosstalk coefficients.")
711 else:
712 self.log.info("Applying crosstalk correction.")
713 crosstalk.subtractCrosstalk(exposure, crosstalkCoeffs=crosstalk.coeffs,
714 minPixelToMask=self.config.minPixelToMask,
715 crosstalkStr=self.config.crosstalkMaskPlane, isTrimmed=isTrimmed,
716 backgroundMethod=self.config.crosstalkBackgroundMethod)
718 if crosstalk.interChip:
719 if crosstalkSources:
720 for detName in crosstalk.interChip:
721 if isinstance(crosstalkSources[0], 'lsst.afw.image.Exposure'):
722 # Received afwImage.Exposure
723 sourceNames = [exp.getDetector().getName() for exp in crosstalkSources]
724 else:
725 # Received dafButler.DeferredDatasetHandle
726 sourceNames = [expRef.get(datasetType='isrOscanCorr').getDetector().getName()
727 for expRef in crosstalkSources]
728 if detName not in sourceNames:
729 self.log.warn("Crosstalk lists %s, not found in sources: %s",
730 detName, sourceNames)
731 continue
732 interChipCoeffs = crosstalk.interChip[detName]
733 sourceExposure = crosstalkSources[sourceNames.index(detName)]
734 crosstalk.subtractCrosstalk(exposure, sourceExposure=sourceExposure,
735 crosstalkCoeffs=interChipCoeffs,
736 minPixelToMask=self.config.minPixelToMask,
737 crosstalkStr=self.config.crosstalkMaskPlane,
738 isTrimmed=isTrimmed,
739 backgroundMethod=self.config.crosstalkBackgroundMethod)
740 else:
741 self.log.warn("Crosstalk contains interChip coefficients, but no sources found!")
744class NullCrosstalkTask(CrosstalkTask):
745 def run(self, exposure, crosstalkSources=None):
746 self.log.info("Not performing any crosstalk correction")