34 from contextlib
import contextmanager
36 from .overscan
import OverscanCorrectionTask, OverscanCorrectionTaskConfig
40 """Make a double Gaussian PSF. 45 FWHM of double Gaussian smoothing kernel. 49 psf : `lsst.meas.algorithms.DoubleGaussianPsf` 50 The created smoothing kernel. 52 ksize = 4*int(fwhm) + 1
53 return measAlg.DoubleGaussianPsf(ksize, ksize, fwhm/(2*math.sqrt(2*math.log(2))))
57 """Make a transposed copy of a masked image. 61 maskedImage : `lsst.afw.image.MaskedImage` 66 transposed : `lsst.afw.image.MaskedImage` 67 The transposed copy of the input image. 69 transposed = maskedImage.Factory(
lsst.geom.Extent2I(maskedImage.getHeight(), maskedImage.getWidth()))
70 transposed.getImage().getArray()[:] = maskedImage.getImage().getArray().T
71 transposed.getMask().getArray()[:] = maskedImage.getMask().getArray().T
72 transposed.getVariance().getArray()[:] = maskedImage.getVariance().getArray().T
77 """Interpolate over defects specified in a defect list. 81 maskedImage : `lsst.afw.image.MaskedImage` 83 defectList : `lsst.meas.algorithms.Defects` 84 List of defects to interpolate over. 86 FWHM of double Gaussian smoothing kernel. 87 fallbackValue : scalar, optional 88 Fallback value if an interpolated value cannot be determined. 89 If None, then the clipped mean of the image is used. 92 if fallbackValue
is None:
93 fallbackValue = afwMath.makeStatistics(maskedImage.getImage(), afwMath.MEANCLIP).getValue()
94 if 'INTRP' not in maskedImage.getMask().getMaskPlaneDict():
95 maskedImage.getMask().addMaskPlane(
'INTRP')
96 measAlg.interpolateOverDefects(maskedImage, psf, defectList, fallbackValue,
True)
101 """Mask pixels based on threshold detection. 105 maskedImage : `lsst.afw.image.MaskedImage` 106 Image to process. Only the mask plane is updated. 109 growFootprints : scalar, optional 110 Number of pixels to grow footprints of detected regions. 111 maskName : str, optional 112 Mask plane name, or list of names to convert 116 defectList : `lsst.meas.algorithms.Defects` 117 Defect list constructed from pixels above the threshold. 120 thresh = afwDetection.Threshold(threshold)
121 fs = afwDetection.FootprintSet(maskedImage, thresh)
123 if growFootprints > 0:
124 fs = afwDetection.FootprintSet(fs, rGrow=growFootprints, isotropic=
False)
125 fpList = fs.getFootprints()
128 mask = maskedImage.getMask()
129 bitmask = mask.getPlaneBitMask(maskName)
130 afwDetection.setMaskFromFootprintList(mask, fpList, bitmask)
132 return measAlg.Defects.fromFootprintList(fpList)
135 def growMasks(mask, radius=0, maskNameList=['BAD'], maskValue="BAD"):
136 """Grow a mask by an amount and add to the requested plane. 140 mask : `lsst.afw.image.Mask` 141 Mask image to process. 143 Amount to grow the mask. 144 maskNameList : `str` or `list` [`str`] 145 Mask names that should be grown. 147 Mask plane to assign the newly masked pixels to. 150 thresh = afwDetection.Threshold(mask.getPlaneBitMask(maskNameList), afwDetection.Threshold.BITMASK)
151 fpSet = afwDetection.FootprintSet(mask, thresh)
152 fpSet = afwDetection.FootprintSet(fpSet, rGrow=radius, isotropic=
False)
153 fpSet.setMask(mask, maskValue)
157 maskNameList=['SAT'], fallbackValue=None):
158 """Interpolate over defects identified by a particular set of mask planes. 162 maskedImage : `lsst.afw.image.MaskedImage` 165 FWHM of double Gaussian smoothing kernel. 166 growSaturatedFootprints : scalar, optional 167 Number of pixels to grow footprints for saturated pixels. 168 maskNameList : `List` of `str`, optional 170 fallbackValue : scalar, optional 171 Value of last resort for interpolation. 173 mask = maskedImage.getMask()
175 if growSaturatedFootprints > 0
and "SAT" in maskNameList:
178 growMasks(mask, radius=growSaturatedFootprints, maskNameList=[
'SAT'], maskValue=
"SAT")
180 thresh = afwDetection.Threshold(mask.getPlaneBitMask(maskNameList), afwDetection.Threshold.BITMASK)
181 fpSet = afwDetection.FootprintSet(mask, thresh)
182 defectList = measAlg.Defects.fromFootprintList(fpSet.getFootprints())
189 def saturationCorrection(maskedImage, saturation, fwhm, growFootprints=1, interpolate=True, maskName='SAT',
191 """Mark saturated pixels and optionally interpolate over them 195 maskedImage : `lsst.afw.image.MaskedImage` 198 Saturation level used as the detection threshold. 200 FWHM of double Gaussian smoothing kernel. 201 growFootprints : scalar, optional 202 Number of pixels to grow footprints of detected regions. 203 interpolate : Bool, optional 204 If True, saturated pixels are interpolated over. 205 maskName : str, optional 207 fallbackValue : scalar, optional 208 Value of last resort for interpolation. 211 maskedImage=maskedImage,
212 threshold=saturation,
213 growFootprints=growFootprints,
223 """Compute number of edge trim pixels to match the calibration data. 225 Use the dimension difference between the raw exposure and the 226 calibration exposure to compute the edge trim pixels. This trim 227 is applied symmetrically, with the same number of pixels masked on 232 rawMaskedImage : `lsst.afw.image.MaskedImage` 234 calibMaskedImage : `lsst.afw.image.MaskedImage` 235 Calibration image to draw new bounding box from. 239 replacementMaskedImage : `lsst.afw.image.MaskedImage` 240 ``rawMaskedImage`` trimmed to the appropriate size 244 Rasied if ``rawMaskedImage`` cannot be symmetrically trimmed to 245 match ``calibMaskedImage``. 247 nx, ny = rawMaskedImage.getBBox().getDimensions() - calibMaskedImage.getBBox().getDimensions()
249 raise RuntimeError(
"Raw and calib maskedImages are trimmed differently in X and Y.")
251 raise RuntimeError(
"Calibration maskedImage is trimmed unevenly in X.")
253 raise RuntimeError(
"Calibration maskedImage is larger than raw data.")
257 replacementMaskedImage = rawMaskedImage[nEdge:-nEdge, nEdge:-nEdge, afwImage.LOCAL]
258 SourceDetectionTask.setEdgeBits(
260 replacementMaskedImage.getBBox(),
261 rawMaskedImage.getMask().getPlaneBitMask(
"EDGE")
264 replacementMaskedImage = rawMaskedImage
266 return replacementMaskedImage
270 """Apply bias correction in place. 274 maskedImage : `lsst.afw.image.MaskedImage` 275 Image to process. The image is modified by this method. 276 biasMaskedImage : `lsst.afw.image.MaskedImage` 277 Bias image of the same size as ``maskedImage`` 278 trimToFit : `Bool`, optional 279 If True, raw data is symmetrically trimmed to match 285 Raised if ``maskedImage`` and ``biasMaskedImage`` do not have 292 if maskedImage.getBBox(afwImage.LOCAL) != biasMaskedImage.getBBox(afwImage.LOCAL):
293 raise RuntimeError(
"maskedImage bbox %s != biasMaskedImage bbox %s" %
294 (maskedImage.getBBox(afwImage.LOCAL), biasMaskedImage.getBBox(afwImage.LOCAL)))
295 maskedImage -= biasMaskedImage
298 def darkCorrection(maskedImage, darkMaskedImage, expScale, darkScale, invert=False, trimToFit=False):
299 """Apply dark correction in place. 303 maskedImage : `lsst.afw.image.MaskedImage` 304 Image to process. The image is modified by this method. 305 darkMaskedImage : `lsst.afw.image.MaskedImage` 306 Dark image of the same size as ``maskedImage``. 308 Dark exposure time for ``maskedImage``. 310 Dark exposure time for ``darkMaskedImage``. 311 invert : `Bool`, optional 312 If True, re-add the dark to an already corrected image. 313 trimToFit : `Bool`, optional 314 If True, raw data is symmetrically trimmed to match 320 Raised if ``maskedImage`` and ``darkMaskedImage`` do not have 325 The dark correction is applied by calculating: 326 maskedImage -= dark * expScaling / darkScaling 331 if maskedImage.getBBox(afwImage.LOCAL) != darkMaskedImage.getBBox(afwImage.LOCAL):
332 raise RuntimeError(
"maskedImage bbox %s != darkMaskedImage bbox %s" %
333 (maskedImage.getBBox(afwImage.LOCAL), darkMaskedImage.getBBox(afwImage.LOCAL)))
335 scale = expScale / darkScale
337 maskedImage.scaledMinus(scale, darkMaskedImage)
339 maskedImage.scaledPlus(scale, darkMaskedImage)
343 """Set the variance plane based on the image plane. 347 maskedImage : `lsst.afw.image.MaskedImage` 348 Image to process. The variance plane is modified. 350 The amplifier gain in electrons/ADU. 352 The amplifier read nmoise in ADU/pixel. 354 var = maskedImage.getVariance()
355 var[:] = maskedImage.getImage()
360 def flatCorrection(maskedImage, flatMaskedImage, scalingType, userScale=1.0, invert=False, trimToFit=False):
361 """Apply flat correction in place. 365 maskedImage : `lsst.afw.image.MaskedImage` 366 Image to process. The image is modified. 367 flatMaskedImage : `lsst.afw.image.MaskedImage` 368 Flat image of the same size as ``maskedImage`` 370 Flat scale computation method. Allowed values are 'MEAN', 372 userScale : scalar, optional 373 Scale to use if ``scalingType``='USER'. 374 invert : `Bool`, optional 375 If True, unflatten an already flattened image. 376 trimToFit : `Bool`, optional 377 If True, raw data is symmetrically trimmed to match 383 Raised if ``maskedImage`` and ``flatMaskedImage`` do not have 384 the same size or if ``scalingType`` is not an allowed value. 389 if maskedImage.getBBox(afwImage.LOCAL) != flatMaskedImage.getBBox(afwImage.LOCAL):
390 raise RuntimeError(
"maskedImage bbox %s != flatMaskedImage bbox %s" %
391 (maskedImage.getBBox(afwImage.LOCAL), flatMaskedImage.getBBox(afwImage.LOCAL)))
396 if scalingType
in (
'MEAN',
'MEDIAN'):
397 scalingType = afwMath.stringToStatisticsProperty(scalingType)
398 flatScale = afwMath.makeStatistics(flatMaskedImage.image, scalingType).getValue()
399 elif scalingType ==
'USER':
400 flatScale = userScale
402 raise RuntimeError(
'%s : %s not implemented' % (
"flatCorrection", scalingType))
405 maskedImage.scaledDivides(1.0/flatScale, flatMaskedImage)
407 maskedImage.scaledMultiplies(1.0/flatScale, flatMaskedImage)
411 """Apply illumination correction in place. 415 maskedImage : `lsst.afw.image.MaskedImage` 416 Image to process. The image is modified. 417 illumMaskedImage : `lsst.afw.image.MaskedImage` 418 Illumination correction image of the same size as ``maskedImage``. 420 Scale factor for the illumination correction. 421 trimToFit : `Bool`, optional 422 If True, raw data is symmetrically trimmed to match 428 Raised if ``maskedImage`` and ``illumMaskedImage`` do not have 434 if maskedImage.getBBox(afwImage.LOCAL) != illumMaskedImage.getBBox(afwImage.LOCAL):
435 raise RuntimeError(
"maskedImage bbox %s != illumMaskedImage bbox %s" %
436 (maskedImage.getBBox(afwImage.LOCAL), illumMaskedImage.getBBox(afwImage.LOCAL)))
438 maskedImage.scaledDivides(1.0/illumScale, illumMaskedImage)
441 def overscanCorrection(ampMaskedImage, overscanImage, fitType='MEDIAN', order=1, collapseRej=3.0,
442 statControl=None, overscanIsInt=True):
443 """Apply overscan correction in place. 447 ampMaskedImage : `lsst.afw.image.MaskedImage` 448 Image of amplifier to correct; modified. 449 overscanImage : `lsst.afw.image.Image` or `lsst.afw.image.MaskedImage` 450 Image of overscan; modified. 452 Type of fit for overscan correction. May be one of: 454 - ``MEAN``: use mean of overscan. 455 - ``MEANCLIP``: use clipped mean of overscan. 456 - ``MEDIAN``: use median of overscan. 457 - ``MEDIAN_PER_ROW``: use median per row of overscan. 458 - ``POLY``: fit with ordinary polynomial. 459 - ``CHEB``: fit with Chebyshev polynomial. 460 - ``LEG``: fit with Legendre polynomial. 461 - ``NATURAL_SPLINE``: fit with natural spline. 462 - ``CUBIC_SPLINE``: fit with cubic spline. 463 - ``AKIMA_SPLINE``: fit with Akima spline. 466 Polynomial order or number of spline knots; ignored unless 467 ``fitType`` indicates a polynomial or spline. 468 statControl : `lsst.afw.math.StatisticsControl` 469 Statistics control object. In particular, we pay attention to numSigmaClip 470 overscanIsInt : `bool` 471 Treat the overscan region as consisting of integers, even if it's been 472 converted to float. E.g. handle ties properly. 476 result : `lsst.pipe.base.Struct` 477 Result struct with components: 479 - ``imageFit``: Value(s) removed from image (scalar or 480 `lsst.afw.image.Image`) 481 - ``overscanFit``: Value(s) removed from overscan (scalar or 482 `lsst.afw.image.Image`) 483 - ``overscanImage``: Overscan corrected overscan region 484 (`lsst.afw.image.Image`) 488 Raised if ``fitType`` is not an allowed value. 492 The ``ampMaskedImage`` and ``overscanImage`` are modified, with the fit 493 subtracted. Note that the ``overscanImage`` should not be a subimage of 494 the ``ampMaskedImage``, to avoid being subtracted twice. 496 Debug plots are available for the SPLINE fitTypes by setting the 497 `debug.display` for `name` == "lsst.ip.isr.isrFunctions". These 498 plots show the scatter plot of the overscan data (collapsed along 499 the perpendicular dimension) as a function of position on the CCD 500 (normalized between +/-1). 502 ampImage = ampMaskedImage.getImage()
506 config.fitType = fitType
510 config.numSigmaClip = collapseRej
512 config.overscanIsInt =
True 515 return overscanTask.run(ampImage, overscanImage)
519 """Apply brighter fatter correction in place for the image. 523 exposure : `lsst.afw.image.Exposure` 524 Exposure to have brighter-fatter correction applied. Modified 526 kernel : `numpy.ndarray` 527 Brighter-fatter kernel to apply. 529 Number of correction iterations to run. 531 Convergence threshold in terms of the sum of absolute 532 deviations between an iteration and the previous one. 534 If True, then the exposure values are scaled by the gain prior 536 gains : `dict` [`str`, `float`] 537 A dictionary, keyed by amplifier name, of the gains to use. 538 If gains is None, the nominal gains in the amplifier object are used. 543 Final difference between iterations achieved in correction. 545 Number of iterations used to calculate correction. 549 This correction takes a kernel that has been derived from flat 550 field images to redistribute the charge. The gradient of the 551 kernel is the deflection field due to the accumulated charge. 553 Given the original image I(x) and the kernel K(x) we can compute 554 the corrected image Ic(x) using the following equation: 556 Ic(x) = I(x) + 0.5*d/dx(I(x)*d/dx(int( dy*K(x-y)*I(y)))) 558 To evaluate the derivative term we expand it as follows: 560 0.5 * ( d/dx(I(x))*d/dx(int(dy*K(x-y)*I(y))) + I(x)*d^2/dx^2(int(dy* K(x-y)*I(y))) ) 562 Because we use the measured counts instead of the incident counts 563 we apply the correction iteratively to reconstruct the original 564 counts and the correction. We stop iterating when the summed 565 difference between the current corrected image and the one from 566 the previous iteration is below the threshold. We do not require 567 convergence because the number of iterations is too large a 568 computational cost. How we define the threshold still needs to be 569 evaluated, the current default was shown to work reasonably well 570 on a small set of images. For more information on the method see 571 DocuShare Document-19407. 573 The edges as defined by the kernel are not corrected because they 574 have spurious values due to the convolution. 576 image = exposure.getMaskedImage().getImage()
579 with
gainContext(exposure, image, applyGain, gains):
581 kLx = numpy.shape(kernel)[0]
582 kLy = numpy.shape(kernel)[1]
583 kernelImage = afwImage.ImageD(kLx, kLy)
584 kernelImage.getArray()[:, :] = kernel
585 tempImage = image.clone()
587 nanIndex = numpy.isnan(tempImage.getArray())
588 tempImage.getArray()[nanIndex] = 0.
590 outImage = afwImage.ImageF(image.getDimensions())
591 corr = numpy.zeros_like(image.getArray())
592 prev_image = numpy.zeros_like(image.getArray())
593 convCntrl = afwMath.ConvolutionControl(
False,
True, 1)
594 fixedKernel = afwMath.FixedKernel(kernelImage)
604 for iteration
in range(maxIter):
606 afwMath.convolve(outImage, tempImage, fixedKernel, convCntrl)
607 tmpArray = tempImage.getArray()
608 outArray = outImage.getArray()
610 with numpy.errstate(invalid=
"ignore", over=
"ignore"):
612 gradTmp = numpy.gradient(tmpArray[startY:endY, startX:endX])
613 gradOut = numpy.gradient(outArray[startY:endY, startX:endX])
614 first = (gradTmp[0]*gradOut[0] + gradTmp[1]*gradOut[1])[1:-1, 1:-1]
617 diffOut20 = numpy.diff(outArray, 2, 0)[startY:endY, startX + 1:endX - 1]
618 diffOut21 = numpy.diff(outArray, 2, 1)[startY + 1:endY - 1, startX:endX]
619 second = tmpArray[startY + 1:endY - 1, startX + 1:endX - 1]*(diffOut20 + diffOut21)
621 corr[startY + 1:endY - 1, startX + 1:endX - 1] = 0.5*(first + second)
623 tmpArray[:, :] = image.getArray()[:, :]
624 tmpArray[nanIndex] = 0.
625 tmpArray[startY:endY, startX:endX] += corr[startY:endY, startX:endX]
628 diff = numpy.sum(numpy.abs(prev_image - tmpArray))
632 prev_image[:, :] = tmpArray[:, :]
634 image.getArray()[startY + 1:endY - 1, startX + 1:endX - 1] += \
635 corr[startY + 1:endY - 1, startX + 1:endX - 1]
637 return diff, iteration
642 """Context manager that applies and removes gain. 646 exp : `lsst.afw.image.Exposure` 647 Exposure to apply/remove gain. 648 image : `lsst.afw.image.Image` 649 Image to apply/remove gain. 651 If True, apply and remove the amplifier gain. 652 gains : `dict` [`str`, `float`] 653 A dictionary, keyed by amplifier name, of the gains to use. 654 If gains is None, the nominal gains in the amplifier object are used. 658 exp : `lsst.afw.image.Exposure` 659 Exposure with the gain applied. 663 if gains
and apply
is True:
664 ampNames = [amp.getName()
for amp
in exp.getDetector()]
665 for ampName
in ampNames:
666 if ampName
not in gains.keys():
667 raise RuntimeError(f
"Gains provided to gain context, but no entry found for amp {ampName}")
670 ccd = exp.getDetector()
672 sim = image.Factory(image, amp.getBBox())
674 gain = gains[amp.getName()]
683 ccd = exp.getDetector()
685 sim = image.Factory(image, amp.getBBox())
687 gain = gains[amp.getName()]
694 sensorTransmission=None, atmosphereTransmission=None):
695 """Attach a TransmissionCurve to an Exposure, given separate curves for 696 different components. 700 exposure : `lsst.afw.image.Exposure` 701 Exposure object to modify by attaching the product of all given 702 ``TransmissionCurves`` in post-assembly trimmed detector coordinates. 703 Must have a valid ``Detector`` attached that matches the detector 704 associated with sensorTransmission. 705 opticsTransmission : `lsst.afw.image.TransmissionCurve` 706 A ``TransmissionCurve`` that represents the throughput of the optics, 707 to be evaluated in focal-plane coordinates. 708 filterTransmission : `lsst.afw.image.TransmissionCurve` 709 A ``TransmissionCurve`` that represents the throughput of the filter 710 itself, to be evaluated in focal-plane coordinates. 711 sensorTransmission : `lsst.afw.image.TransmissionCurve` 712 A ``TransmissionCurve`` that represents the throughput of the sensor 713 itself, to be evaluated in post-assembly trimmed detector coordinates. 714 atmosphereTransmission : `lsst.afw.image.TransmissionCurve` 715 A ``TransmissionCurve`` that represents the throughput of the 716 atmosphere, assumed to be spatially constant. 720 combined : `lsst.afw.image.TransmissionCurve` 721 The TransmissionCurve attached to the exposure. 725 All ``TransmissionCurve`` arguments are optional; if none are provided, the 726 attached ``TransmissionCurve`` will have unit transmission everywhere. 728 combined = afwImage.TransmissionCurve.makeIdentity()
729 if atmosphereTransmission
is not None:
730 combined *= atmosphereTransmission
731 if opticsTransmission
is not None:
732 combined *= opticsTransmission
733 if filterTransmission
is not None:
734 combined *= filterTransmission
735 detector = exposure.getDetector()
736 fpToPix = detector.getTransform(fromSys=camGeom.FOCAL_PLANE,
737 toSys=camGeom.PIXELS)
738 combined = combined.transformedBy(fpToPix)
739 if sensorTransmission
is not None:
740 combined *= sensorTransmission
741 exposure.getInfo().setTransmissionCurve(combined)
746 """Scale an exposure by the amplifier gains. 750 exposure : `lsst.afw.image.Exposure` 751 Exposure to process. The image is modified. 752 normalizeGains : `Bool`, optional 753 If True, then amplifiers are scaled to force the median of 754 each amplifier to equal the median of those medians. 756 ccd = exposure.getDetector()
757 ccdImage = exposure.getMaskedImage()
761 sim = ccdImage.Factory(ccdImage, amp.getBBox())
765 medians.append(numpy.median(sim.getImage().getArray()))
768 median = numpy.median(numpy.array(medians))
769 for index, amp
in enumerate(ccd):
770 sim = ccdImage.Factory(ccdImage, amp.getBBox())
771 if medians[index] != 0.0:
772 sim *= median/medians[index]
776 """Grow the saturation trails by an amount dependent on the width of the trail. 780 mask : `lsst.afw.image.Mask` 781 Mask which will have the saturated areas grown. 785 for i
in range(1, 6):
787 for i
in range(6, 8):
789 for i
in range(8, 10):
793 if extraGrowMax <= 0:
796 saturatedBit = mask.getPlaneBitMask(
"SAT")
798 xmin, ymin = mask.getBBox().getMin()
799 width = mask.getWidth()
801 thresh = afwDetection.Threshold(saturatedBit, afwDetection.Threshold.BITMASK)
802 fpList = afwDetection.FootprintSet(mask, thresh).getFootprints()
805 for s
in fp.getSpans():
806 x0, x1 = s.getX0(), s.getX1()
808 extraGrow = extraGrowDict.get(x1 - x0 + 1, extraGrowMax)
811 x0 -= xmin + extraGrow
812 x1 -= xmin - extraGrow
819 mask.array[y, x0:x1+1] |= saturatedBit
823 """Set all BAD areas of the chip to the average of the rest of the exposure 827 exposure : `lsst.afw.image.Exposure` 828 Exposure to mask. The exposure mask is modified. 829 badStatistic : `str`, optional 830 Statistic to use to generate the replacement value from the 831 image data. Allowed values are 'MEDIAN' or 'MEANCLIP'. 835 badPixelCount : scalar 836 Number of bad pixels masked. 837 badPixelValue : scalar 838 Value substituted for bad pixels. 843 Raised if `badStatistic` is not an allowed value. 845 if badStatistic ==
"MEDIAN":
846 statistic = afwMath.MEDIAN
847 elif badStatistic ==
"MEANCLIP":
848 statistic = afwMath.MEANCLIP
850 raise RuntimeError(
"Impossible method %s of bad region correction" % badStatistic)
852 mi = exposure.getMaskedImage()
854 BAD = mask.getPlaneBitMask(
"BAD")
855 INTRP = mask.getPlaneBitMask(
"INTRP")
857 sctrl = afwMath.StatisticsControl()
858 sctrl.setAndMask(BAD)
859 value = afwMath.makeStatistics(mi, statistic, sctrl).getValue()
861 maskArray = mask.getArray()
862 imageArray = mi.getImage().getArray()
863 badPixels = numpy.logical_and((maskArray & BAD) > 0, (maskArray & INTRP) == 0)
864 imageArray[:] = numpy.where(badPixels, value, imageArray)
866 return badPixels.sum(), value
def saturationCorrection(maskedImage, saturation, fwhm, growFootprints=1, interpolate=True, maskName='SAT', fallbackValue=None)
def setBadRegions(exposure, badStatistic="MEDIAN")
def interpolateFromMask(maskedImage, fwhm, growSaturatedFootprints=1, maskNameList=['SAT'], fallbackValue=None)
def interpolateDefectList(maskedImage, defectList, fwhm, fallbackValue=None)
def brighterFatterCorrection(exposure, kernel, maxIter, threshold, applyGain, gains=None)
def transposeMaskedImage(maskedImage)
def trimToMatchCalibBBox(rawMaskedImage, calibMaskedImage)
def attachTransmissionCurve(exposure, opticsTransmission=None, filterTransmission=None, sensorTransmission=None, atmosphereTransmission=None)
def flatCorrection(maskedImage, flatMaskedImage, scalingType, userScale=1.0, invert=False, trimToFit=False)
def makeThresholdMask(maskedImage, threshold, growFootprints=1, maskName='SAT')
def growMasks(mask, radius=0, maskNameList=['BAD'], maskValue="BAD")
def widenSaturationTrails(mask)
def applyGains(exposure, normalizeGains=False)
def darkCorrection(maskedImage, darkMaskedImage, expScale, darkScale, invert=False, trimToFit=False)
def biasCorrection(maskedImage, biasMaskedImage, trimToFit=False)
def gainContext(exp, image, apply, gains=None)
def illuminationCorrection(maskedImage, illumMaskedImage, illumScale, trimToFit=True)
def updateVariance(maskedImage, gain, readNoise)
def overscanCorrection(ampMaskedImage, overscanImage, fitType='MEDIAN', order=1, collapseRej=3.0, statControl=None, overscanIsInt=True)