LSST Applications g07dc498a13+7851b72aa9,g1409bbee79+7851b72aa9,g1a7e361dbc+7851b72aa9,g1fd858c14a+a4e18a0dda,g33399d78f5+a0324bbf49,g35bb328faa+e55fef2c71,g3bd4b5ce2c+8524b1c0c8,g53246c7159+e55fef2c71,g579b87e3d2+a58ba40925,g60b5630c4e+7b4465799a,g78460c75b0+8427c4cc8f,g78619a8342+5517f7db9e,g786e29fd12+307f82e6af,g8534526c7b+8e1c6b434f,g89139ef638+7851b72aa9,g8b49a6ea8e+7b4465799a,g8ffcb69f3d+0065d7bbc8,g9125e01d80+e55fef2c71,g97b8272a79+a8c4cb337e,g989de1cb63+7851b72aa9,g9f33ca652e+747bd1f1f9,gaaedd4e678+7851b72aa9,gabe3b4be73+9c0c3c7524,gb1101e3267+c03a154bbb,gb58c049af0+28045f66fd,gc1fe0db326+7b4465799a,gca43fec769+e55fef2c71,gce7788e931+99adca4f64,gcf25f946ba+a0324bbf49,gd397e13551+18f805d5e0,gd6cbbdb0b4+f6e5445f66,gde0f65d7ad+78b6ec8427,ge278dab8ac+b4c2c8faf7,geab183fbe5+7b4465799a,gecb8035dfe+1f480bec5e,gf58bf46354+e55fef2c71,gf92a8ffd38+e7bc33f3ea,gfe7187db8c+38a2c5c626,w.2025.03
LSST Data Management Base Package
Loading...
Searching...
No Matches
lsst.ip.isr.isrStatistics.IsrStatisticsTask Class Reference
Inheritance diagram for lsst.ip.isr.isrStatistics.IsrStatisticsTask:

Public Member Functions

 __init__ (self, statControl=None, **kwargs)
 
 run (self, inputExp, ptc=None, overscanResults=None, **kwargs)
 
 measureCti (self, inputExp, overscans, gains)
 
 measureBanding (self, inputExp, overscans)
 
 measureProjectionStatistics (self, inputExp, overscans)
 
 copyCalibDistributionStatistics (self, inputExp, **kwargs)
 
 measureBiasShifts (self, inputExp, overscanResults)
 
 measureAmpCorrelations (self, inputExp, overscanResults)
 
 measureDivisaderoStatistics (self, inputExp, **kwargs)
 

Static Public Member Functions

 makeKernel (kernelSize)
 

Public Attributes

 statControl
 
 statType
 

Static Public Attributes

 ConfigClass = IsrStatisticsTaskConfig
 

Protected Member Functions

 _scan_for_shifts (self, overscanData)
 
 _satisfies_flatness (self, shiftRow, shiftPeak, overscanData)
 

Static Protected Attributes

str _DefaultName = "isrStatistics"
 

Detailed Description

Task to measure arbitrary statistics on ISR processed exposures.

The goal is to wrap a number of optional measurements that are
useful for calibration production and detector stability.

Definition at line 213 of file isrStatistics.py.

Constructor & Destructor Documentation

◆ __init__()

lsst.ip.isr.isrStatistics.IsrStatisticsTask.__init__ ( self,
statControl = None,
** kwargs )

Definition at line 222 of file isrStatistics.py.

222 def __init__(self, statControl=None, **kwargs):
223 super().__init__(**kwargs)
224 self.statControl = afwMath.StatisticsControl(self.config.nSigmaClip, self.config.nIter,
225 afwImage.Mask.getPlaneBitMask(self.config.badMask))
226 self.statType = afwMath.stringToStatisticsProperty(self.config.stat)
227
Pass parameters to a Statistics object.
Definition Statistics.h:83
Property stringToStatisticsProperty(std::string const property)
Conversion function to switch a string to a Property (see Statistics.h)

Member Function Documentation

◆ _satisfies_flatness()

lsst.ip.isr.isrStatistics.IsrStatisticsTask._satisfies_flatness ( self,
shiftRow,
shiftPeak,
overscanData )
protected
Determine if a region is flat.

Parameters
----------
shiftRow : `int`
    Row with possible peak.
shiftPeak : `float`
    Value at the possible peak.
overscanData : `list` [`float`]
    Overscan data used to fit around the possible peak.

Returns
-------
isFlat : `bool`
    Indicates if the region is flat, and so the peak is valid.

Definition at line 737 of file isrStatistics.py.

737 def _satisfies_flatness(self, shiftRow, shiftPeak, overscanData):
738 """Determine if a region is flat.
739
740 Parameters
741 ----------
742 shiftRow : `int`
743 Row with possible peak.
744 shiftPeak : `float`
745 Value at the possible peak.
746 overscanData : `list` [`float`]
747 Overscan data used to fit around the possible peak.
748
749 Returns
750 -------
751 isFlat : `bool`
752 Indicates if the region is flat, and so the peak is valid.
753 """
754 prerange = np.arange(shiftRow - self.config.biasShiftWindow, shiftRow)
755 postrange = np.arange(shiftRow, shiftRow + self.config.biasShiftWindow)
756
757 preFit = linregress(prerange, overscanData[prerange])
758 postFit = linregress(postrange, overscanData[postrange])
759
760 if shiftPeak > 0:
761 preTrend = (2*preFit[0]*len(prerange) < shiftPeak)
762 postTrend = (2*postFit[0]*len(postrange) < shiftPeak)
763 else:
764 preTrend = (2*preFit[0]*len(prerange) > shiftPeak)
765 postTrend = (2*postFit[0]*len(postrange) > shiftPeak)
766
767 return (preTrend and postTrend)
768

◆ _scan_for_shifts()

lsst.ip.isr.isrStatistics.IsrStatisticsTask._scan_for_shifts ( self,
overscanData )
protected
Scan overscan data for shifts.

Parameters
----------
overscanData : `list` [`float`]
     Overscan data to search for shifts.

Returns
-------
noise : `float`
    Noise estimated from Butterworth filtered overscan data.
peaks : `list` [`float`, `float`, `int`, `int`]
    Shift peak information, containing the convolved peak
    value, the raw peak value, and the lower and upper bounds
    of the region checked.

Definition at line 693 of file isrStatistics.py.

693 def _scan_for_shifts(self, overscanData):
694 """Scan overscan data for shifts.
695
696 Parameters
697 ----------
698 overscanData : `list` [`float`]
699 Overscan data to search for shifts.
700
701 Returns
702 -------
703 noise : `float`
704 Noise estimated from Butterworth filtered overscan data.
705 peaks : `list` [`float`, `float`, `int`, `int`]
706 Shift peak information, containing the convolved peak
707 value, the raw peak value, and the lower and upper bounds
708 of the region checked.
709 """
710 numerator, denominator = butter(self.config.biasShiftFilterOrder,
711 self.config.biasShiftCutoff,
712 btype="high", analog=False)
713 noise = np.std(filtfilt(numerator, denominator, overscanData))
714 kernel = np.concatenate([np.arange(self.config.biasShiftWindow),
715 np.arange(-self.config.biasShiftWindow + 1, 0)])
716 kernel = kernel/np.sum(kernel[:self.config.biasShiftWindow])
717
718 convolved = np.convolve(overscanData, kernel, mode="valid")
719 convolved = np.pad(convolved, (self.config.biasShiftWindow - 1, self.config.biasShiftWindow))
720
721 shift_check = np.abs(convolved)/noise
722 shift_mask = shift_check > self.config.biasShiftThreshold
723 shift_mask[:self.config.biasShiftRowSkip] = False
724
725 shift_regions = np.flatnonzero(np.diff(np.r_[np.int8(0),
726 shift_mask.view(np.int8),
727 np.int8(0)])).reshape(-1, 2)
728 shift_peaks = []
729 for region in shift_regions:
730 region_peak = np.argmax(shift_check[region[0]:region[1]]) + region[0]
731 if self._satisfies_flatness(region_peak, convolved[region_peak], overscanData):
732 shift_peaks.append(
733 [float(convolved[region_peak]), float(region_peak),
734 int(region[0]), int(region[1])])
735 return noise, shift_peaks
736

◆ copyCalibDistributionStatistics()

lsst.ip.isr.isrStatistics.IsrStatisticsTask.copyCalibDistributionStatistics ( self,
inputExp,
** kwargs )
Copy calibration statistics for this exposure.

Parameters
----------
inputExp : `lsst.afw.image.Exposure`
    The exposure being processed.
**kwargs :
    Keyword arguments with calibrations.

Returns
-------
outputStats : `dict` [`str`, [`dict` [`str`, `float`]]]
    Dictionary of measurements, keyed by amplifier name and
    statistics segment.

Definition at line 595 of file isrStatistics.py.

595 def copyCalibDistributionStatistics(self, inputExp, **kwargs):
596 """Copy calibration statistics for this exposure.
597
598 Parameters
599 ----------
600 inputExp : `lsst.afw.image.Exposure`
601 The exposure being processed.
602 **kwargs :
603 Keyword arguments with calibrations.
604
605 Returns
606 -------
607 outputStats : `dict` [`str`, [`dict` [`str`, `float`]]]
608 Dictionary of measurements, keyed by amplifier name and
609 statistics segment.
610 """
611 outputStats = {}
612
613 # Amp level elements
614 for amp in inputExp.getDetector():
615 ampStats = {}
616
617 for calibType in ("bias", "dark", "flat"):
618 if kwargs.get(calibType, None) is not None:
619 metadata = kwargs[calibType].getMetadata()
620 for pct in self.config.expectedDistributionLevels:
621 key = f"LSST CALIB {calibType.upper()} {amp.getName()} DISTRIBUTION {pct}-PCT"
622 ampStats[key] = metadata.get(key, np.nan)
623
624 for calibType in ("defects"):
625 if kwargs.get(calibType, None) is not None:
626 metadata = kwargs[calibType].getMetadata()
627 for key in (f"LSST CALIB {calibType.upper()} {amp.getName()} N_HOT",
628 f"LSST CALIB {calibType.upper()} {amp.getName()} N_COLD"):
629 ampStats[key] = metadata.get(key, np.nan)
630 outputStats[amp.getName()] = ampStats
631
632 # Detector level elements
633 for calibType in ("defects"):
634 if kwargs.get(calibType, None) is not None:
635 metadata = kwargs[calibType].getMetadata()
636 for key in (f"LSST CALIB {calibType.upper()} N_BAD_COLUMNS"):
637 outputStats["detector"][key] = metadata.get(key, np.nan)
638
639 return outputStats
640

◆ makeKernel()

lsst.ip.isr.isrStatistics.IsrStatisticsTask.makeKernel ( kernelSize)
static
Make a boxcar smoothing kernel.

Parameters
----------
kernelSize : `int`
    Size of the kernel in pixels.

Returns
-------
kernel : `np.array`
    Kernel for boxcar smoothing.

Definition at line 441 of file isrStatistics.py.

441 def makeKernel(kernelSize):
442 """Make a boxcar smoothing kernel.
443
444 Parameters
445 ----------
446 kernelSize : `int`
447 Size of the kernel in pixels.
448
449 Returns
450 -------
451 kernel : `np.array`
452 Kernel for boxcar smoothing.
453 """
454 if kernelSize > 0:
455 kernel = np.full(kernelSize, 1.0 / kernelSize)
456 else:
457 kernel = np.array([1.0])
458 return kernel
459

◆ measureAmpCorrelations()

lsst.ip.isr.isrStatistics.IsrStatisticsTask.measureAmpCorrelations ( self,
inputExp,
overscanResults )
Measure correlations between amplifier segments.

Parameters
----------
inputExp : `lsst.afw.image.Exposure`
    Exposure to measure.
overscans : `list` [`lsst.pipe.base.Struct`]
    List of overscan results.  Expected fields are:

    ``imageFit``
        Value or fit subtracted from the amplifier image data
        (scalar or `lsst.afw.image.Image`).
    ``overscanFit``
        Value or fit subtracted from the overscan image data
        (scalar or `lsst.afw.image.Image`).
    ``overscanImage``
        Image of the overscan region with the overscan
        correction applied (`lsst.afw.image.Image`). This
        quantity is used to estimate the amplifier read noise
        empirically.

Returns
-------
outputStats : `dict` [`str`, [`dict` [`str`, `float`]]]
    Dictionary of measurements, keyed by amplifier name and
    statistics segment.

Notes
-----
Based on eo_pipe implementation:
https://github.com/lsst-camera-dh/eo_pipe/blob/main/python/lsst/eo/pipe/raft_level_correlations.py  # noqa: E501 W505

Definition at line 769 of file isrStatistics.py.

769 def measureAmpCorrelations(self, inputExp, overscanResults):
770 """Measure correlations between amplifier segments.
771
772 Parameters
773 ----------
774 inputExp : `lsst.afw.image.Exposure`
775 Exposure to measure.
776 overscans : `list` [`lsst.pipe.base.Struct`]
777 List of overscan results. Expected fields are:
778
779 ``imageFit``
780 Value or fit subtracted from the amplifier image data
781 (scalar or `lsst.afw.image.Image`).
782 ``overscanFit``
783 Value or fit subtracted from the overscan image data
784 (scalar or `lsst.afw.image.Image`).
785 ``overscanImage``
786 Image of the overscan region with the overscan
787 correction applied (`lsst.afw.image.Image`). This
788 quantity is used to estimate the amplifier read noise
789 empirically.
790
791 Returns
792 -------
793 outputStats : `dict` [`str`, [`dict` [`str`, `float`]]]
794 Dictionary of measurements, keyed by amplifier name and
795 statistics segment.
796
797 Notes
798 -----
799 Based on eo_pipe implementation:
800 https://github.com/lsst-camera-dh/eo_pipe/blob/main/python/lsst/eo/pipe/raft_level_correlations.py # noqa: E501 W505
801 """
802 detector = inputExp.getDetector()
803 rows = []
804
805 for ampId, overscan in enumerate(overscanResults):
806 rawOverscan = overscan.overscanImage.image.array + overscan.overscanFit
807 rawOverscan = rawOverscan.ravel()
808
809 ampImage = inputExp[detector[ampId].getBBox()]
810 ampImage = ampImage.image.array.ravel()
811
812 for ampId2, overscan2 in enumerate(overscanResults):
813 osC = 1.0
814 imC = 1.0
815 if ampId2 != ampId:
816 rawOverscan2 = overscan2.overscanImage.image.array + overscan2.overscanFit
817 rawOverscan2 = rawOverscan2.ravel()
818
819 osC = np.corrcoef(rawOverscan, rawOverscan2)[0, 1]
820
821 ampImage2 = inputExp[detector[ampId2].getBBox()]
822 ampImage2 = ampImage2.image.array.ravel()
823
824 imC = np.corrcoef(ampImage, ampImage2)[0, 1]
825 rows.append(
826 {'detector': detector.getId(),
827 'detectorComp': detector.getId(),
828 'ampName': detector[ampId].getName(),
829 'ampComp': detector[ampId2].getName(),
830 'imageCorr': float(imC),
831 'overscanCorr': float(osC),
832 }
833 )
834 return rows
835

◆ measureBanding()

lsst.ip.isr.isrStatistics.IsrStatisticsTask.measureBanding ( self,
inputExp,
overscans )
Task to measure banding statistics.

Parameters
----------
inputExp : `lsst.afw.image.Exposure`
    Exposure to measure.
overscans : `list` [`lsst.pipe.base.Struct`]
    List of overscan results.  Expected fields are:

    ``imageFit``
        Value or fit subtracted from the amplifier image data
        (scalar or `lsst.afw.image.Image`).
    ``overscanFit``
        Value or fit subtracted from the overscan image data
        (scalar or `lsst.afw.image.Image`).
    ``overscanImage``
        Image of the overscan region with the overscan
        correction applied (`lsst.afw.image.Image`). This
        quantity is used to estimate the amplifier read noise
        empirically.

Returns
-------
outputStats : `dict` [`str`, [`dict` [`str`, `float`]]]
    Dictionary of measurements, keyed by amplifier name and
    statistics segment.

Definition at line 460 of file isrStatistics.py.

460 def measureBanding(self, inputExp, overscans):
461 """Task to measure banding statistics.
462
463 Parameters
464 ----------
465 inputExp : `lsst.afw.image.Exposure`
466 Exposure to measure.
467 overscans : `list` [`lsst.pipe.base.Struct`]
468 List of overscan results. Expected fields are:
469
470 ``imageFit``
471 Value or fit subtracted from the amplifier image data
472 (scalar or `lsst.afw.image.Image`).
473 ``overscanFit``
474 Value or fit subtracted from the overscan image data
475 (scalar or `lsst.afw.image.Image`).
476 ``overscanImage``
477 Image of the overscan region with the overscan
478 correction applied (`lsst.afw.image.Image`). This
479 quantity is used to estimate the amplifier read noise
480 empirically.
481
482 Returns
483 -------
484 outputStats : `dict` [`str`, [`dict` [`str`, `float`]]]
485 Dictionary of measurements, keyed by amplifier name and
486 statistics segment.
487 """
488 outputStats = {}
489
490 detector = inputExp.getDetector()
491 kernel = self.makeKernel(self.config.bandingKernelSize)
492
493 outputStats["AMP_BANDING"] = []
494 for amp, overscanData in zip(detector.getAmplifiers(), overscans):
495 overscanFit = np.array(overscanData.overscanFit)
496 overscanArray = overscanData.overscanImage.image.array
497 rawOverscan = np.mean(overscanArray + overscanFit, axis=1)
498
499 smoothedOverscan = np.convolve(rawOverscan, kernel, mode="valid")
500
501 low, high = np.quantile(smoothedOverscan, [self.config.bandingFractionLow,
502 self.config.bandingFractionHigh])
503 outputStats["AMP_BANDING"].append(float(high - low))
504
505 if self.config.bandingUseHalfDetector:
506 fullLength = len(outputStats["AMP_BANDING"])
507 outputStats["DET_BANDING"] = float(np.nanmedian(outputStats["AMP_BANDING"][0:fullLength//2]))
508 else:
509 outputStats["DET_BANDING"] = float(np.nanmedian(outputStats["AMP_BANDING"]))
510
511 return outputStats
512

◆ measureBiasShifts()

lsst.ip.isr.isrStatistics.IsrStatisticsTask.measureBiasShifts ( self,
inputExp,
overscanResults )
Measure number of bias shifts from overscan data.

Parameters
----------
inputExp : `lsst.afw.image.Exposure`
    Exposure to measure.
overscans : `list` [`lsst.pipe.base.Struct`]
    List of overscan results.  Expected fields are:

    ``imageFit``
        Value or fit subtracted from the amplifier image data
        (scalar or `lsst.afw.image.Image`).
    ``overscanFit``
        Value or fit subtracted from the overscan image data
        (scalar or `lsst.afw.image.Image`).
    ``overscanImage``
        Image of the overscan region with the overscan
        correction applied (`lsst.afw.image.Image`). This
        quantity is used to estimate the amplifier read noise
        empirically.

Returns
-------
outputStats : `dict` [`str`, [`dict` [`str`, `float`]]]
    Dictionary of measurements, keyed by amplifier name and
    statistics segment.

Notes
-----
Based on eop_pipe implementation:
https://github.com/lsst-camera-dh/eo_pipe/blob/main/python/lsst/eo/pipe/biasShiftsTask.py  # noqa: E501 W505

Definition at line 641 of file isrStatistics.py.

641 def measureBiasShifts(self, inputExp, overscanResults):
642 """Measure number of bias shifts from overscan data.
643
644 Parameters
645 ----------
646 inputExp : `lsst.afw.image.Exposure`
647 Exposure to measure.
648 overscans : `list` [`lsst.pipe.base.Struct`]
649 List of overscan results. Expected fields are:
650
651 ``imageFit``
652 Value or fit subtracted from the amplifier image data
653 (scalar or `lsst.afw.image.Image`).
654 ``overscanFit``
655 Value or fit subtracted from the overscan image data
656 (scalar or `lsst.afw.image.Image`).
657 ``overscanImage``
658 Image of the overscan region with the overscan
659 correction applied (`lsst.afw.image.Image`). This
660 quantity is used to estimate the amplifier read noise
661 empirically.
662
663 Returns
664 -------
665 outputStats : `dict` [`str`, [`dict` [`str`, `float`]]]
666 Dictionary of measurements, keyed by amplifier name and
667 statistics segment.
668
669 Notes
670 -----
671 Based on eop_pipe implementation:
672 https://github.com/lsst-camera-dh/eo_pipe/blob/main/python/lsst/eo/pipe/biasShiftsTask.py # noqa: E501 W505
673 """
674 outputStats = {}
675
676 detector = inputExp.getDetector()
677 for amp, overscans in zip(detector, overscanResults):
678 ampStats = {}
679 # Add fit back to data
680 rawOverscan = overscans.overscanImage.image.array + overscans.overscanFit
681
682 # Collapse array, skipping first three columns
683 rawOverscan = np.mean(rawOverscan[:, self.config.biasShiftColumnSkip:], axis=1)
684
685 # Scan for shifts
686 noise, shift_peaks = self._scan_for_shifts(rawOverscan)
687 ampStats["LOCAL_NOISE"] = float(noise)
688 ampStats["BIAS_SHIFTS"] = shift_peaks
689
690 outputStats[amp.getName()] = ampStats
691 return outputStats
692

◆ measureCti()

lsst.ip.isr.isrStatistics.IsrStatisticsTask.measureCti ( self,
inputExp,
overscans,
gains )
Task to measure CTI statistics.

Parameters
----------
inputExp : `lsst.afw.image.Exposure`
    Exposure to measure.
overscans : `list` [`lsst.pipe.base.Struct`]
    List of overscan results (expects base units of adu).
    Expected fields are:

    ``imageFit``
        Value or fit subtracted from the amplifier image data
        (scalar or `lsst.afw.image.Image`).
    ``overscanFit``
        Value or fit subtracted from the overscan image data
        (scalar or `lsst.afw.image.Image`).
    ``overscanImage``
        Image of the overscan region with the overscan
        correction applied (`lsst.afw.image.Image`). This
        quantity is used to estimate the amplifier read noise
        empirically.
gains : `dict` [`str` `float`]
    Dictionary of per-amplifier gains, indexed by amplifier name.

Returns
-------
outputStats : `dict` [`str`, [`dict` [`str`, `float`]]]
    Dictionary of measurements, keyed by amplifier name and
    statistics segment. Everything in units based on electron.

Definition at line 320 of file isrStatistics.py.

320 def measureCti(self, inputExp, overscans, gains):
321 """Task to measure CTI statistics.
322
323 Parameters
324 ----------
325 inputExp : `lsst.afw.image.Exposure`
326 Exposure to measure.
327 overscans : `list` [`lsst.pipe.base.Struct`]
328 List of overscan results (expects base units of adu).
329 Expected fields are:
330
331 ``imageFit``
332 Value or fit subtracted from the amplifier image data
333 (scalar or `lsst.afw.image.Image`).
334 ``overscanFit``
335 Value or fit subtracted from the overscan image data
336 (scalar or `lsst.afw.image.Image`).
337 ``overscanImage``
338 Image of the overscan region with the overscan
339 correction applied (`lsst.afw.image.Image`). This
340 quantity is used to estimate the amplifier read noise
341 empirically.
342 gains : `dict` [`str` `float`]
343 Dictionary of per-amplifier gains, indexed by amplifier name.
344
345 Returns
346 -------
347 outputStats : `dict` [`str`, [`dict` [`str`, `float`]]]
348 Dictionary of measurements, keyed by amplifier name and
349 statistics segment. Everything in units based on electron.
350 """
351 outputStats = {}
352
353 detector = inputExp.getDetector()
354 image = inputExp.image
355
356 # It only makes sense to measure CTI in electron units.
357 # Make it so.
358 imageUnits = inputExp.getMetadata().get("LSST ISR UNITS")
359 applyGain = False
360 if imageUnits == "adu":
361 applyGain = True
362
363 # Check if the image is trimmed.
364 isTrimmed = isTrimmedExposure(inputExp)
365
366 # Ensure we have the same number of overscans as amplifiers.
367 assert len(overscans) == len(detector.getAmplifiers())
368
369 with gainContext(inputExp, image, applyGain, gains, isTrimmed=isTrimmed):
370 for ampIter, amp in enumerate(detector.getAmplifiers()):
371 ampStats = {}
372 readoutCorner = amp.getReadoutCorner()
373
374 # Full data region.
375 dataRegion = image[amp.getBBox() if isTrimmed else amp.getRawDataBBox()]
376
377 # Get the mean of the image
378 ampStats["IMAGE_MEAN"] = afwMath.makeStatistics(dataRegion, self.statType,
379 self.statControl).getValue()
380
381 # First and last image columns.
382 pixelA = afwMath.makeStatistics(dataRegion.array[:, 0],
383 self.statType,
384 self.statControl).getValue()
385 pixelZ = afwMath.makeStatistics(dataRegion.array[:, -1],
386 self.statType,
387 self.statControl).getValue()
388
389 # We want these relative to the readout corner. If that's
390 # on the right side, we need to swap them.
391 if readoutCorner in (ReadoutCorner.LR, ReadoutCorner.UR):
392 ampStats["FIRST_MEAN"] = pixelZ
393 ampStats["LAST_MEAN"] = pixelA
394 else:
395 ampStats["FIRST_MEAN"] = pixelA
396 ampStats["LAST_MEAN"] = pixelZ
397
398 # Measure the columns of the overscan.
399 if overscans[ampIter] is None:
400 # The amplifier is likely entirely bad, and needs to
401 # be skipped.
402 self.log.warning("No overscan information available for ISR statistics for amp %s.",
403 amp.getName())
404 nCols = amp.getRawSerialOverscanBBox().getWidth()
405 ampStats["OVERSCAN_COLUMNS"] = np.full((nCols, ), np.nan)
406 ampStats["OVERSCAN_VALUES"] = np.full((nCols, ), np.nan)
407 else:
408 overscanImage = overscans[ampIter].overscanImage
409
410 columns = []
411 values = []
412 for column in range(0, overscanImage.getWidth()):
413 # If overscan.doParallelOverscan=True, the
414 # overscanImage will contain both the serial
415 # and parallel overscan regions.
416 # Only the serial CTI correction has been
417 # implemented, so we must select only the
418 # serial overscan rows for a given column.
419 nRows = amp.getRawSerialOverscanBBox().getHeight()
420 osMean = afwMath.makeStatistics(overscanImage.image.array[:nRows, column],
421 self.statType, self.statControl).getValue()
422 columns.append(column)
423 # The overscan input is always in adu, but it only
424 # makes sense to measure CTI in electron units.
425 values.append(osMean * gains[amp.getName()])
426
427 # We want these relative to the readout corner. If that's
428 # on the right side, we need to swap them.
429 if readoutCorner in (ReadoutCorner.LR, ReadoutCorner.UR):
430 ampStats["OVERSCAN_COLUMNS"] = list(reversed(columns))
431 ampStats["OVERSCAN_VALUES"] = list(reversed(values))
432 else:
433 ampStats["OVERSCAN_COLUMNS"] = columns
434 ampStats["OVERSCAN_VALUES"] = values
435
436 outputStats[amp.getName()] = ampStats
437
438 return outputStats
439
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

◆ measureDivisaderoStatistics()

lsst.ip.isr.isrStatistics.IsrStatisticsTask.measureDivisaderoStatistics ( self,
inputExp,
** kwargs )
Measure Max Divisadero Tearing effect per amp.

Parameters
----------
inputExp : `lsst.afw.image.Exposure`
    Exposure to measure. Usually a flat.
**kwargs :
    The flat will be selected from here.

Returns
-------
outputStats : `dict` [`str`, [`dict` [`str`, `float`]]]
    Dictionary of measurements, keyed by amplifier name and
    statistics segment.
    Measurements include

    - DIVISADERO_PROFILE: Robust mean of rows between
      divisaderoProjection<Maximum|Minumum> on readout edge of ccd
      normalized by a linear fit to the same rows.
    - DIVISADERO_MAX_PAIR: Tuple of maximum of the absolute values of
      the DIVISADERO_PROFILE, for number of pixels (specified by
      divisaderoNumImpactPixels on left and right side of amp.
    - DIVISADERO_MAX: Maximum of the absolute values of the
      the DIVISADERO_PROFILE, for the divisaderoNumImpactPixels on
      boundaries of neighboring amps (including the pixels in those
      neighborboring amps).

Definition at line 836 of file isrStatistics.py.

836 def measureDivisaderoStatistics(self, inputExp, **kwargs):
837 """Measure Max Divisadero Tearing effect per amp.
838
839 Parameters
840 ----------
841 inputExp : `lsst.afw.image.Exposure`
842 Exposure to measure. Usually a flat.
843 **kwargs :
844 The flat will be selected from here.
845
846 Returns
847 -------
848 outputStats : `dict` [`str`, [`dict` [`str`, `float`]]]
849 Dictionary of measurements, keyed by amplifier name and
850 statistics segment.
851 Measurements include
852
853 - DIVISADERO_PROFILE: Robust mean of rows between
854 divisaderoProjection<Maximum|Minumum> on readout edge of ccd
855 normalized by a linear fit to the same rows.
856 - DIVISADERO_MAX_PAIR: Tuple of maximum of the absolute values of
857 the DIVISADERO_PROFILE, for number of pixels (specified by
858 divisaderoNumImpactPixels on left and right side of amp.
859 - DIVISADERO_MAX: Maximum of the absolute values of the
860 the DIVISADERO_PROFILE, for the divisaderoNumImpactPixels on
861 boundaries of neighboring amps (including the pixels in those
862 neighborboring amps).
863 """
864 outputStats = {}
865
866 for amp in inputExp.getDetector():
867 # Copy unneeded if we do not ever modify the array by flipping
868 ampArray = inputExp.image[amp.getBBox()].array
869 # slice the outer top or bottom of the amp: the readout side
870 if amp.getReadoutCorner().name in ('UL', 'UR'):
871 minRow = amp.getBBox().getHeight() - self.config.divisaderoProjectionMaximum
872 maxRow = amp.getBBox().getHeight() - self.config.divisaderoProjectionMinimum
873 else:
874 minRow = self.config.divisaderoProjectionMinimum
875 maxRow = self.config.divisaderoProjectionMaximum
876
877 segment = slice(minRow, maxRow)
878 projection, _, _ = astropy.stats.sigma_clipped_stats(ampArray[segment, :], axis=0)
879
880 ampStats = {}
881 projection /= np.median(projection)
882 columns = np.arange(len(projection))
883
884 segment = slice(self.config.divisaderoEdgePixels, -self.config.divisaderoEdgePixels)
885 model = np.polyfit(columns[segment], projection[segment], 1)
886 modelProjection = model[0] * columns + model[1]
887 divisaderoProfile = projection / modelProjection
888
889 # look for max at the edges:
890 leftMax = np.nanmax(np.abs(divisaderoProfile[0:self.config.divisaderoNumImpactPixels] - 1.0))
891 rightMax = np.nanmax(np.abs(divisaderoProfile[-self.config.divisaderoNumImpactPixels:] - 1.0))
892
893 ampStats['DIVISADERO_PROFILE'] = np.array(divisaderoProfile).tolist()
894 ampStats['DIVISADERO_MAX_PAIR'] = [leftMax, rightMax]
895 outputStats[amp.getName()] = ampStats
896
897 detector = inputExp.getDetector()
898 xCenters = [amp.getBBox().getCenterX() for amp in detector]
899 yCenters = [amp.getBBox().getCenterY() for amp in detector]
900 xIndices = np.ceil(xCenters / np.min(xCenters) / 2).astype(int) - 1
901 yIndices = np.ceil(yCenters / np.min(yCenters) / 2).astype(int) - 1
902 ampIds = np.zeros((len(set(yIndices)), len(set(xIndices))), dtype=int)
903 for ampId, xIndex, yIndex in zip(np.arange(len(detector)), xIndices, yIndices):
904 ampIds[yIndex, xIndex] = ampId
905
906 # Loop over amps again because the DIVISIDERO_MAX will be the max
907 # of the profile on its boundary with its neighboring amps
908 for i, amp in enumerate(detector):
909 y, x = np.where(ampIds == i)
910 end = ampIds.shape[1] - 1
911 xInd = x[0]
912 yInd = y[0]
913 thisAmpsPair = outputStats[amp.getName()]['DIVISADERO_MAX_PAIR']
914
915 if x == 0:
916 # leftmost amp: take the max of your right side and
917 myMax = thisAmpsPair[1]
918 # your neighbor's left side
919 neighborMax = outputStats[detector[ampIds[yInd, 1]].getName()]['DIVISADERO_MAX_PAIR'][0]
920 elif x == end:
921 # rightmost amp: take the max of your left side and
922 myMax = thisAmpsPair[0]
923 # your neighbor's right side
924 neighborMax = outputStats[detector[ampIds[yInd, end - 1]].getName()]['DIVISADERO_MAX_PAIR'][1]
925 else:
926 # Middle amp: take the max of both your own sides and the
927 myMax = max(thisAmpsPair)
928 leftName = detector[ampIds[yInd, max(xInd - 1, 0)]].getName()
929 rightName = detector[ampIds[yInd, min(xInd + 1, ampIds.shape[1] - 1)]].getName()
930 # right side of the neighbor to your left
931 # and left side of your neighbor to your right
932 neighborMax = max(outputStats[leftName]['DIVISADERO_MAX_PAIR'][1],
933 outputStats[rightName]['DIVISADERO_MAX_PAIR'][0])
934
935 divisaderoMax = max([myMax, neighborMax])
936 outputStats[amp.getName()]['DIVISADERO_MAX'] = divisaderoMax
937
938 return outputStats
int min
int max

◆ measureProjectionStatistics()

lsst.ip.isr.isrStatistics.IsrStatisticsTask.measureProjectionStatistics ( self,
inputExp,
overscans )
Task to measure metrics from image slicing.

Parameters
----------
inputExp : `lsst.afw.image.Exposure`
    Exposure to measure.
overscans : `list` [`lsst.pipe.base.Struct`]
    List of overscan results.  Expected fields are:

    ``imageFit``
        Value or fit subtracted from the amplifier image data
        (scalar or `lsst.afw.image.Image`).
    ``overscanFit``
        Value or fit subtracted from the overscan image data
        (scalar or `lsst.afw.image.Image`).
    ``overscanImage``
        Image of the overscan region with the overscan
        correction applied (`lsst.afw.image.Image`). This
        quantity is used to estimate the amplifier read noise
        empirically.

Returns
-------
outputStats : `dict` [`str`, [`dict` [`str`, `float`]]]
    Dictionary of measurements, keyed by amplifier name and
    statistics segment.

Definition at line 513 of file isrStatistics.py.

513 def measureProjectionStatistics(self, inputExp, overscans):
514 """Task to measure metrics from image slicing.
515
516 Parameters
517 ----------
518 inputExp : `lsst.afw.image.Exposure`
519 Exposure to measure.
520 overscans : `list` [`lsst.pipe.base.Struct`]
521 List of overscan results. Expected fields are:
522
523 ``imageFit``
524 Value or fit subtracted from the amplifier image data
525 (scalar or `lsst.afw.image.Image`).
526 ``overscanFit``
527 Value or fit subtracted from the overscan image data
528 (scalar or `lsst.afw.image.Image`).
529 ``overscanImage``
530 Image of the overscan region with the overscan
531 correction applied (`lsst.afw.image.Image`). This
532 quantity is used to estimate the amplifier read noise
533 empirically.
534
535 Returns
536 -------
537 outputStats : `dict` [`str`, [`dict` [`str`, `float`]]]
538 Dictionary of measurements, keyed by amplifier name and
539 statistics segment.
540 """
541 outputStats = {}
542
543 detector = inputExp.getDetector()
544 kernel = self.makeKernel(self.config.projectionKernelSize)
545
546 outputStats["AMP_VPROJECTION"] = {}
547 outputStats["AMP_HPROJECTION"] = {}
548 convolveMode = "valid"
549 if self.config.doProjectionFft:
550 outputStats["AMP_VFFT_REAL"] = {}
551 outputStats["AMP_VFFT_IMAG"] = {}
552 outputStats["AMP_HFFT_REAL"] = {}
553 outputStats["AMP_HFFT_IMAG"] = {}
554 convolveMode = "same"
555
556 for amp in detector.getAmplifiers():
557 ampArray = inputExp.image[amp.getBBox()].array
558
559 horizontalProjection = np.mean(ampArray, axis=0)
560 verticalProjection = np.mean(ampArray, axis=1)
561
562 horizontalProjection = np.convolve(horizontalProjection, kernel, mode=convolveMode)
563 verticalProjection = np.convolve(verticalProjection, kernel, mode=convolveMode)
564
565 outputStats["AMP_HPROJECTION"][amp.getName()] = horizontalProjection.tolist()
566 outputStats["AMP_VPROJECTION"][amp.getName()] = verticalProjection.tolist()
567
568 if self.config.doProjectionFft:
569 horizontalWindow = np.ones_like(horizontalProjection)
570 verticalWindow = np.ones_like(verticalProjection)
571 if self.config.projectionFftWindow == "NONE":
572 pass
573 elif self.config.projectionFftWindow == "HAMMING":
574 horizontalWindow = hamming(len(horizontalProjection))
575 verticalWindow = hamming(len(verticalProjection))
576 elif self.config.projectionFftWindow == "HANN":
577 horizontalWindow = hann(len(horizontalProjection))
578 verticalWindow = hann(len(verticalProjection))
579 elif self.config.projectionFftWindow == "GAUSSIAN":
580 horizontalWindow = gaussian(len(horizontalProjection))
581 verticalWindow = gaussian(len(verticalProjection))
582 else:
583 raise RuntimeError(f"Invalid window function: {self.config.projectionFftWindow}")
584
585 horizontalFFT = np.fft.rfft(np.multiply(horizontalProjection, horizontalWindow))
586 verticalFFT = np.fft.rfft(np.multiply(verticalProjection, verticalWindow))
587
588 outputStats["AMP_HFFT_REAL"][amp.getName()] = np.real(horizontalFFT).tolist()
589 outputStats["AMP_HFFT_IMAG"][amp.getName()] = np.imag(horizontalFFT).tolist()
590 outputStats["AMP_VFFT_REAL"][amp.getName()] = np.real(verticalFFT).tolist()
591 outputStats["AMP_VFFT_IMAG"][amp.getName()] = np.imag(verticalFFT).tolist()
592
593 return outputStats
594

◆ run()

lsst.ip.isr.isrStatistics.IsrStatisticsTask.run ( self,
inputExp,
ptc = None,
overscanResults = None,
** kwargs )
Task to run arbitrary statistics.

The statistics should be measured by individual methods, and
add to the dictionary in the return struct.

Parameters
----------
inputExp : `lsst.afw.image.Exposure`
    The exposure to measure.
ptc : `lsst.ip.isr.PtcDataset`, optional
    A PTC object containing gains to use.
overscanResults : `list` [`lsst.pipe.base.Struct`], optional
    List of overscan results.  Expected fields are:

    ``imageFit``
        Value or fit subtracted from the amplifier image data
        (scalar or `lsst.afw.image.Image`).
    ``overscanFit``
        Value or fit subtracted from the overscan image data
        (scalar or `lsst.afw.image.Image`).
    ``overscanImage``
        Image of the overscan region with the overscan
        correction applied (`lsst.afw.image.Image`). This
        quantity is used to estimate the amplifier read noise
        empirically.
**kwargs :
     Keyword arguments.  Calibrations being passed in should
     have an entry here.

Returns
-------
resultStruct : `lsst.pipe.base.Struct`
    Contains the measured statistics as a dict stored in a
    field named ``results``.

Raises
------
RuntimeError
    Raised if the amplifier gains could not be found.

Definition at line 228 of file isrStatistics.py.

228 def run(self, inputExp, ptc=None, overscanResults=None, **kwargs):
229 """Task to run arbitrary statistics.
230
231 The statistics should be measured by individual methods, and
232 add to the dictionary in the return struct.
233
234 Parameters
235 ----------
236 inputExp : `lsst.afw.image.Exposure`
237 The exposure to measure.
238 ptc : `lsst.ip.isr.PtcDataset`, optional
239 A PTC object containing gains to use.
240 overscanResults : `list` [`lsst.pipe.base.Struct`], optional
241 List of overscan results. Expected fields are:
242
243 ``imageFit``
244 Value or fit subtracted from the amplifier image data
245 (scalar or `lsst.afw.image.Image`).
246 ``overscanFit``
247 Value or fit subtracted from the overscan image data
248 (scalar or `lsst.afw.image.Image`).
249 ``overscanImage``
250 Image of the overscan region with the overscan
251 correction applied (`lsst.afw.image.Image`). This
252 quantity is used to estimate the amplifier read noise
253 empirically.
254 **kwargs :
255 Keyword arguments. Calibrations being passed in should
256 have an entry here.
257
258 Returns
259 -------
260 resultStruct : `lsst.pipe.base.Struct`
261 Contains the measured statistics as a dict stored in a
262 field named ``results``.
263
264 Raises
265 ------
266 RuntimeError
267 Raised if the amplifier gains could not be found.
268 """
269 # Find gains.
270 detector = inputExp.getDetector()
271 if ptc is not None:
272 gains = ptc.gain
273 elif detector is not None:
274 gains = {amp.getName(): amp.getGain() for amp in detector.getAmplifiers()}
275 else:
276 raise RuntimeError("No source of gains provided.")
277
278 ctiResults = None
279 if self.config.doCtiStatistics:
280 ctiResults = self.measureCti(inputExp, overscanResults, gains)
281
282 bandingResults = None
283 if self.config.doBandingStatistics:
284 bandingResults = self.measureBanding(inputExp, overscanResults)
285
286 projectionResults = None
287 if self.config.doProjectionStatistics:
288 projectionResults = self.measureProjectionStatistics(inputExp, overscanResults)
289
290 divisaderoResults = None
291 if self.config.doDivisaderoStatistics:
292 divisaderoResults = self.measureDivisaderoStatistics(inputExp, **kwargs)
293
294 calibDistributionResults = None
295 if self.config.doCopyCalibDistributionStatistics:
296 calibDistributionResults = self.copyCalibDistributionStatistics(inputExp, **kwargs)
297
298 biasShiftResults = None
299 if self.config.doBiasShiftStatistics:
300 biasShiftResults = self.measureBiasShifts(inputExp, overscanResults)
301
302 ampCorrelationResults = None
303 if self.config.doAmplifierCorrelationStatistics:
304 ampCorrelationResults = self.measureAmpCorrelations(inputExp, overscanResults)
305
306 mjd = inputExp.getMetadata().get("MJD", None)
307
308 return pipeBase.Struct(
309 results={"CTI": ctiResults,
310 "BANDING": bandingResults,
311 "PROJECTION": projectionResults,
312 "CALIBDIST": calibDistributionResults,
313 "BIASSHIFT": biasShiftResults,
314 "AMPCORR": ampCorrelationResults,
315 "MJD": mjd,
316 'DIVISADERO': divisaderoResults,
317 },
318 )
319

Member Data Documentation

◆ _DefaultName

str lsst.ip.isr.isrStatistics.IsrStatisticsTask._DefaultName = "isrStatistics"
staticprotected

Definition at line 220 of file isrStatistics.py.

◆ ConfigClass

lsst.ip.isr.isrStatistics.IsrStatisticsTask.ConfigClass = IsrStatisticsTaskConfig
static

Definition at line 219 of file isrStatistics.py.

◆ statControl

lsst.ip.isr.isrStatistics.IsrStatisticsTask.statControl

Definition at line 224 of file isrStatistics.py.

◆ statType

lsst.ip.isr.isrStatistics.IsrStatisticsTask.statType

Definition at line 226 of file isrStatistics.py.


The documentation for this class was generated from the following file: