LSSTApplications  18.1.0
LSSTDataManagementBasePackage
Functions
lsst.ip.isr.isrFunctions Namespace Reference

Functions

def createPsf (fwhm)
 
def transposeMaskedImage (maskedImage)
 
def interpolateDefectList (maskedImage, defectList, fwhm, fallbackValue=None)
 
def defectListFromFootprintList (fpList)
 
def transposeDefectList (defectList)
 
def maskPixelsFromDefectList (maskedImage, defectList, maskName='BAD')
 
def getDefectListFromMask (maskedImage, maskName)
 
def makeThresholdMask (maskedImage, threshold, growFootprints=1, maskName='SAT')
 
def interpolateFromMask (maskedImage, fwhm, growSaturatedFootprints=1, maskNameList=['SAT'], fallbackValue=None)
 
def saturationCorrection (maskedImage, saturation, fwhm, growFootprints=1, interpolate=True, maskName='SAT', fallbackValue=None)
 
def trimToMatchCalibBBox (rawMaskedImage, calibMaskedImage)
 
def biasCorrection (maskedImage, biasMaskedImage, trimToFit=False)
 
def darkCorrection (maskedImage, darkMaskedImage, expScale, darkScale, invert=False, trimToFit=False)
 
def updateVariance (maskedImage, gain, readNoise)
 
def flatCorrection (maskedImage, flatMaskedImage, scalingType, userScale=1.0, invert=False, trimToFit=False)
 
def illuminationCorrection (maskedImage, illumMaskedImage, illumScale)
 
def overscanCorrection (ampMaskedImage, overscanImage, fitType='MEDIAN', order=1, collapseRej=3.0, statControl=None, overscanIsInt=True)
 
def brighterFatterCorrection (exposure, kernel, maxIter, threshold, applyGain)
 
def gainContext (exp, image, apply)
 
def attachTransmissionCurve (exposure, opticsTransmission=None, filterTransmission=None, sensorTransmission=None, atmosphereTransmission=None)
 
def addDistortionModel (exposure, camera)
 Update the WCS in exposure with a distortion model based on camera geometry. More...
 
def applyGains (exposure, normalizeGains=False)
 
def widenSaturationTrails (mask)
 
def setBadRegions (exposure, badStatistic="MEDIAN")
 

Function Documentation

◆ addDistortionModel()

def lsst.ip.isr.isrFunctions.addDistortionModel (   exposure,
  camera 
)

Update the WCS in exposure with a distortion model based on camera geometry.

Parameters

exposure : lsst.afw.image.Exposure Exposure to process. Must contain a Detector and WCS. The exposure is modified. camera : lsst.afw.cameraGeom.Camera Camera geometry.

Raises

RuntimeError Raised if exposure is lacking a Detector or WCS, or if camera is None.

Notes

Add a model for optical distortion based on geometry found in camera and the exposure's detector. The raw input exposure is assumed have a TAN WCS that has no compensation for optical distortion. Two other possibilities are:

  • The raw input exposure already has a model for optical distortion, as is the case for raw DECam data. In that case you should set config.doAddDistortionModel False.
  • The raw input exposure has a model for distortion, but it has known deficiencies severe enough to be worth fixing (e.g. because they cause problems for fitting a better WCS). In that case you should override this method with a version suitable for your raw data.

Definition at line 919 of file isrFunctions.py.

919 def addDistortionModel(exposure, camera):
920  """!Update the WCS in exposure with a distortion model based on camera
921  geometry.
922 
923  Parameters
924  ----------
925  exposure : `lsst.afw.image.Exposure`
926  Exposure to process. Must contain a Detector and WCS. The
927  exposure is modified.
928  camera : `lsst.afw.cameraGeom.Camera`
929  Camera geometry.
930 
931  Raises
932  ------
933  RuntimeError
934  Raised if ``exposure`` is lacking a Detector or WCS, or if
935  ``camera`` is None.
936  Notes
937  -----
938  Add a model for optical distortion based on geometry found in ``camera``
939  and the ``exposure``'s detector. The raw input exposure is assumed
940  have a TAN WCS that has no compensation for optical distortion.
941  Two other possibilities are:
942  - The raw input exposure already has a model for optical distortion,
943  as is the case for raw DECam data.
944  In that case you should set config.doAddDistortionModel False.
945  - The raw input exposure has a model for distortion, but it has known
946  deficiencies severe enough to be worth fixing (e.g. because they
947  cause problems for fitting a better WCS). In that case you should
948  override this method with a version suitable for your raw data.
949 
950  """
951  wcs = exposure.getWcs()
952  if wcs is None:
953  raise RuntimeError("exposure has no WCS")
954  if camera is None:
955  raise RuntimeError("camera is None")
956  detector = exposure.getDetector()
957  if detector is None:
958  raise RuntimeError("exposure has no Detector")
959  pixelToFocalPlane = detector.getTransform(camGeom.PIXELS, camGeom.FOCAL_PLANE)
960  focalPlaneToFieldAngle = camera.getTransformMap().getTransform(camGeom.FOCAL_PLANE,
961  camGeom.FIELD_ANGLE)
962  distortedWcs = makeDistortedTanWcs(wcs, pixelToFocalPlane, focalPlaneToFieldAngle)
963  exposure.setWcs(distortedWcs)
964 
965 
def addDistortionModel(exposure, camera)
Update the WCS in exposure with a distortion model based on camera geometry.
def makeDistortedTanWcs(tanWcs, pixelToFocalPlane, focalPlaneToFieldAngle)

◆ applyGains()

def lsst.ip.isr.isrFunctions.applyGains (   exposure,
  normalizeGains = False 
)
Scale an exposure by the amplifier gains.

Parameters
----------
exposure : `lsst.afw.image.Exposure`
    Exposure to process.  The image is modified.
normalizeGains : `Bool`, optional
    If True, then amplifiers are scaled to force the median of
    each amplifier to equal the median of those medians.

Definition at line 966 of file isrFunctions.py.

966 def applyGains(exposure, normalizeGains=False):
967  """Scale an exposure by the amplifier gains.
968 
969  Parameters
970  ----------
971  exposure : `lsst.afw.image.Exposure`
972  Exposure to process. The image is modified.
973  normalizeGains : `Bool`, optional
974  If True, then amplifiers are scaled to force the median of
975  each amplifier to equal the median of those medians.
976  """
977  ccd = exposure.getDetector()
978  ccdImage = exposure.getMaskedImage()
979 
980  medians = []
981  for amp in ccd:
982  sim = ccdImage.Factory(ccdImage, amp.getBBox())
983  sim *= amp.getGain()
984 
985  if normalizeGains:
986  medians.append(numpy.median(sim.getImage().getArray()))
987 
988  if normalizeGains:
989  median = numpy.median(numpy.array(medians))
990  for index, amp in enumerate(ccd):
991  sim = ccdImage.Factory(ccdImage, amp.getBBox())
992  if medians[index] != 0.0:
993  sim *= median/medians[index]
994 
995 
def applyGains(exposure, normalizeGains=False)

◆ attachTransmissionCurve()

def lsst.ip.isr.isrFunctions.attachTransmissionCurve (   exposure,
  opticsTransmission = None,
  filterTransmission = None,
  sensorTransmission = None,
  atmosphereTransmission = None 
)
Attach a TransmissionCurve to an Exposure, given separate curves for
different components.

Parameters
----------
exposure : `lsst.afw.image.Exposure`
    Exposure object to modify by attaching the product of all given
    ``TransmissionCurves`` in post-assembly trimmed detector coordinates.
    Must have a valid ``Detector`` attached that matches the detector
    associated with sensorTransmission.
opticsTransmission : `lsst.afw.image.TransmissionCurve`
    A ``TransmissionCurve`` that represents the throughput of the optics,
    to be evaluated in focal-plane coordinates.
filterTransmission : `lsst.afw.image.TransmissionCurve`
    A ``TransmissionCurve`` that represents the throughput of the filter
    itself, to be evaluated in focal-plane coordinates.
sensorTransmission : `lsst.afw.image.TransmissionCurve`
    A ``TransmissionCurve`` that represents the throughput of the sensor
    itself, to be evaluated in post-assembly trimmed detector coordinates.
atmosphereTransmission : `lsst.afw.image.TransmissionCurve`
    A ``TransmissionCurve`` that represents the throughput of the
    atmosphere, assumed to be spatially constant.

Returns
-------
combined : `lsst.afw.image.TransmissionCurve`
    The TransmissionCurve attached to the exposure.

Notes
-----
All ``TransmissionCurve`` arguments are optional; if none are provided, the
attached ``TransmissionCurve`` will have unit transmission everywhere.

Definition at line 868 of file isrFunctions.py.

868  sensorTransmission=None, atmosphereTransmission=None):
869  """Attach a TransmissionCurve to an Exposure, given separate curves for
870  different components.
871 
872  Parameters
873  ----------
874  exposure : `lsst.afw.image.Exposure`
875  Exposure object to modify by attaching the product of all given
876  ``TransmissionCurves`` in post-assembly trimmed detector coordinates.
877  Must have a valid ``Detector`` attached that matches the detector
878  associated with sensorTransmission.
879  opticsTransmission : `lsst.afw.image.TransmissionCurve`
880  A ``TransmissionCurve`` that represents the throughput of the optics,
881  to be evaluated in focal-plane coordinates.
882  filterTransmission : `lsst.afw.image.TransmissionCurve`
883  A ``TransmissionCurve`` that represents the throughput of the filter
884  itself, to be evaluated in focal-plane coordinates.
885  sensorTransmission : `lsst.afw.image.TransmissionCurve`
886  A ``TransmissionCurve`` that represents the throughput of the sensor
887  itself, to be evaluated in post-assembly trimmed detector coordinates.
888  atmosphereTransmission : `lsst.afw.image.TransmissionCurve`
889  A ``TransmissionCurve`` that represents the throughput of the
890  atmosphere, assumed to be spatially constant.
891 
892  Returns
893  -------
894  combined : `lsst.afw.image.TransmissionCurve`
895  The TransmissionCurve attached to the exposure.
896 
897  Notes
898  -----
899  All ``TransmissionCurve`` arguments are optional; if none are provided, the
900  attached ``TransmissionCurve`` will have unit transmission everywhere.
901  """
902  combined = afwImage.TransmissionCurve.makeIdentity()
903  if atmosphereTransmission is not None:
904  combined *= atmosphereTransmission
905  if opticsTransmission is not None:
906  combined *= opticsTransmission
907  if filterTransmission is not None:
908  combined *= filterTransmission
909  detector = exposure.getDetector()
910  fpToPix = detector.getTransform(fromSys=camGeom.FOCAL_PLANE,
911  toSys=camGeom.PIXELS)
912  combined = combined.transformedBy(fpToPix)
913  if sensorTransmission is not None:
914  combined *= sensorTransmission
915  exposure.getInfo().setTransmissionCurve(combined)
916  return combined
917 
918 

◆ biasCorrection()

def lsst.ip.isr.isrFunctions.biasCorrection (   maskedImage,
  biasMaskedImage,
  trimToFit = False 
)
Apply bias correction in place.

Parameters
----------
maskedImage : `lsst.afw.image.MaskedImage`
   Image to process.  The image is modified by this method.
biasMaskedImage : `lsst.afw.image.MaskedImage`
    Bias image of the same size as ``maskedImage``
trimToFit : `Bool`, optional
    If True, raw data is symmetrically trimmed to match
    calibration size.

Raises
------
RuntimeError
    Raised if ``maskedImage`` and ``biasMaskedImage`` do not have
    the same size.

Definition at line 327 of file isrFunctions.py.

327 def biasCorrection(maskedImage, biasMaskedImage, trimToFit=False):
328  """Apply bias correction in place.
329 
330  Parameters
331  ----------
332  maskedImage : `lsst.afw.image.MaskedImage`
333  Image to process. The image is modified by this method.
334  biasMaskedImage : `lsst.afw.image.MaskedImage`
335  Bias image of the same size as ``maskedImage``
336  trimToFit : `Bool`, optional
337  If True, raw data is symmetrically trimmed to match
338  calibration size.
339 
340  Raises
341  ------
342  RuntimeError
343  Raised if ``maskedImage`` and ``biasMaskedImage`` do not have
344  the same size.
345 
346  """
347  if trimToFit:
348  maskedImage = trimToMatchCalibBBox(maskedImage, biasMaskedImage)
349 
350  if maskedImage.getBBox(afwImage.LOCAL) != biasMaskedImage.getBBox(afwImage.LOCAL):
351  raise RuntimeError("maskedImage bbox %s != biasMaskedImage bbox %s" %
352  (maskedImage.getBBox(afwImage.LOCAL), biasMaskedImage.getBBox(afwImage.LOCAL)))
353  maskedImage -= biasMaskedImage
354 
355 
def trimToMatchCalibBBox(rawMaskedImage, calibMaskedImage)
def biasCorrection(maskedImage, biasMaskedImage, trimToFit=False)

◆ brighterFatterCorrection()

def lsst.ip.isr.isrFunctions.brighterFatterCorrection (   exposure,
  kernel,
  maxIter,
  threshold,
  applyGain 
)
Apply brighter fatter correction in place for the image.

Parameters
----------
exposure : `lsst.afw.image.Exposure`
    Exposure to have brighter-fatter correction applied.  Modified
    by this method.
kernel : `numpy.ndarray`
    Brighter-fatter kernel to apply.
maxIter : scalar
    Number of correction iterations to run.
threshold : scalar
    Convergence threshold in terms of the sum of absolute
    deviations between an iteration and the previous one.
applyGain : `Bool`
    If True, then the exposure values are scaled by the gain prior
    to correction.

Notes
-----
This correction takes a kernel that has been derived from flat
field images to redistribute the charge.  The gradient of the
kernel is the deflection field due to the accumulated charge.

Given the original image I(x) and the kernel K(x) we can compute
the corrected image Ic(x) using the following equation:

Ic(x) = I(x) + 0.5*d/dx(I(x)*d/dx(int( dy*K(x-y)*I(y))))

To evaluate the derivative term we expand it as follows:

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))) )

Because we use the measured counts instead of the incident counts
we apply the correction iteratively to reconstruct the original
counts and the correction.  We stop iterating when the summed
difference between the current corrected image and the one from
the previous iteration is below the threshold.  We do not require
convergence because the number of iterations is too large a
computational cost.  How we define the threshold still needs to be
evaluated, the current default was shown to work reasonably well
on a small set of images.  For more information on the method see
DocuShare Document-19407.

The edges as defined by the kernel are not corrected because they
have spurious values due to the convolution.

Definition at line 718 of file isrFunctions.py.

718 def brighterFatterCorrection(exposure, kernel, maxIter, threshold, applyGain):
719  """Apply brighter fatter correction in place for the image.
720 
721  Parameters
722  ----------
723  exposure : `lsst.afw.image.Exposure`
724  Exposure to have brighter-fatter correction applied. Modified
725  by this method.
726  kernel : `numpy.ndarray`
727  Brighter-fatter kernel to apply.
728  maxIter : scalar
729  Number of correction iterations to run.
730  threshold : scalar
731  Convergence threshold in terms of the sum of absolute
732  deviations between an iteration and the previous one.
733  applyGain : `Bool`
734  If True, then the exposure values are scaled by the gain prior
735  to correction.
736 
737  Notes
738  -----
739  This correction takes a kernel that has been derived from flat
740  field images to redistribute the charge. The gradient of the
741  kernel is the deflection field due to the accumulated charge.
742 
743  Given the original image I(x) and the kernel K(x) we can compute
744  the corrected image Ic(x) using the following equation:
745 
746  Ic(x) = I(x) + 0.5*d/dx(I(x)*d/dx(int( dy*K(x-y)*I(y))))
747 
748  To evaluate the derivative term we expand it as follows:
749 
750  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))) )
751 
752  Because we use the measured counts instead of the incident counts
753  we apply the correction iteratively to reconstruct the original
754  counts and the correction. We stop iterating when the summed
755  difference between the current corrected image and the one from
756  the previous iteration is below the threshold. We do not require
757  convergence because the number of iterations is too large a
758  computational cost. How we define the threshold still needs to be
759  evaluated, the current default was shown to work reasonably well
760  on a small set of images. For more information on the method see
761  DocuShare Document-19407.
762 
763  The edges as defined by the kernel are not corrected because they
764  have spurious values due to the convolution.
765  """
766  image = exposure.getMaskedImage().getImage()
767 
768  # The image needs to be units of electrons/holes
769  with gainContext(exposure, image, applyGain):
770 
771  kLx = numpy.shape(kernel)[0]
772  kLy = numpy.shape(kernel)[1]
773  kernelImage = afwImage.ImageD(kLx, kLy)
774  kernelImage.getArray()[:, :] = kernel
775  tempImage = image.clone()
776 
777  nanIndex = numpy.isnan(tempImage.getArray())
778  tempImage.getArray()[nanIndex] = 0.
779 
780  outImage = afwImage.ImageF(image.getDimensions())
781  corr = numpy.zeros_like(image.getArray())
782  prev_image = numpy.zeros_like(image.getArray())
783  convCntrl = afwMath.ConvolutionControl(False, True, 1)
784  fixedKernel = afwMath.FixedKernel(kernelImage)
785 
786  # Define boundary by convolution region. The region that the correction will be
787  # calculated for is one fewer in each dimension because of the second derivative terms.
788  # NOTE: these need to use integer math, as we're using start:end as numpy index ranges.
789  startX = kLx//2
790  endX = -kLx//2
791  startY = kLy//2
792  endY = -kLy//2
793 
794  for iteration in range(maxIter):
795 
796  afwMath.convolve(outImage, tempImage, fixedKernel, convCntrl)
797  tmpArray = tempImage.getArray()
798  outArray = outImage.getArray()
799 
800  with numpy.errstate(invalid="ignore", over="ignore"):
801  # First derivative term
802  gradTmp = numpy.gradient(tmpArray[startY:endY, startX:endX])
803  gradOut = numpy.gradient(outArray[startY:endY, startX:endX])
804  first = (gradTmp[0]*gradOut[0] + gradTmp[1]*gradOut[1])[1:-1, 1:-1]
805 
806  # Second derivative term
807  diffOut20 = numpy.diff(outArray, 2, 0)[startY:endY, startX + 1:endX - 1]
808  diffOut21 = numpy.diff(outArray, 2, 1)[startY + 1:endY - 1, startX:endX]
809  second = tmpArray[startY + 1:endY - 1, startX + 1:endX - 1]*(diffOut20 + diffOut21)
810 
811  corr[startY + 1:endY - 1, startX + 1:endX - 1] = 0.5*(first + second)
812 
813  tmpArray[:, :] = image.getArray()[:, :]
814  tmpArray[nanIndex] = 0.
815  tmpArray[startY:endY, startX:endX] += corr[startY:endY, startX:endX]
816 
817  if iteration > 0:
818  diff = numpy.sum(numpy.abs(prev_image - tmpArray))
819 
820  if diff < threshold:
821  break
822  prev_image[:, :] = tmpArray[:, :]
823 
824  # if iteration == maxIter - 1:
825  # self.log.warn("Brighter fatter correction did not converge,
826  # final difference %f" % diff)
827 
828  # self.log.info("Finished brighter fatter in %d iterations" % (iteration + 1))
829  image.getArray()[startY + 1:endY - 1, startX + 1:endX - 1] += \
830  corr[startY + 1:endY - 1, startX + 1:endX - 1]
831 
832 
833 @contextmanager
def brighterFatterCorrection(exposure, kernel, maxIter, threshold, applyGain)
Parameters to control convolution.
Definition: ConvolveImage.h:50
void convolve(OutImageT &convolvedImage, InImageT const &inImage, KernelT const &kernel, bool doNormalize, bool doCopyEdge=false)
Old, deprecated version of convolve.
def gainContext(exp, image, apply)
A kernel created from an Image.
Definition: Kernel.h:509

◆ createPsf()

def lsst.ip.isr.isrFunctions.createPsf (   fwhm)
Make a double Gaussian PSF.

Parameters
----------
fwhm : scalar
    FWHM of double Gaussian smoothing kernel.

Returns
-------
psf : `lsst.meas.algorithms.DoubleGaussianPsf`
    The created smoothing kernel.

Definition at line 41 of file isrFunctions.py.

41 def createPsf(fwhm):
42  """Make a double Gaussian PSF.
43 
44  Parameters
45  ----------
46  fwhm : scalar
47  FWHM of double Gaussian smoothing kernel.
48 
49  Returns
50  -------
51  psf : `lsst.meas.algorithms.DoubleGaussianPsf`
52  The created smoothing kernel.
53  """
54  ksize = 4*int(fwhm) + 1
55  return measAlg.DoubleGaussianPsf(ksize, ksize, fwhm/(2*math.sqrt(2*math.log(2))))
56 
57 

◆ darkCorrection()

def lsst.ip.isr.isrFunctions.darkCorrection (   maskedImage,
  darkMaskedImage,
  expScale,
  darkScale,
  invert = False,
  trimToFit = False 
)
Apply dark correction in place.

Parameters
----------
maskedImage : `lsst.afw.image.MaskedImage`
   Image to process.  The image is modified by this method.
darkMaskedImage : `lsst.afw.image.MaskedImage`
    Dark image of the same size as ``maskedImage``.
expScale : scalar
    Dark exposure time for ``maskedImage``.
darkScale : scalar
    Dark exposure time for ``darkMaskedImage``.
invert : `Bool`, optional
    If True, re-add the dark to an already corrected image.
trimToFit : `Bool`, optional
    If True, raw data is symmetrically trimmed to match
    calibration size.

Raises
------
RuntimeError
    Raised if ``maskedImage`` and ``darkMaskedImage`` do not have
    the same size.

Notes
-----
The dark correction is applied by calculating:
    maskedImage -= dark * expScaling / darkScaling

Definition at line 356 of file isrFunctions.py.

356 def darkCorrection(maskedImage, darkMaskedImage, expScale, darkScale, invert=False, trimToFit=False):
357  """Apply dark correction in place.
358 
359  Parameters
360  ----------
361  maskedImage : `lsst.afw.image.MaskedImage`
362  Image to process. The image is modified by this method.
363  darkMaskedImage : `lsst.afw.image.MaskedImage`
364  Dark image of the same size as ``maskedImage``.
365  expScale : scalar
366  Dark exposure time for ``maskedImage``.
367  darkScale : scalar
368  Dark exposure time for ``darkMaskedImage``.
369  invert : `Bool`, optional
370  If True, re-add the dark to an already corrected image.
371  trimToFit : `Bool`, optional
372  If True, raw data is symmetrically trimmed to match
373  calibration size.
374 
375  Raises
376  ------
377  RuntimeError
378  Raised if ``maskedImage`` and ``darkMaskedImage`` do not have
379  the same size.
380 
381  Notes
382  -----
383  The dark correction is applied by calculating:
384  maskedImage -= dark * expScaling / darkScaling
385  """
386  if trimToFit:
387  maskedImage = trimToMatchCalibBBox(maskedImage, darkMaskedImage)
388 
389  if maskedImage.getBBox(afwImage.LOCAL) != darkMaskedImage.getBBox(afwImage.LOCAL):
390  raise RuntimeError("maskedImage bbox %s != darkMaskedImage bbox %s" %
391  (maskedImage.getBBox(afwImage.LOCAL), darkMaskedImage.getBBox(afwImage.LOCAL)))
392 
393  scale = expScale / darkScale
394  if not invert:
395  maskedImage.scaledMinus(scale, darkMaskedImage)
396  else:
397  maskedImage.scaledPlus(scale, darkMaskedImage)
398 
399 
def trimToMatchCalibBBox(rawMaskedImage, calibMaskedImage)
def darkCorrection(maskedImage, darkMaskedImage, expScale, darkScale, invert=False, trimToFit=False)

◆ defectListFromFootprintList()

def lsst.ip.isr.isrFunctions.defectListFromFootprintList (   fpList)
Compute a defect list from a footprint list, optionally growing the footprints.

Parameters
----------
fpList : `list` of `lsst.afw.detection.Footprint`
    Footprint list to process.

Returns
-------
defectList : `lsst.meas.algorithms.Defects`
    List of defects.

Definition at line 104 of file isrFunctions.py.

104 def defectListFromFootprintList(fpList):
105  """Compute a defect list from a footprint list, optionally growing the footprints.
106 
107  Parameters
108  ----------
109  fpList : `list` of `lsst.afw.detection.Footprint`
110  Footprint list to process.
111 
112  Returns
113  -------
114  defectList : `lsst.meas.algorithms.Defects`
115  List of defects.
116  """
117  return measAlg.Defects.fromFootprintList(fpList)
118 
119 @deprecated(reason="Replaced by Defects.transpose() (will be removed after v18)",
120  category=FutureWarning)
def defectListFromFootprintList(fpList)

◆ flatCorrection()

def lsst.ip.isr.isrFunctions.flatCorrection (   maskedImage,
  flatMaskedImage,
  scalingType,
  userScale = 1.0,
  invert = False,
  trimToFit = False 
)
Apply flat correction in place.

Parameters
----------
maskedImage : `lsst.afw.image.MaskedImage`
    Image to process.  The image is modified.
flatMaskedImage : `lsst.afw.image.MaskedImage`
    Flat image of the same size as ``maskedImage``
scalingType : str
    Flat scale computation method.  Allowed values are 'MEAN',
    'MEDIAN', or 'USER'.
userScale : scalar, optional
    Scale to use if ``scalingType``='USER'.
invert : `Bool`, optional
    If True, unflatten an already flattened image.
trimToFit : `Bool`, optional
    If True, raw data is symmetrically trimmed to match
    calibration size.

Raises
------
RuntimeError
    Raised if ``maskedImage`` and ``flatMaskedImage`` do not have
    the same size or if ``scalingType`` is not an allowed value.

Definition at line 418 of file isrFunctions.py.

418 def flatCorrection(maskedImage, flatMaskedImage, scalingType, userScale=1.0, invert=False, trimToFit=False):
419  """Apply flat correction in place.
420 
421  Parameters
422  ----------
423  maskedImage : `lsst.afw.image.MaskedImage`
424  Image to process. The image is modified.
425  flatMaskedImage : `lsst.afw.image.MaskedImage`
426  Flat image of the same size as ``maskedImage``
427  scalingType : str
428  Flat scale computation method. Allowed values are 'MEAN',
429  'MEDIAN', or 'USER'.
430  userScale : scalar, optional
431  Scale to use if ``scalingType``='USER'.
432  invert : `Bool`, optional
433  If True, unflatten an already flattened image.
434  trimToFit : `Bool`, optional
435  If True, raw data is symmetrically trimmed to match
436  calibration size.
437 
438  Raises
439  ------
440  RuntimeError
441  Raised if ``maskedImage`` and ``flatMaskedImage`` do not have
442  the same size or if ``scalingType`` is not an allowed value.
443  """
444  if trimToFit:
445  maskedImage = trimToMatchCalibBBox(maskedImage, flatMaskedImage)
446 
447  if maskedImage.getBBox(afwImage.LOCAL) != flatMaskedImage.getBBox(afwImage.LOCAL):
448  raise RuntimeError("maskedImage bbox %s != flatMaskedImage bbox %s" %
449  (maskedImage.getBBox(afwImage.LOCAL), flatMaskedImage.getBBox(afwImage.LOCAL)))
450 
451  # Figure out scale from the data
452  # Ideally the flats are normalized by the calibration product pipeline, but this allows some flexibility
453  # in the case that the flat is created by some other mechanism.
454  if scalingType in ('MEAN', 'MEDIAN'):
455  scalingType = afwMath.stringToStatisticsProperty(scalingType)
456  flatScale = afwMath.makeStatistics(flatMaskedImage.image, scalingType).getValue()
457  elif scalingType == 'USER':
458  flatScale = userScale
459  else:
460  raise RuntimeError('%s : %s not implemented' % ("flatCorrection", scalingType))
461 
462  if not invert:
463  maskedImage.scaledDivides(1.0/flatScale, flatMaskedImage)
464  else:
465  maskedImage.scaledMultiplies(1.0/flatScale, flatMaskedImage)
466 
467 
def trimToMatchCalibBBox(rawMaskedImage, calibMaskedImage)
Statistics makeStatistics(lsst::afw::math::MaskedVector< EntryT > const &mv, std::vector< WeightPixel > const &vweights, int const flags, StatisticsControl const &sctrl=StatisticsControl())
The makeStatistics() overload to handle lsst::afw::math::MaskedVector<>
Definition: Statistics.h:520
Property stringToStatisticsProperty(std::string const property)
Conversion function to switch a string to a Property (see Statistics.h)
Definition: Statistics.cc:747
def flatCorrection(maskedImage, flatMaskedImage, scalingType, userScale=1.0, invert=False, trimToFit=False)

◆ gainContext()

def lsst.ip.isr.isrFunctions.gainContext (   exp,
  image,
  apply 
)
Context manager that applies and removes gain.

Parameters
----------
exp : `lsst.afw.image.Exposure`
    Exposure to apply/remove gain.
image : `lsst.afw.image.Image`
    Image to apply/remove gain.
apply : `Bool`
    If True, apply and remove the amplifier gain.

Yields
------
exp : `lsst.afw.image.Exposure`
    Exposure with the gain applied.

Definition at line 834 of file isrFunctions.py.

834 def gainContext(exp, image, apply):
835  """Context manager that applies and removes gain.
836 
837  Parameters
838  ----------
839  exp : `lsst.afw.image.Exposure`
840  Exposure to apply/remove gain.
841  image : `lsst.afw.image.Image`
842  Image to apply/remove gain.
843  apply : `Bool`
844  If True, apply and remove the amplifier gain.
845 
846  Yields
847  ------
848  exp : `lsst.afw.image.Exposure`
849  Exposure with the gain applied.
850  """
851  if apply:
852  ccd = exp.getDetector()
853  for amp in ccd:
854  sim = image.Factory(image, amp.getBBox())
855  sim *= amp.getGain()
856 
857  try:
858  yield exp
859  finally:
860  if apply:
861  ccd = exp.getDetector()
862  for amp in ccd:
863  sim = image.Factory(image, amp.getBBox())
864  sim /= amp.getGain()
865 
866 
def gainContext(exp, image, apply)

◆ getDefectListFromMask()

def lsst.ip.isr.isrFunctions.getDefectListFromMask (   maskedImage,
  maskName 
)
Compute a defect list from a specified mask plane.

Parameters
----------
maskedImage : `lsst.afw.image.MaskedImage`
    Image to process.
maskName : `str` or `list`
    Mask plane name, or list of names to convert.

Returns
-------
defectList : `lsst.meas.algorithms.Defects`
    Defect list constructed from masked pixels.

Definition at line 158 of file isrFunctions.py.

158 def getDefectListFromMask(maskedImage, maskName):
159  """Compute a defect list from a specified mask plane.
160 
161  Parameters
162  ----------
163  maskedImage : `lsst.afw.image.MaskedImage`
164  Image to process.
165  maskName : `str` or `list`
166  Mask plane name, or list of names to convert.
167 
168  Returns
169  -------
170  defectList : `lsst.meas.algorithms.Defects`
171  Defect list constructed from masked pixels.
172  """
173  return measAlg.Defects.fromMask(maskedImage, maskName)
174 
175 
def getDefectListFromMask(maskedImage, maskName)

◆ illuminationCorrection()

def lsst.ip.isr.isrFunctions.illuminationCorrection (   maskedImage,
  illumMaskedImage,
  illumScale 
)
Apply illumination correction in place.

Parameters
----------
maskedImage : `lsst.afw.image.MaskedImage`
    Image to process.  The image is modified.
illumMaskedImage : `lsst.afw.image.MaskedImage`
    Illumination correction image of the same size as ``maskedImage``.
illumScale : scalar
    Scale factor for the illumination correction.

Raises
------
RuntimeError
    Raised if ``maskedImage`` and ``illumMaskedImage`` do not have
    the same size.

Definition at line 468 of file isrFunctions.py.

468 def illuminationCorrection(maskedImage, illumMaskedImage, illumScale):
469  """Apply illumination correction in place.
470 
471  Parameters
472  ----------
473  maskedImage : `lsst.afw.image.MaskedImage`
474  Image to process. The image is modified.
475  illumMaskedImage : `lsst.afw.image.MaskedImage`
476  Illumination correction image of the same size as ``maskedImage``.
477  illumScale : scalar
478  Scale factor for the illumination correction.
479 
480  Raises
481  ------
482  RuntimeError
483  Raised if ``maskedImage`` and ``illumMaskedImage`` do not have
484  the same size.
485  """
486  if maskedImage.getBBox(afwImage.LOCAL) != illumMaskedImage.getBBox(afwImage.LOCAL):
487  raise RuntimeError("maskedImage bbox %s != illumMaskedImage bbox %s" %
488  (maskedImage.getBBox(afwImage.LOCAL), illumMaskedImage.getBBox(afwImage.LOCAL)))
489 
490  maskedImage.scaledDivides(1./illumScale, illumMaskedImage)
491 
492 
def illuminationCorrection(maskedImage, illumMaskedImage, illumScale)

◆ interpolateDefectList()

def lsst.ip.isr.isrFunctions.interpolateDefectList (   maskedImage,
  defectList,
  fwhm,
  fallbackValue = None 
)
Interpolate over defects specified in a defect list.

Parameters
----------
maskedImage : `lsst.afw.image.MaskedImage`
    Image to process.
defectList : `lsst.meas.algorithms.Defects`
    List of defects to interpolate over.
fwhm : scalar
    FWHM of double Gaussian smoothing kernel.
fallbackValue : scalar, optional
    Fallback value if an interpolated value cannot be determined.
    If None, then the clipped mean of the image is used.

Definition at line 78 of file isrFunctions.py.

78 def interpolateDefectList(maskedImage, defectList, fwhm, fallbackValue=None):
79  """Interpolate over defects specified in a defect list.
80 
81  Parameters
82  ----------
83  maskedImage : `lsst.afw.image.MaskedImage`
84  Image to process.
85  defectList : `lsst.meas.algorithms.Defects`
86  List of defects to interpolate over.
87  fwhm : scalar
88  FWHM of double Gaussian smoothing kernel.
89  fallbackValue : scalar, optional
90  Fallback value if an interpolated value cannot be determined.
91  If None, then the clipped mean of the image is used.
92  """
93  psf = createPsf(fwhm)
94  if fallbackValue is None:
95  fallbackValue = afwMath.makeStatistics(maskedImage.getImage(), afwMath.MEANCLIP).getValue()
96  if 'INTRP' not in maskedImage.getMask().getMaskPlaneDict():
97  maskedImage.getMask().addMaskPlane('INTRP')
98  measAlg.interpolateOverDefects(maskedImage, psf, defectList, fallbackValue, True)
99  return maskedImage
100 
101 
102 @deprecated(reason="Replaced by Defects.fromFootPrintList() (will be removed after v18)",
103  category=FutureWarning)
def interpolateDefectList(maskedImage, defectList, fwhm, fallbackValue=None)
Definition: isrFunctions.py:78
Statistics makeStatistics(lsst::afw::math::MaskedVector< EntryT > const &mv, std::vector< WeightPixel > const &vweights, int const flags, StatisticsControl const &sctrl=StatisticsControl())
The makeStatistics() overload to handle lsst::afw::math::MaskedVector<>
Definition: Statistics.h:520

◆ interpolateFromMask()

def lsst.ip.isr.isrFunctions.interpolateFromMask (   maskedImage,
  fwhm,
  growSaturatedFootprints = 1,
  maskNameList = ['SAT'],
  fallbackValue = None 
)
Interpolate over defects identified by a particular set of mask planes.

Parameters
----------
maskedImage : `lsst.afw.image.MaskedImage`
    Image to process.
fwhm : scalar
    FWHM of double Gaussian smoothing kernel.
growSaturatedFootprints : scalar, optional
    Number of pixels to grow footprints for saturated pixels.
maskNameList : `List` of `str`, optional
    Mask plane name.
fallbackValue : scalar, optional
    Value of last resort for interpolation.

Definition at line 212 of file isrFunctions.py.

212  maskNameList=['SAT'], fallbackValue=None):
213  """Interpolate over defects identified by a particular set of mask planes.
214 
215  Parameters
216  ----------
217  maskedImage : `lsst.afw.image.MaskedImage`
218  Image to process.
219  fwhm : scalar
220  FWHM of double Gaussian smoothing kernel.
221  growSaturatedFootprints : scalar, optional
222  Number of pixels to grow footprints for saturated pixels.
223  maskNameList : `List` of `str`, optional
224  Mask plane name.
225  fallbackValue : scalar, optional
226  Value of last resort for interpolation.
227  """
228  mask = maskedImage.getMask()
229 
230  if growSaturatedFootprints > 0 and "SAT" in maskNameList:
231  thresh = afwDetection.Threshold(mask.getPlaneBitMask("SAT"), afwDetection.Threshold.BITMASK)
232  fpSet = afwDetection.FootprintSet(mask, thresh)
233  # If we are interpolating over an area larger than the original masked region, we need
234  # to expand the original mask bit to the full area to explain why we interpolated there.
235  fpSet = afwDetection.FootprintSet(fpSet, rGrow=growSaturatedFootprints, isotropic=False)
236  fpSet.setMask(mask, "SAT")
237 
238  thresh = afwDetection.Threshold(mask.getPlaneBitMask(maskNameList), afwDetection.Threshold.BITMASK)
239  fpSet = afwDetection.FootprintSet(mask, thresh)
240  defectList = measAlg.Defects.fromFootprintList(fpSet.getFootprints())
241 
242  interpolateDefectList(maskedImage, defectList, fwhm, fallbackValue=fallbackValue)
243 
244  return maskedImage
245 
246 
def interpolateDefectList(maskedImage, defectList, fwhm, fallbackValue=None)
Definition: isrFunctions.py:78

◆ makeThresholdMask()

def lsst.ip.isr.isrFunctions.makeThresholdMask (   maskedImage,
  threshold,
  growFootprints = 1,
  maskName = 'SAT' 
)
Mask pixels based on threshold detection.

Parameters
----------
maskedImage : `lsst.afw.image.MaskedImage`
    Image to process.  Only the mask plane is updated.
threshold : scalar
    Detection threshold.
growFootprints : scalar, optional
    Number of pixels to grow footprints of detected regions.
maskName : str, optional
    Mask plane name, or list of names to convert

Returns
-------
defectList : `lsst.meas.algorithms.Defects`
    Defect list constructed from pixels above the threshold.

Definition at line 176 of file isrFunctions.py.

176 def makeThresholdMask(maskedImage, threshold, growFootprints=1, maskName='SAT'):
177  """Mask pixels based on threshold detection.
178 
179  Parameters
180  ----------
181  maskedImage : `lsst.afw.image.MaskedImage`
182  Image to process. Only the mask plane is updated.
183  threshold : scalar
184  Detection threshold.
185  growFootprints : scalar, optional
186  Number of pixels to grow footprints of detected regions.
187  maskName : str, optional
188  Mask plane name, or list of names to convert
189 
190  Returns
191  -------
192  defectList : `lsst.meas.algorithms.Defects`
193  Defect list constructed from pixels above the threshold.
194  """
195  # find saturated regions
196  thresh = afwDetection.Threshold(threshold)
197  fs = afwDetection.FootprintSet(maskedImage, thresh)
198 
199  if growFootprints > 0:
200  fs = afwDetection.FootprintSet(fs, rGrow=growFootprints, isotropic=False)
201  fpList = fs.getFootprints()
202 
203  # set mask
204  mask = maskedImage.getMask()
205  bitmask = mask.getPlaneBitMask(maskName)
206  afwDetection.setMaskFromFootprintList(mask, fpList, bitmask)
207 
208  return measAlg.Defects.fromFootprintList(fpList)
209 
210 
def makeThresholdMask(maskedImage, threshold, growFootprints=1, maskName='SAT')

◆ maskPixelsFromDefectList()

def lsst.ip.isr.isrFunctions.maskPixelsFromDefectList (   maskedImage,
  defectList,
  maskName = 'BAD' 
)
Set mask plane based on a defect list.

Parameters
----------
maskedImage : `lsst.afw.image.MaskedImage`
    Image to process.  Only the mask plane is updated.
defectList : `lsst.meas.algorithms.Defects`
    Defect list to mask.
maskName : str, optional
    Mask plane name to use.

Definition at line 141 of file isrFunctions.py.

141 def maskPixelsFromDefectList(maskedImage, defectList, maskName='BAD'):
142  """Set mask plane based on a defect list.
143 
144  Parameters
145  ----------
146  maskedImage : `lsst.afw.image.MaskedImage`
147  Image to process. Only the mask plane is updated.
148  defectList : `lsst.meas.algorithms.Defects`
149  Defect list to mask.
150  maskName : str, optional
151  Mask plane name to use.
152  """
153  return lsst.meas.algorithms.Defects(defectList).maskPixels(maskedImage, maskName=maskName)
154 
155 
156 @deprecated(reason="Replaced by Defects.fromMask() (will be removed after v18)",
157  category=FutureWarning)
def maskPixelsFromDefectList(maskedImage, defectList, maskName='BAD')

◆ overscanCorrection()

def lsst.ip.isr.isrFunctions.overscanCorrection (   ampMaskedImage,
  overscanImage,
  fitType = 'MEDIAN',
  order = 1,
  collapseRej = 3.0,
  statControl = None,
  overscanIsInt = True 
)
Apply overscan correction in place.

Parameters
----------
ampMaskedImage : `lsst.afw.image.MaskedImage`
    Image of amplifier to correct; modified.
overscanImage : `lsst.afw.image.Image` or `lsst.afw.image.MaskedImage`
    Image of overscan; modified.
fitType : `str`
    Type of fit for overscan correction. May be one of:

    - ``MEAN``: use mean of overscan.
    - ``MEANCLIP``: use clipped mean of overscan.
    - ``MEDIAN``: use median of overscan.
    - ``POLY``: fit with ordinary polynomial.
    - ``CHEB``: fit with Chebyshev polynomial.
    - ``LEG``: fit with Legendre polynomial.
    - ``NATURAL_SPLINE``: fit with natural spline.
    - ``CUBIC_SPLINE``: fit with cubic spline.
    - ``AKIMA_SPLINE``: fit with Akima spline.

order : `int`
    Polynomial order or number of spline knots; ignored unless
    ``fitType`` indicates a polynomial or spline.
statControl : `lsst.afw.math.StatisticsControl`
    Statistics control object.  In particular, we pay attention to numSigmaClip
overscanIsInt : `bool`
    Treat the overscan region as consisting of integers, even if it's been
    converted to float.  E.g. handle ties properly.

Returns
-------
result : `lsst.pipe.base.Struct`
    Result struct with components:

    - ``imageFit``: Value(s) removed from image (scalar or
        `lsst.afw.image.Image`)
    - ``overscanFit``: Value(s) removed from overscan (scalar or
        `lsst.afw.image.Image`)
    - ``overscanImage``: Overscan corrected overscan region
        (`lsst.afw.image.Image`)
Raises
------
pexExcept.Exception
    Raised if ``fitType`` is not an allowed value.

Notes
-----
The ``ampMaskedImage`` and ``overscanImage`` are modified, with the fit
subtracted. Note that the ``overscanImage`` should not be a subimage of
the ``ampMaskedImage``, to avoid being subtracted twice.

Debug plots are available for the SPLINE fitTypes by setting the
`debug.display` for `name` == "lsst.ip.isr.isrFunctions".  These
plots show the scatter plot of the overscan data (collapsed along
the perpendicular dimension) as a function of position on the CCD
(normalized between +/-1).

Definition at line 494 of file isrFunctions.py.

494  statControl=None, overscanIsInt=True):
495  """Apply overscan correction in place.
496 
497  Parameters
498  ----------
499  ampMaskedImage : `lsst.afw.image.MaskedImage`
500  Image of amplifier to correct; modified.
501  overscanImage : `lsst.afw.image.Image` or `lsst.afw.image.MaskedImage`
502  Image of overscan; modified.
503  fitType : `str`
504  Type of fit for overscan correction. May be one of:
505 
506  - ``MEAN``: use mean of overscan.
507  - ``MEANCLIP``: use clipped mean of overscan.
508  - ``MEDIAN``: use median of overscan.
509  - ``POLY``: fit with ordinary polynomial.
510  - ``CHEB``: fit with Chebyshev polynomial.
511  - ``LEG``: fit with Legendre polynomial.
512  - ``NATURAL_SPLINE``: fit with natural spline.
513  - ``CUBIC_SPLINE``: fit with cubic spline.
514  - ``AKIMA_SPLINE``: fit with Akima spline.
515 
516  order : `int`
517  Polynomial order or number of spline knots; ignored unless
518  ``fitType`` indicates a polynomial or spline.
519  statControl : `lsst.afw.math.StatisticsControl`
520  Statistics control object. In particular, we pay attention to numSigmaClip
521  overscanIsInt : `bool`
522  Treat the overscan region as consisting of integers, even if it's been
523  converted to float. E.g. handle ties properly.
524 
525  Returns
526  -------
527  result : `lsst.pipe.base.Struct`
528  Result struct with components:
529 
530  - ``imageFit``: Value(s) removed from image (scalar or
531  `lsst.afw.image.Image`)
532  - ``overscanFit``: Value(s) removed from overscan (scalar or
533  `lsst.afw.image.Image`)
534  - ``overscanImage``: Overscan corrected overscan region
535  (`lsst.afw.image.Image`)
536  Raises
537  ------
538  pexExcept.Exception
539  Raised if ``fitType`` is not an allowed value.
540 
541  Notes
542  -----
543  The ``ampMaskedImage`` and ``overscanImage`` are modified, with the fit
544  subtracted. Note that the ``overscanImage`` should not be a subimage of
545  the ``ampMaskedImage``, to avoid being subtracted twice.
546 
547  Debug plots are available for the SPLINE fitTypes by setting the
548  `debug.display` for `name` == "lsst.ip.isr.isrFunctions". These
549  plots show the scatter plot of the overscan data (collapsed along
550  the perpendicular dimension) as a function of position on the CCD
551  (normalized between +/-1).
552  """
553  ampImage = ampMaskedImage.getImage()
554  if statControl is None:
555  statControl = afwMath.StatisticsControl()
556 
557  numSigmaClip = statControl.getNumSigmaClip()
558 
559  if fitType in ('MEAN', 'MEANCLIP'):
560  fitType = afwMath.stringToStatisticsProperty(fitType)
561  offImage = afwMath.makeStatistics(overscanImage, fitType, statControl).getValue()
562  overscanFit = offImage
563  elif fitType in ('MEDIAN',):
564  if overscanIsInt:
565  # we need an image with integer pixels to handle ties properly
566  if hasattr(overscanImage, "image"):
567  imageI = overscanImage.image.convertI()
568  overscanImageI = afwImage.MaskedImageI(imageI, overscanImage.mask, overscanImage.variance)
569  else:
570  overscanImageI = overscanImage.convertI()
571  else:
572  overscanImageI = overscanImage
573 
574  fitType = afwMath.stringToStatisticsProperty(fitType)
575  offImage = afwMath.makeStatistics(overscanImageI, fitType, statControl).getValue()
576  overscanFit = offImage
577 
578  if overscanIsInt:
579  del overscanImageI
580  elif fitType in ('POLY', 'CHEB', 'LEG', 'NATURAL_SPLINE', 'CUBIC_SPLINE', 'AKIMA_SPLINE'):
581  if hasattr(overscanImage, "getImage"):
582  biasArray = overscanImage.getImage().getArray()
583  biasArray = numpy.ma.masked_where(overscanImage.getMask().getArray() & statControl.getAndMask(),
584  biasArray)
585  else:
586  biasArray = overscanImage.getArray()
587  # Fit along the long axis, so collapse along each short row and fit the resulting array
588  shortInd = numpy.argmin(biasArray.shape)
589  if shortInd == 0:
590  # Convert to some 'standard' representation to make things easier
591  biasArray = numpy.transpose(biasArray)
592 
593  # Do a single round of clipping to weed out CR hits and signal leaking into the overscan
594  percentiles = numpy.percentile(biasArray, [25.0, 50.0, 75.0], axis=1)
595  medianBiasArr = percentiles[1]
596  stdevBiasArr = 0.74*(percentiles[2] - percentiles[0]) # robust stdev
597  diff = numpy.abs(biasArray - medianBiasArr[:, numpy.newaxis])
598  biasMaskedArr = numpy.ma.masked_where(diff > numSigmaClip*stdevBiasArr[:, numpy.newaxis], biasArray)
599  collapsed = numpy.mean(biasMaskedArr, axis=1)
600  if collapsed.mask.sum() > 0:
601  collapsed.data[collapsed.mask] = numpy.mean(biasArray.data[collapsed.mask], axis=1)
602  del biasArray, percentiles, stdevBiasArr, diff, biasMaskedArr
603 
604  if shortInd == 0:
605  collapsed = numpy.transpose(collapsed)
606 
607  num = len(collapsed)
608  indices = 2.0*numpy.arange(num)/float(num) - 1.0
609 
610  if fitType in ('POLY', 'CHEB', 'LEG'):
611  # A numpy polynomial
612  poly = numpy.polynomial
613  fitter, evaler = {"POLY": (poly.polynomial.polyfit, poly.polynomial.polyval),
614  "CHEB": (poly.chebyshev.chebfit, poly.chebyshev.chebval),
615  "LEG": (poly.legendre.legfit, poly.legendre.legval),
616  }[fitType]
617 
618  coeffs = fitter(indices, collapsed, order)
619  fitBiasArr = evaler(indices, coeffs)
620  elif 'SPLINE' in fitType:
621  # An afw interpolation
622  numBins = order
623  #
624  # numpy.histogram needs a real array for the mask, but numpy.ma "optimises" the case
625  # no-values-are-masked by replacing the mask array by a scalar, numpy.ma.nomask
626  #
627  # Issue DM-415
628  #
629  collapsedMask = collapsed.mask
630  try:
631  if collapsedMask == numpy.ma.nomask:
632  collapsedMask = numpy.array(len(collapsed)*[numpy.ma.nomask])
633  except ValueError: # If collapsedMask is an array the test fails [needs .all()]
634  pass
635 
636  numPerBin, binEdges = numpy.histogram(indices, bins=numBins,
637  weights=1-collapsedMask.astype(int))
638  # Binning is just a histogram, with weights equal to the values.
639  # Use a similar trick to get the bin centers (this deals with different numbers per bin).
640  with numpy.errstate(invalid="ignore"): # suppress NAN warnings
641  values = numpy.histogram(indices, bins=numBins,
642  weights=collapsed.data*~collapsedMask)[0]/numPerBin
643  binCenters = numpy.histogram(indices, bins=numBins,
644  weights=indices*~collapsedMask)[0]/numPerBin
645  interp = afwMath.makeInterpolate(binCenters.astype(float)[numPerBin > 0],
646  values.astype(float)[numPerBin > 0],
648  fitBiasArr = numpy.array([interp.interpolate(i) for i in indices])
649 
650  import lsstDebug
651  if lsstDebug.Info(__name__).display:
652  import matplotlib.pyplot as plot
653  figure = plot.figure(1)
654  figure.clear()
655  axes = figure.add_axes((0.1, 0.1, 0.8, 0.8))
656  axes.plot(indices[~collapsedMask], collapsed[~collapsedMask], 'k+')
657  if collapsedMask.sum() > 0:
658  axes.plot(indices[collapsedMask], collapsed.data[collapsedMask], 'b+')
659  axes.plot(indices, fitBiasArr, 'r-')
660  plot.xlabel("centered/scaled position along overscan region")
661  plot.ylabel("pixel value/fit value")
662  figure.show()
663  prompt = "Press Enter or c to continue [chp]... "
664  while True:
665  ans = input(prompt).lower()
666  if ans in ("", "c",):
667  break
668  if ans in ("p",):
669  import pdb
670  pdb.set_trace()
671  elif ans in ("h", ):
672  print("h[elp] c[ontinue] p[db]")
673  plot.close()
674 
675  offImage = ampImage.Factory(ampImage.getDimensions())
676  offArray = offImage.getArray()
677  overscanFit = afwImage.ImageF(overscanImage.getDimensions())
678  overscanArray = overscanFit.getArray()
679  if shortInd == 1:
680  offArray[:, :] = fitBiasArr[:, numpy.newaxis]
681  overscanArray[:, :] = fitBiasArr[:, numpy.newaxis]
682  else:
683  offArray[:, :] = fitBiasArr[numpy.newaxis, :]
684  overscanArray[:, :] = fitBiasArr[numpy.newaxis, :]
685 
686  # We don't trust any extrapolation: mask those pixels as SUSPECT
687  # This will occur when the top and or bottom edges of the overscan
688  # contain saturated values. The values will be extrapolated from
689  # the surrounding pixels, but we cannot entirely trust the value of
690  # the extrapolation, and will mark the image mask plane to flag the
691  # image as such.
692  mask = ampMaskedImage.getMask()
693  maskArray = mask.getArray() if shortInd == 1 else mask.getArray().transpose()
694  suspect = mask.getPlaneBitMask("SUSPECT")
695  try:
696  if collapsed.mask == numpy.ma.nomask:
697  # There is no mask, so the whole array is fine
698  pass
699  except ValueError: # If collapsed.mask is an array the test fails [needs .all()]
700  for low in range(num):
701  if not collapsed.mask[low]:
702  break
703  if low > 0:
704  maskArray[:low, :] |= suspect
705  for high in range(1, num):
706  if not collapsed.mask[-high]:
707  break
708  if high > 1:
709  maskArray[-high:, :] |= suspect
710 
711  else:
712  raise pexExcept.Exception('%s : %s an invalid overscan type' % ("overscanCorrection", fitType))
713  ampImage -= offImage
714  overscanImage -= overscanFit
715  return Struct(imageFit=offImage, overscanFit=overscanFit, overscanImage=overscanImage)
716 
717 
Interpolate::Style stringToInterpStyle(std::string const &style)
Conversion function to switch a string to an Interpolate::Style.
Definition: Interpolate.cc:257
Provides consistent interface for LSST exceptions.
Definition: Exception.h:107
Statistics makeStatistics(lsst::afw::math::MaskedVector< EntryT > const &mv, std::vector< WeightPixel > const &vweights, int const flags, StatisticsControl const &sctrl=StatisticsControl())
The makeStatistics() overload to handle lsst::afw::math::MaskedVector<>
Definition: Statistics.h:520
Property stringToStatisticsProperty(std::string const property)
Conversion function to switch a string to a Property (see Statistics.h)
Definition: Statistics.cc:747
Pass parameters to a Statistics object.
Definition: Statistics.h:93
std::shared_ptr< Interpolate > makeInterpolate(ndarray::Array< double const, 1 > const &x, ndarray::Array< double const, 1 > const &y, Interpolate::Style const style=Interpolate::AKIMA_SPLINE)
Definition: Interpolate.cc:353

◆ saturationCorrection()

def lsst.ip.isr.isrFunctions.saturationCorrection (   maskedImage,
  saturation,
  fwhm,
  growFootprints = 1,
  interpolate = True,
  maskName = 'SAT',
  fallbackValue = None 
)
Mark saturated pixels and optionally interpolate over them

Parameters
----------
maskedImage : `lsst.afw.image.MaskedImage`
    Image to process.
saturation  : scalar
    Saturation level used as the detection threshold.
fwhm : scalar
    FWHM of double Gaussian smoothing kernel.
growFootprints : scalar, optional
    Number of pixels to grow footprints of detected regions.
interpolate : Bool, optional
    If True, saturated pixels are interpolated over.
maskName : str, optional
    Mask plane name.
fallbackValue : scalar, optional
    Value of last resort for interpolation.

Definition at line 248 of file isrFunctions.py.

248  fallbackValue=None):
249  """Mark saturated pixels and optionally interpolate over them
250 
251  Parameters
252  ----------
253  maskedImage : `lsst.afw.image.MaskedImage`
254  Image to process.
255  saturation : scalar
256  Saturation level used as the detection threshold.
257  fwhm : scalar
258  FWHM of double Gaussian smoothing kernel.
259  growFootprints : scalar, optional
260  Number of pixels to grow footprints of detected regions.
261  interpolate : Bool, optional
262  If True, saturated pixels are interpolated over.
263  maskName : str, optional
264  Mask plane name.
265  fallbackValue : scalar, optional
266  Value of last resort for interpolation.
267  """
268  defectList = makeThresholdMask(
269  maskedImage=maskedImage,
270  threshold=saturation,
271  growFootprints=growFootprints,
272  maskName=maskName,
273  )
274  if interpolate:
275  interpolateDefectList(maskedImage, defectList, fwhm, fallbackValue=fallbackValue)
276 
277  return maskedImage
278 
279 
def interpolateDefectList(maskedImage, defectList, fwhm, fallbackValue=None)
Definition: isrFunctions.py:78
def makeThresholdMask(maskedImage, threshold, growFootprints=1, maskName='SAT')

◆ setBadRegions()

def lsst.ip.isr.isrFunctions.setBadRegions (   exposure,
  badStatistic = "MEDIAN" 
)
Set all BAD areas of the chip to the average of the rest of the exposure

Parameters
----------
exposure : `lsst.afw.image.Exposure`
    Exposure to mask.  The exposure mask is modified.
badStatistic : `str`, optional
    Statistic to use to generate the replacement value from the
    image data.  Allowed values are 'MEDIAN' or 'MEANCLIP'.

Returns
-------
badPixelCount : scalar
    Number of bad pixels masked.
badPixelValue : scalar
    Value substituted for bad pixels.

Raises
------
RuntimeError
    Raised if `badStatistic` is not an allowed value.

Definition at line 1043 of file isrFunctions.py.

1043 def setBadRegions(exposure, badStatistic="MEDIAN"):
1044  """Set all BAD areas of the chip to the average of the rest of the exposure
1045 
1046  Parameters
1047  ----------
1048  exposure : `lsst.afw.image.Exposure`
1049  Exposure to mask. The exposure mask is modified.
1050  badStatistic : `str`, optional
1051  Statistic to use to generate the replacement value from the
1052  image data. Allowed values are 'MEDIAN' or 'MEANCLIP'.
1053 
1054  Returns
1055  -------
1056  badPixelCount : scalar
1057  Number of bad pixels masked.
1058  badPixelValue : scalar
1059  Value substituted for bad pixels.
1060 
1061  Raises
1062  ------
1063  RuntimeError
1064  Raised if `badStatistic` is not an allowed value.
1065  """
1066  if badStatistic == "MEDIAN":
1067  statistic = afwMath.MEDIAN
1068  elif badStatistic == "MEANCLIP":
1069  statistic = afwMath.MEANCLIP
1070  else:
1071  raise RuntimeError("Impossible method %s of bad region correction" % badStatistic)
1072 
1073  mi = exposure.getMaskedImage()
1074  mask = mi.getMask()
1075  BAD = mask.getPlaneBitMask("BAD")
1076  INTRP = mask.getPlaneBitMask("INTRP")
1077 
1078  sctrl = afwMath.StatisticsControl()
1079  sctrl.setAndMask(BAD)
1080  value = afwMath.makeStatistics(mi, statistic, sctrl).getValue()
1081 
1082  maskArray = mask.getArray()
1083  imageArray = mi.getImage().getArray()
1084  badPixels = numpy.logical_and((maskArray & BAD) > 0, (maskArray & INTRP) == 0)
1085  imageArray[:] = numpy.where(badPixels, value, imageArray)
1086 
1087  return badPixels.sum(), value
1088 
def setBadRegions(exposure, badStatistic="MEDIAN")
Statistics makeStatistics(lsst::afw::math::MaskedVector< EntryT > const &mv, std::vector< WeightPixel > const &vweights, int const flags, StatisticsControl const &sctrl=StatisticsControl())
The makeStatistics() overload to handle lsst::afw::math::MaskedVector<>
Definition: Statistics.h:520
Pass parameters to a Statistics object.
Definition: Statistics.h:93

◆ transposeDefectList()

def lsst.ip.isr.isrFunctions.transposeDefectList (   defectList)
Make a transposed copy of a defect list.

Parameters
----------
defectList : `lsst.meas.algorithms.Defects`
    Input list of defects.

Returns
-------
retDefectList : `lsst.meas.algorithms.Defects`
    Transposed list of defects.

Definition at line 121 of file isrFunctions.py.

121 def transposeDefectList(defectList):
122  """Make a transposed copy of a defect list.
123 
124  Parameters
125  ----------
126  defectList : `lsst.meas.algorithms.Defects`
127  Input list of defects.
128 
129  Returns
130  -------
131  retDefectList : `lsst.meas.algorithms.Defects`
132  Transposed list of defects.
133  """
134  if isinstance(defectList, measAlg.Defects):
135  return defectList.transpose()
136  return measAlg.Defects(defectList).transpose()
137 
138 
139 @deprecated(reason="Replaced by Defects.maskPixels() (will be removed after v18)",
140  category=FutureWarning)
def transposeDefectList(defectList)

◆ transposeMaskedImage()

def lsst.ip.isr.isrFunctions.transposeMaskedImage (   maskedImage)
Make a transposed copy of a masked image.

Parameters
----------
maskedImage : `lsst.afw.image.MaskedImage`
    Image to process.

Returns
-------
transposed : `lsst.afw.image.MaskedImage`
    The transposed copy of the input image.

Definition at line 58 of file isrFunctions.py.

58 def transposeMaskedImage(maskedImage):
59  """Make a transposed copy of a masked image.
60 
61  Parameters
62  ----------
63  maskedImage : `lsst.afw.image.MaskedImage`
64  Image to process.
65 
66  Returns
67  -------
68  transposed : `lsst.afw.image.MaskedImage`
69  The transposed copy of the input image.
70  """
71  transposed = maskedImage.Factory(lsst.geom.Extent2I(maskedImage.getHeight(), maskedImage.getWidth()))
72  transposed.getImage().getArray()[:] = maskedImage.getImage().getArray().T
73  transposed.getMask().getArray()[:] = maskedImage.getMask().getArray().T
74  transposed.getVariance().getArray()[:] = maskedImage.getVariance().getArray().T
75  return transposed
76 
77 
def transposeMaskedImage(maskedImage)
Definition: isrFunctions.py:58

◆ trimToMatchCalibBBox()

def lsst.ip.isr.isrFunctions.trimToMatchCalibBBox (   rawMaskedImage,
  calibMaskedImage 
)
Compute number of edge trim pixels to match the calibration data.

Use the dimension difference between the raw exposure and the
calibration exposure to compute the edge trim pixels.  This trim
is applied symmetrically, with the same number of pixels masked on
each side.

Parameters
----------
rawMaskedImage : `lsst.afw.image.MaskedImage`
    Image to trim.
calibMaskedImage : `lsst.afw.image.MaskedImage`
    Calibration image to draw new bounding box from.

Returns
-------
replacementMaskedImage : `lsst.afw.image.MaskedImage`
    ``rawMaskedImage`` trimmed to the appropriate size
Raises
------
RuntimeError
   Rasied if ``rawMaskedImage`` cannot be symmetrically trimmed to
   match ``calibMaskedImage``.

Definition at line 280 of file isrFunctions.py.

280 def trimToMatchCalibBBox(rawMaskedImage, calibMaskedImage):
281  """Compute number of edge trim pixels to match the calibration data.
282 
283  Use the dimension difference between the raw exposure and the
284  calibration exposure to compute the edge trim pixels. This trim
285  is applied symmetrically, with the same number of pixels masked on
286  each side.
287 
288  Parameters
289  ----------
290  rawMaskedImage : `lsst.afw.image.MaskedImage`
291  Image to trim.
292  calibMaskedImage : `lsst.afw.image.MaskedImage`
293  Calibration image to draw new bounding box from.
294 
295  Returns
296  -------
297  replacementMaskedImage : `lsst.afw.image.MaskedImage`
298  ``rawMaskedImage`` trimmed to the appropriate size
299  Raises
300  ------
301  RuntimeError
302  Rasied if ``rawMaskedImage`` cannot be symmetrically trimmed to
303  match ``calibMaskedImage``.
304  """
305  nx, ny = rawMaskedImage.getBBox().getDimensions() - calibMaskedImage.getBBox().getDimensions()
306  if nx != ny:
307  raise RuntimeError("Raw and calib maskedImages are trimmed differently in X and Y.")
308  if nx % 2 != 0:
309  raise RuntimeError("Calibration maskedImage is trimmed unevenly in X.")
310  if nx < 0:
311  raise RuntimeError("Calibration maskedImage is larger than raw data.")
312 
313  nEdge = nx//2
314  if nEdge > 0:
315  replacementMaskedImage = rawMaskedImage[nEdge:-nEdge, nEdge:-nEdge, afwImage.LOCAL]
316  SourceDetectionTask.setEdgeBits(
317  rawMaskedImage,
318  replacementMaskedImage.getBBox(),
319  rawMaskedImage.getMask().getPlaneBitMask("EDGE")
320  )
321  else:
322  replacementMaskedImage = rawMaskedImage
323 
324  return replacementMaskedImage
325 
326 
def trimToMatchCalibBBox(rawMaskedImage, calibMaskedImage)

◆ updateVariance()

def lsst.ip.isr.isrFunctions.updateVariance (   maskedImage,
  gain,
  readNoise 
)
Set the variance plane based on the image plane.

Parameters
----------
maskedImage : `lsst.afw.image.MaskedImage`
    Image to process.  The variance plane is modified.
gain : scalar
    The amplifier gain in electrons/ADU.
readNoise : scalar
    The amplifier read nmoise in ADU/pixel.

Definition at line 400 of file isrFunctions.py.

400 def updateVariance(maskedImage, gain, readNoise):
401  """Set the variance plane based on the image plane.
402 
403  Parameters
404  ----------
405  maskedImage : `lsst.afw.image.MaskedImage`
406  Image to process. The variance plane is modified.
407  gain : scalar
408  The amplifier gain in electrons/ADU.
409  readNoise : scalar
410  The amplifier read nmoise in ADU/pixel.
411  """
412  var = maskedImage.getVariance()
413  var[:] = maskedImage.getImage()
414  var /= gain
415  var += readNoise**2
416 
417 
def updateVariance(maskedImage, gain, readNoise)

◆ widenSaturationTrails()

def lsst.ip.isr.isrFunctions.widenSaturationTrails (   mask)
Grow the saturation trails by an amount dependent on the width of the trail.

Parameters
----------
mask : `lsst.afw.image.Mask`
    Mask which will have the saturated areas grown.

Definition at line 996 of file isrFunctions.py.

996 def widenSaturationTrails(mask):
997  """Grow the saturation trails by an amount dependent on the width of the trail.
998 
999  Parameters
1000  ----------
1001  mask : `lsst.afw.image.Mask`
1002  Mask which will have the saturated areas grown.
1003  """
1004 
1005  extraGrowDict = {}
1006  for i in range(1, 6):
1007  extraGrowDict[i] = 0
1008  for i in range(6, 8):
1009  extraGrowDict[i] = 1
1010  for i in range(8, 10):
1011  extraGrowDict[i] = 3
1012  extraGrowMax = 4
1013 
1014  if extraGrowMax <= 0:
1015  return
1016 
1017  saturatedBit = mask.getPlaneBitMask("SAT")
1018 
1019  xmin, ymin = mask.getBBox().getMin()
1020  width = mask.getWidth()
1021 
1022  thresh = afwDetection.Threshold(saturatedBit, afwDetection.Threshold.BITMASK)
1023  fpList = afwDetection.FootprintSet(mask, thresh).getFootprints()
1024 
1025  for fp in fpList:
1026  for s in fp.getSpans():
1027  x0, x1 = s.getX0(), s.getX1()
1028 
1029  extraGrow = extraGrowDict.get(x1 - x0 + 1, extraGrowMax)
1030  if extraGrow > 0:
1031  y = s.getY() - ymin
1032  x0 -= xmin + extraGrow
1033  x1 -= xmin - extraGrow
1034 
1035  if x0 < 0:
1036  x0 = 0
1037  if x1 >= width - 1:
1038  x1 = width - 1
1039 
1040  mask.array[y, x0:x1+1] |= saturatedBit
1041 
1042