LSST Applications g0fba68d861+bb7a7cfa1f,g1ec0fe41b4+f536777771,g1fd858c14a+470a99fdf4,g216c3ac8a7+0d4d80193f,g35bb328faa+fcb1d3bbc8,g4d2262a081+23bd310d1b,g53246c7159+fcb1d3bbc8,g56a49b3a55+369644a549,g5a012ec0e7+3632fc3ff3,g60b5630c4e+3bfb9058a5,g67b6fd64d1+ed4b5058f4,g78460c75b0+2f9a1b4bcd,g786e29fd12+cf7ec2a62a,g8180f54f50+60bd39f3b6,g8352419a5c+fcb1d3bbc8,g87d29937c9+57a68d035f,g8852436030+4699110379,g89139ef638+ed4b5058f4,g9125e01d80+fcb1d3bbc8,g94187f82dc+3bfb9058a5,g989de1cb63+ed4b5058f4,g9ccd5d7f00+b7cae620c0,g9d31334357+3bfb9058a5,g9f33ca652e+00883ace41,gabe3b4be73+1e0a283bba,gabf8522325+fa80ff7197,gb1101e3267+27b24065a3,gb58c049af0+f03b321e39,gb89ab40317+ed4b5058f4,gc0af124501+708fe67c54,gcf25f946ba+4699110379,gd6cbbdb0b4+bb83cc51f8,gde0f65d7ad+acd5afb0eb,ge1ad929117+3bfb9058a5,ge278dab8ac+d65b3c2b70,ge410e46f29+ed4b5058f4,gf5e32f922b+fcb1d3bbc8,gf67bdafdda+ed4b5058f4,w.2025.17
LSST Data Management Base Package
Loading...
Searching...
No Matches
lsst.ip.diffim.subtractImages.AlardLuptonPreconvolveSubtractTask Class Reference
Inheritance diagram for lsst.ip.diffim.subtractImages.AlardLuptonPreconvolveSubtractTask:
lsst.ip.diffim.subtractImages.AlardLuptonSubtractTask

Public Member Functions

 run (self, template, science, sources, visitSummary=None)
 
 runPreconvolve (self, template, science, matchedScience, selectSources, preConvKernel)
 
 computeImageMetrics (self, science, difference, stars)
 
 runConvolveTemplate (self, template, science, selectSources)
 
 runConvolveScience (self, template, science, selectSources)
 
 finalize (self, template, science, difference, kernel, templateMatched=True, preConvMode=False, preConvKernel=None, spatiallyVarying=False)
 
 updateMasks (self, template, science)
 

Public Attributes

 convolutionControl = lsst.afw.math.ConvolutionControl()
 
 templatePsfSize = getPsfFwhm(template.psf)
 
 sciencePsfSize = getPsfFwhm(science.psf)
 
 log = self.runConvolveScience(template, science, selectSources)
 

Static Public Attributes

 ConfigClass = AlardLuptonSubtractConfig
 

Protected Member Functions

 _applyExternalCalibrations (self, exposure, visitSummary)
 
 _calculateMagLim (self, exposure, nsigma=5.0, fallbackPsfSize=None)
 
 _convolveExposure (self, exposure, kernel, convolutionControl, bbox=None, psf=None, photoCalib=None, interpolateBadMaskPlanes=False)
 
 _sourceSelector (self, sources, mask)
 
 _prepareInputs (self, template, science, visitSummary=None)
 
 _clearMask (self, mask, clearMaskPlanes=None)
 

Static Protected Member Functions

 _validateExposures (template, science)
 
 _checkMask (mask, sources, excludeMaskPlanes, checkAdjacent=True)
 
 _renameMaskPlanes (mask, maskPlane, newMaskPlane)
 

Static Protected Attributes

str _DefaultName = "alardLuptonSubtract"
 

Detailed Description

Subtract a template from a science image, convolving the science image
before computing the kernel, and also convolving the template before
subtraction.

Definition at line 1071 of file subtractImages.py.

Member Function Documentation

◆ _applyExternalCalibrations()

lsst.ip.diffim.subtractImages.AlardLuptonSubtractTask._applyExternalCalibrations ( self,
exposure,
visitSummary )
protectedinherited
Replace calibrations (psf, and ApCorrMap) on this exposure with
external ones.".

Parameters
----------
exposure : `lsst.afw.image.exposure.Exposure`
    Input exposure to adjust calibrations.
visitSummary : `lsst.afw.table.ExposureCatalog`
    Exposure catalog with external calibrations to be applied. Catalog
    uses the detector id for the catalog id, sorted on id for fast
    lookup.

Returns
-------
exposure : `lsst.afw.image.exposure.Exposure`
    Exposure with adjusted calibrations.

Definition at line 290 of file subtractImages.py.

290 def _applyExternalCalibrations(self, exposure, visitSummary):
291 """Replace calibrations (psf, and ApCorrMap) on this exposure with
292 external ones.".
293
294 Parameters
295 ----------
296 exposure : `lsst.afw.image.exposure.Exposure`
297 Input exposure to adjust calibrations.
298 visitSummary : `lsst.afw.table.ExposureCatalog`
299 Exposure catalog with external calibrations to be applied. Catalog
300 uses the detector id for the catalog id, sorted on id for fast
301 lookup.
302
303 Returns
304 -------
305 exposure : `lsst.afw.image.exposure.Exposure`
306 Exposure with adjusted calibrations.
307 """
308 detectorId = exposure.info.getDetector().getId()
309
310 row = visitSummary.find(detectorId)
311 if row is None:
312 self.log.warning("Detector id %s not found in external calibrations catalog; "
313 "Using original calibrations.", detectorId)
314 else:
315 psf = row.getPsf()
316 apCorrMap = row.getApCorrMap()
317 if psf is None:
318 self.log.warning("Detector id %s has None for psf in "
319 "external calibrations catalog; Using original psf and aperture correction.",
320 detectorId)
321 elif apCorrMap is None:
322 self.log.warning("Detector id %s has None for apCorrMap in "
323 "external calibrations catalog; Using original psf and aperture correction.",
324 detectorId)
325 else:
326 exposure.setPsf(psf)
327 exposure.info.setApCorrMap(apCorrMap)
328
329 return exposure
330

◆ _calculateMagLim()

lsst.ip.diffim.subtractImages.AlardLuptonSubtractTask._calculateMagLim ( self,
exposure,
nsigma = 5.0,
fallbackPsfSize = None )
protectedinherited
Calculate an exposure's limiting magnitude.

This method uses the photometric zeropoint together with the
PSF size from the average position of the exposure.

Parameters
----------
exposure : `lsst.afw.image.Exposure`
    The target exposure to calculate the limiting magnitude for.
nsigma : `float`, optional
    The detection threshold in sigma.
fallbackPsfSize : `float`, optional
    PSF FWHM to use in the event the exposure PSF cannot be retrieved.

Returns
-------
maglim : `astropy.units.Quantity`
    The limiting magnitude of the exposure, or np.nan.

Definition at line 729 of file subtractImages.py.

729 def _calculateMagLim(self, exposure, nsigma=5.0, fallbackPsfSize=None):
730 """Calculate an exposure's limiting magnitude.
731
732 This method uses the photometric zeropoint together with the
733 PSF size from the average position of the exposure.
734
735 Parameters
736 ----------
737 exposure : `lsst.afw.image.Exposure`
738 The target exposure to calculate the limiting magnitude for.
739 nsigma : `float`, optional
740 The detection threshold in sigma.
741 fallbackPsfSize : `float`, optional
742 PSF FWHM to use in the event the exposure PSF cannot be retrieved.
743
744 Returns
745 -------
746 maglim : `astropy.units.Quantity`
747 The limiting magnitude of the exposure, or np.nan.
748 """
749 try:
750 psf = exposure.getPsf()
751 psf_shape = psf.computeShape(psf.getAveragePosition())
752 except (lsst.pex.exceptions.InvalidParameterError, afwDetection.InvalidPsfError):
753 if fallbackPsfSize is not None:
754 self.log.info("Unable to evaluate PSF, using fallback FWHM %f", fallbackPsfSize)
755 psf_area = np.pi*(fallbackPsfSize/2)**2
756 zeropoint = exposure.getPhotoCalib().instFluxToMagnitude(1)
757 maglim = zeropoint - 2.5*np.log10(nsigma*np.sqrt(psf_area))
758 else:
759 self.log.info("Unable to evaluate PSF, setting maglim to nan")
760 maglim = np.nan
761 else:
762 # Get a more accurate area than `psf_shape.getArea()` via moments
763 psf_area = np.pi*np.sqrt(psf_shape.getIxx()*psf_shape.getIyy())
764 zeropoint = exposure.getPhotoCalib().instFluxToMagnitude(1)
765 maglim = zeropoint - 2.5*np.log10(nsigma*np.sqrt(psf_area))
766 finally:
767 return maglim
768
Reports invalid arguments.
Definition Runtime.h:66

◆ _checkMask()

lsst.ip.diffim.subtractImages.AlardLuptonSubtractTask._checkMask ( mask,
sources,
excludeMaskPlanes,
checkAdjacent = True )
staticprotectedinherited
Exclude sources that are located on or adjacent to masked pixels.

Parameters
----------
mask : `lsst.afw.image.Mask`
    The image mask plane to use to reject sources
    based on the location of their centroid on the ccd.
sources : `lsst.afw.table.SourceCatalog`
    The source catalog to evaluate.
excludeMaskPlanes : `list` of `str`
    List of the names of the mask planes to exclude.

Returns
-------
flags : `numpy.ndarray` of `bool`
    Array indicating whether each source in the catalog should be
    kept (True) or rejected (False) based on the value of the
    mask plane at its location.

Definition at line 898 of file subtractImages.py.

898 def _checkMask(mask, sources, excludeMaskPlanes, checkAdjacent=True):
899 """Exclude sources that are located on or adjacent to masked pixels.
900
901 Parameters
902 ----------
903 mask : `lsst.afw.image.Mask`
904 The image mask plane to use to reject sources
905 based on the location of their centroid on the ccd.
906 sources : `lsst.afw.table.SourceCatalog`
907 The source catalog to evaluate.
908 excludeMaskPlanes : `list` of `str`
909 List of the names of the mask planes to exclude.
910
911 Returns
912 -------
913 flags : `numpy.ndarray` of `bool`
914 Array indicating whether each source in the catalog should be
915 kept (True) or rejected (False) based on the value of the
916 mask plane at its location.
917 """
918 setExcludeMaskPlanes = [
919 maskPlane for maskPlane in excludeMaskPlanes if maskPlane in mask.getMaskPlaneDict()
920 ]
921
922 excludePixelMask = mask.getPlaneBitMask(setExcludeMaskPlanes)
923
924 xv = (np.rint(sources.getX() - mask.getX0())).astype(int)
925 yv = (np.rint(sources.getY() - mask.getY0())).astype(int)
926
927 flags = np.ones(len(sources), dtype=bool)
928 if checkAdjacent:
929 pixRange = (0, -1, 1)
930 else:
931 pixRange = (0,)
932 for j in pixRange:
933 for i in pixRange:
934 mv = mask.array[yv + j, xv + i]
935 flags *= np.bitwise_and(mv, excludePixelMask) == 0
936 return flags
937

◆ _clearMask()

lsst.ip.diffim.subtractImages.AlardLuptonSubtractTask._clearMask ( self,
mask,
clearMaskPlanes = None )
protectedinherited
Clear the mask plane of an exposure.

Parameters
----------
mask : `lsst.afw.image.Mask`
    The mask plane to erase, which will be modified in place.
clearMaskPlanes : `list` of `str`, optional
    Erase the specified mask planes.
    If not supplied, the entire mask will be erased.

Definition at line 1043 of file subtractImages.py.

1043 def _clearMask(self, mask, clearMaskPlanes=None):
1044 """Clear the mask plane of an exposure.
1045
1046 Parameters
1047 ----------
1048 mask : `lsst.afw.image.Mask`
1049 The mask plane to erase, which will be modified in place.
1050 clearMaskPlanes : `list` of `str`, optional
1051 Erase the specified mask planes.
1052 If not supplied, the entire mask will be erased.
1053 """
1054 if clearMaskPlanes is None:
1055 clearMaskPlanes = list(mask.getMaskPlaneDict().keys())
1056
1057 bitMaskToClear = mask.getPlaneBitMask(clearMaskPlanes)
1058 mask &= ~bitMaskToClear
1059
1060

◆ _convolveExposure()

lsst.ip.diffim.subtractImages.AlardLuptonSubtractTask._convolveExposure ( self,
exposure,
kernel,
convolutionControl,
bbox = None,
psf = None,
photoCalib = None,
interpolateBadMaskPlanes = False )
protectedinherited
Convolve an exposure with the given kernel.

Parameters
----------
exposure : `lsst.afw.Exposure`
    exposure to convolve.
kernel : `lsst.afw.math.LinearCombinationKernel`
    PSF matching kernel computed in the ``makeKernel`` subtask.
convolutionControl : `lsst.afw.math.ConvolutionControl`
    Configuration for convolve algorithm.
bbox : `lsst.geom.Box2I`, optional
    Bounding box to trim the convolved exposure to.
psf : `lsst.afw.detection.Psf`, optional
    Point spread function (PSF) to set for the convolved exposure.
photoCalib : `lsst.afw.image.PhotoCalib`, optional
    Photometric calibration of the convolved exposure.

Returns
-------
convolvedExp : `lsst.afw.Exposure`
    The convolved image.

Definition at line 799 of file subtractImages.py.

804 ):
805 """Convolve an exposure with the given kernel.
806
807 Parameters
808 ----------
809 exposure : `lsst.afw.Exposure`
810 exposure to convolve.
811 kernel : `lsst.afw.math.LinearCombinationKernel`
812 PSF matching kernel computed in the ``makeKernel`` subtask.
813 convolutionControl : `lsst.afw.math.ConvolutionControl`
814 Configuration for convolve algorithm.
815 bbox : `lsst.geom.Box2I`, optional
816 Bounding box to trim the convolved exposure to.
817 psf : `lsst.afw.detection.Psf`, optional
818 Point spread function (PSF) to set for the convolved exposure.
819 photoCalib : `lsst.afw.image.PhotoCalib`, optional
820 Photometric calibration of the convolved exposure.
821
822 Returns
823 -------
824 convolvedExp : `lsst.afw.Exposure`
825 The convolved image.
826 """
827 convolvedExposure = exposure.clone()
828 if psf is not None:
829 convolvedExposure.setPsf(psf)
830 if photoCalib is not None:
831 convolvedExposure.setPhotoCalib(photoCalib)
832 if interpolateBadMaskPlanes and self.config.badMaskPlanes is not None:
833 nInterp = _interpolateImage(convolvedExposure.maskedImage,
834 self.config.badMaskPlanes)
835 self.metadata["nInterpolated"] = nInterp
836 convolvedImage = lsst.afw.image.MaskedImageF(convolvedExposure.getBBox())
837 lsst.afw.math.convolve(convolvedImage, convolvedExposure.maskedImage, kernel, convolutionControl)
838 convolvedExposure.setMaskedImage(convolvedImage)
839 if bbox is None:
840 return convolvedExposure
841 else:
842 return convolvedExposure[bbox]
843
void convolve(OutImageT &convolvedImage, InImageT const &inImage, KernelT const &kernel, ConvolutionControl const &convolutionControl=ConvolutionControl())
Convolve an Image or MaskedImage with a Kernel, setting pixels of an existing output image.

◆ _prepareInputs()

lsst.ip.diffim.subtractImages.AlardLuptonSubtractTask._prepareInputs ( self,
template,
science,
visitSummary = None )
protectedinherited
Perform preparatory calculations common to all Alard&Lupton Tasks.

Parameters
----------
template : `lsst.afw.image.ExposureF`
    Template exposure, warped to match the science exposure. The
    variance plane of the template image is modified in place.
science : `lsst.afw.image.ExposureF`
    Science exposure to subtract from the template. The variance plane
    of the science image is modified in place.
visitSummary : `lsst.afw.table.ExposureCatalog`, optional
    Exposure catalog with external calibrations to be applied.  Catalog
    uses the detector id for the catalog id, sorted on id for fast
    lookup.

Definition at line 938 of file subtractImages.py.

938 def _prepareInputs(self, template, science, visitSummary=None):
939 """Perform preparatory calculations common to all Alard&Lupton Tasks.
940
941 Parameters
942 ----------
943 template : `lsst.afw.image.ExposureF`
944 Template exposure, warped to match the science exposure. The
945 variance plane of the template image is modified in place.
946 science : `lsst.afw.image.ExposureF`
947 Science exposure to subtract from the template. The variance plane
948 of the science image is modified in place.
949 visitSummary : `lsst.afw.table.ExposureCatalog`, optional
950 Exposure catalog with external calibrations to be applied. Catalog
951 uses the detector id for the catalog id, sorted on id for fast
952 lookup.
953 """
954 self._validateExposures(template, science)
955 if visitSummary is not None:
956 self._applyExternalCalibrations(science, visitSummary=visitSummary)
957 templateCoverageFraction = checkTemplateIsSufficient(
958 template[science.getBBox()], self.log,
959 requiredTemplateFraction=self.config.requiredTemplateFraction,
960 exceptionMessage="Not attempting subtraction. To force subtraction,"
961 " set config requiredTemplateFraction=0"
962 )
963 self.metadata["templateCoveragePercent"] = 100*templateCoverageFraction
964
965 if self.config.doScaleVariance:
966 # Scale the variance of the template and science images before
967 # convolution, subtraction, or decorrelation so that they have the
968 # correct ratio.
969 templateVarFactor = self.scaleVariance.run(template.maskedImage)
970 sciVarFactor = self.scaleVariance.run(science.maskedImage)
971 self.log.info("Template variance scaling factor: %.2f", templateVarFactor)
972 self.metadata["scaleTemplateVarianceFactor"] = templateVarFactor
973 self.log.info("Science variance scaling factor: %.2f", sciVarFactor)
974 self.metadata["scaleScienceVarianceFactor"] = sciVarFactor
975
976 # Erase existing detection mask planes.
977 # We don't want the detection mask from the science image
978 self.updateMasks(template, science)
979

◆ _renameMaskPlanes()

lsst.ip.diffim.subtractImages.AlardLuptonSubtractTask._renameMaskPlanes ( mask,
maskPlane,
newMaskPlane )
staticprotectedinherited
Rename a mask plane by adding the new name and copying the data.

Parameters
----------
mask : `lsst.afw.image.Mask`
    The mask image to update in place.
maskPlane : `str`
    The name of the existing mask plane to copy.
newMaskPlane : `str`
    The new name of the mask plane that will be added.
    If the mask plane already exists, it will be updated in place.

Definition at line 1025 of file subtractImages.py.

1025 def _renameMaskPlanes(mask, maskPlane, newMaskPlane):
1026 """Rename a mask plane by adding the new name and copying the data.
1027
1028 Parameters
1029 ----------
1030 mask : `lsst.afw.image.Mask`
1031 The mask image to update in place.
1032 maskPlane : `str`
1033 The name of the existing mask plane to copy.
1034 newMaskPlane : `str`
1035 The new name of the mask plane that will be added.
1036 If the mask plane already exists, it will be updated in place.
1037 """
1038 mask.addMaskPlane(newMaskPlane)
1039 originBitMask = mask.getPlaneBitMask(maskPlane)
1040 destinationBitMask = mask.getPlaneBitMask(newMaskPlane)
1041 mask.array |= ((mask.array & originBitMask) > 0)*destinationBitMask
1042

◆ _sourceSelector()

lsst.ip.diffim.subtractImages.AlardLuptonSubtractTask._sourceSelector ( self,
sources,
mask )
protectedinherited
Select sources from a catalog that meet the selection criteria.

Parameters
----------
sources : `lsst.afw.table.SourceCatalog`
    Input source catalog to select sources from.
mask : `lsst.afw.image.Mask`
    The image mask plane to use to reject sources
    based on their location on the ccd.

Returns
-------
selectSources : `lsst.afw.table.SourceCatalog`
    The input source catalog, with flagged and low signal-to-noise
    sources removed.

Raises
------
RuntimeError
    If there are too few sources to compute the PSF matching kernel
    remaining after source selection.

Definition at line 844 of file subtractImages.py.

844 def _sourceSelector(self, sources, mask):
845 """Select sources from a catalog that meet the selection criteria.
846
847 Parameters
848 ----------
849 sources : `lsst.afw.table.SourceCatalog`
850 Input source catalog to select sources from.
851 mask : `lsst.afw.image.Mask`
852 The image mask plane to use to reject sources
853 based on their location on the ccd.
854
855 Returns
856 -------
857 selectSources : `lsst.afw.table.SourceCatalog`
858 The input source catalog, with flagged and low signal-to-noise
859 sources removed.
860
861 Raises
862 ------
863 RuntimeError
864 If there are too few sources to compute the PSF matching kernel
865 remaining after source selection.
866 """
867
868 selected = self.sourceSelector.selectSources(sources).selected
869 nInitialSelected = np.count_nonzero(selected)
870 selected *= self._checkMask(mask, sources, self.config.excludeMaskPlanes)
871 nSelected = np.count_nonzero(selected)
872 self.log.info("Rejecting %i candidate sources: an excluded template mask plane is set.",
873 nInitialSelected - nSelected)
874 selectSources = sources[selected].copy(deep=True)
875 # Trim selectSources if they exceed ``maxKernelSources``.
876 # Keep the highest signal-to-noise sources of those selected.
877 if (len(selectSources) > self.config.maxKernelSources) & (self.config.maxKernelSources > 0):
878 signalToNoise = selectSources.getPsfInstFlux()/selectSources.getPsfInstFluxErr()
879 indices = np.argsort(signalToNoise)
880 indices = indices[-self.config.maxKernelSources:]
881 selected = np.zeros(len(selectSources), dtype=bool)
882 selected[indices] = True
883 selectSources = selectSources[selected].copy(deep=True)
884
885 self.log.info("%i/%i=%.1f%% of sources selected for PSF matching from the input catalog",
886 len(selectSources), len(sources), 100*len(selectSources)/len(sources))
887 if len(selectSources) < self.config.minKernelSources:
888 self.log.error("Too few sources to calculate the PSF matching kernel: "
889 "%i selected but %i needed for the calculation.",
890 len(selectSources), self.config.minKernelSources)
891 if not self.config.allowKernelSourceDetection:
892 raise RuntimeError("Cannot compute PSF matching kernel: too few sources selected.")
893 self.metadata["nPsfSources"] = len(selectSources)
894
895 return selectSources
896

◆ _validateExposures()

lsst.ip.diffim.subtractImages.AlardLuptonSubtractTask._validateExposures ( template,
science )
staticprotectedinherited
Check that the WCS of the two Exposures match, the template bbox
contains the science bbox, and that the bands match.

Parameters
----------
template : `lsst.afw.image.ExposureF`
    Template exposure, warped to match the science exposure.
science : `lsst.afw.image.ExposureF`
    Science exposure to subtract from the template.

Raises
------
AssertionError
    Raised if the WCS of the template is not equal to the science WCS,
    if the science image is not fully contained in the template
    bounding box, or if the bands do not match.

Definition at line 770 of file subtractImages.py.

770 def _validateExposures(template, science):
771 """Check that the WCS of the two Exposures match, the template bbox
772 contains the science bbox, and that the bands match.
773
774 Parameters
775 ----------
776 template : `lsst.afw.image.ExposureF`
777 Template exposure, warped to match the science exposure.
778 science : `lsst.afw.image.ExposureF`
779 Science exposure to subtract from the template.
780
781 Raises
782 ------
783 AssertionError
784 Raised if the WCS of the template is not equal to the science WCS,
785 if the science image is not fully contained in the template
786 bounding box, or if the bands do not match.
787 """
788 assert template.wcs == science.wcs, \
789 "Template and science exposure WCS are not identical."
790 templateBBox = template.getBBox()
791 scienceBBox = science.getBBox()
792 assert science.filter.bandLabel == template.filter.bandLabel, \
793 "Science and template exposures have different bands: %s, %s" % \
794 (science.filter, template.filter)
795
796 assert templateBBox.contains(scienceBBox), \
797 "Template bbox does not contain all of the science image."
798

◆ computeImageMetrics()

lsst.ip.diffim.subtractImages.AlardLuptonSubtractTask.computeImageMetrics ( self,
science,
difference,
stars )
inherited
Compute quality metrics (saved to the task metadata) on the
difference image, at the locations of detected stars on the science
image. This restricts the metric to locations that should be
well-subtracted.

Parameters
----------
science : `lsst.afw.image.ExposureF`
Science exposure that was subtracted.
difference : `lsst.afw.image.ExposureF`
Result of subtracting template and science.
stars : `lsst.afw.table.SourceCatalog`
Good calibration sources detected on science image; these
footprints are what the metrics are computed on.

Notes
-----
The task metadata does not include docstrings, so descriptions of the
computed metrics are given here:

differenceFootprintRatioMean
Mean of the ratio of the absolute value of the difference image
(with the mean absolute value of the sky regions on the difference
image removed) to the science image, computed in the footprints
of stars detected on the science image (the sums below are of the
pixels in each star or sky footprint):
:math:`\mathrm{mean}_{footprints}((\sum |difference| -
\mathrm{mean}(\sum |difference_{sky}|)) / \sum science)`
differenceFootprintRatioStdev
Standard Deviation across footprints of the above ratio.
differenceFootprintSkyRatioMean
Mean of the ratio of the absolute value of sky source regions on
the difference image to the science image (the sum below is of the
pixels in each sky source footprint):
:math:`\mathrm{mean}_{footprints}(\sum |difference_{sky}| / \sum science_{sky})`
differenceFootprintSkyRatioStdev
Standard Deivation across footprints of the above sky ratio.

Definition at line 471 of file subtractImages.py.

471 def computeImageMetrics(self, science, difference, stars):
472 r"""Compute quality metrics (saved to the task metadata) on the
473 difference image, at the locations of detected stars on the science
474 image. This restricts the metric to locations that should be
475 well-subtracted.
476
477 Parameters
478 ----------
479 science : `lsst.afw.image.ExposureF`
480 Science exposure that was subtracted.
481 difference : `lsst.afw.image.ExposureF`
482 Result of subtracting template and science.
483 stars : `lsst.afw.table.SourceCatalog`
484 Good calibration sources detected on science image; these
485 footprints are what the metrics are computed on.
486
487 Notes
488 -----
489 The task metadata does not include docstrings, so descriptions of the
490 computed metrics are given here:
491
492 differenceFootprintRatioMean
493 Mean of the ratio of the absolute value of the difference image
494 (with the mean absolute value of the sky regions on the difference
495 image removed) to the science image, computed in the footprints
496 of stars detected on the science image (the sums below are of the
497 pixels in each star or sky footprint):
498 :math:`\mathrm{mean}_{footprints}((\sum |difference| -
499 \mathrm{mean}(\sum |difference_{sky}|)) / \sum science)`
500 differenceFootprintRatioStdev
501 Standard Deviation across footprints of the above ratio.
502 differenceFootprintSkyRatioMean
503 Mean of the ratio of the absolute value of sky source regions on
504 the difference image to the science image (the sum below is of the
505 pixels in each sky source footprint):
506 :math:`\mathrm{mean}_{footprints}(\sum |difference_{sky}| / \sum science_{sky})`
507 differenceFootprintSkyRatioStdev
508 Standard Deivation across footprints of the above sky ratio.
509 """
510 def footprint_mean(sources, sky=0):
511 """Compute ratio of the absolute value of the diffim to the science
512 image, within each source footprint, subtracting the sky from the
513 diffim values if provided.
514 """
515 n = len(sources)
516 science_footprints = np.zeros(n)
517 difference_footprints = np.zeros(n)
518 ratio = np.zeros(n)
519 for i, record in enumerate(sources):
520 footprint = record.getFootprint()
521 heavy = lsst.afw.detection.makeHeavyFootprint(footprint, science.maskedImage)
522 heavy_diff = lsst.afw.detection.makeHeavyFootprint(footprint, difference.maskedImage)
523 science_footprints[i] = abs(heavy.getImageArray()).mean()
524 difference_footprints[i] = abs(heavy_diff.getImageArray()).mean()
525 ratio[i] = abs((difference_footprints[i] - sky)/science_footprints[i])
526 return science_footprints, difference_footprints, ratio
527
528 sky = stars["sky_source"]
529 sky_science, sky_difference, sky_ratio = footprint_mean(stars[sky])
530 science_footprints, difference_footprints, ratio = footprint_mean(stars[~sky], sky_difference.mean())
531
532 self.metadata["differenceFootprintRatioMean"] = ratio.mean()
533 self.metadata["differenceFootprintRatioStdev"] = ratio.std()
534 self.metadata["differenceFootprintSkyRatioMean"] = sky_ratio.mean()
535 self.metadata["differenceFootprintSkyRatioStdev"] = sky_ratio.std()
536 self.log.info("Mean, stdev of ratio of difference to science "
537 "pixels in star footprints: %5.4f, %5.4f",
538 self.metadata["differenceFootprintRatioMean"],
539 self.metadata["differenceFootprintRatioStdev"])
540
HeavyFootprint< ImagePixelT, MaskPixelT, VariancePixelT > makeHeavyFootprint(Footprint const &foot, lsst::afw::image::MaskedImage< ImagePixelT, MaskPixelT, VariancePixelT > const &img, HeavyFootprintCtrl const *ctrl=nullptr)
Create a HeavyFootprint with footprint defined by the given Footprint and pixel values from the given...

◆ finalize()

lsst.ip.diffim.subtractImages.AlardLuptonSubtractTask.finalize ( self,
template,
science,
difference,
kernel,
templateMatched = True,
preConvMode = False,
preConvKernel = None,
spatiallyVarying = False )
inherited
Decorrelate the difference image to undo the noise correlations
caused by convolution.

Parameters
----------
template : `lsst.afw.image.ExposureF`
    Template exposure, warped to match the science exposure.
science : `lsst.afw.image.ExposureF`
    Science exposure to subtract from the template.
difference : `lsst.afw.image.ExposureF`
    Result of subtracting template and science.
kernel : `lsst.afw.math.Kernel`
    An (optionally spatially-varying) PSF matching kernel
templateMatched : `bool`, optional
    Was the template PSF-matched to the science image?
preConvMode : `bool`, optional
    Was the science image preconvolved with its own PSF
    before PSF matching the template?
preConvKernel : `lsst.afw.detection.Psf`, optional
    If not `None`, then the science image was pre-convolved with
    (the reflection of) this kernel. Must be normalized to sum to 1.
spatiallyVarying : `bool`, optional
    Compute the decorrelation kernel spatially varying across the image?

Returns
-------
correctedExposure : `lsst.afw.image.ExposureF`
    The decorrelated image difference.

Definition at line 680 of file subtractImages.py.

684 spatiallyVarying=False):
685 """Decorrelate the difference image to undo the noise correlations
686 caused by convolution.
687
688 Parameters
689 ----------
690 template : `lsst.afw.image.ExposureF`
691 Template exposure, warped to match the science exposure.
692 science : `lsst.afw.image.ExposureF`
693 Science exposure to subtract from the template.
694 difference : `lsst.afw.image.ExposureF`
695 Result of subtracting template and science.
696 kernel : `lsst.afw.math.Kernel`
697 An (optionally spatially-varying) PSF matching kernel
698 templateMatched : `bool`, optional
699 Was the template PSF-matched to the science image?
700 preConvMode : `bool`, optional
701 Was the science image preconvolved with its own PSF
702 before PSF matching the template?
703 preConvKernel : `lsst.afw.detection.Psf`, optional
704 If not `None`, then the science image was pre-convolved with
705 (the reflection of) this kernel. Must be normalized to sum to 1.
706 spatiallyVarying : `bool`, optional
707 Compute the decorrelation kernel spatially varying across the image?
708
709 Returns
710 -------
711 correctedExposure : `lsst.afw.image.ExposureF`
712 The decorrelated image difference.
713 """
714 if self.config.doDecorrelation:
715 self.log.info("Decorrelating image difference.")
716 # We have cleared the template mask plane, so copy the mask plane of
717 # the image difference so that we can calculate correct statistics
718 # during decorrelation
719 correctedExposure = self.decorrelate.run(science, template[science.getBBox()], difference, kernel,
720 templateMatched=templateMatched,
721 preConvMode=preConvMode,
722 preConvKernel=preConvKernel,
723 spatiallyVarying=spatiallyVarying).correctedExposure
724 else:
725 self.log.info("NOT decorrelating image difference.")
726 correctedExposure = difference
727 return correctedExposure
728

◆ run()

lsst.ip.diffim.subtractImages.AlardLuptonPreconvolveSubtractTask.run ( self,
template,
science,
sources,
visitSummary = None )
Preconvolve the science image with its own PSF,
convolve the template image with a PSF-matching kernel and subtract
from the preconvolved science image.

Parameters
----------
template : `lsst.afw.image.ExposureF`
    The template image, which has previously been warped to the science
    image. The template bbox will be padded by a few pixels compared to
    the science bbox.
science : `lsst.afw.image.ExposureF`
    The science exposure.
sources : `lsst.afw.table.SourceCatalog`
    Identified sources on the science exposure. This catalog is used to
    select sources in order to perform the AL PSF matching on stamp
    images around them.
visitSummary : `lsst.afw.table.ExposureCatalog`, optional
    Exposure catalog with complete external calibrations. Catalog uses
    the detector id for the catalog id, sorted on id for fast lookup.

Returns
-------
results : `lsst.pipe.base.Struct`
    ``scoreExposure`` : `lsst.afw.image.ExposureF`
        Result of subtracting the convolved template and science
        images. Attached PSF is that of the original science image.
    ``matchedTemplate`` : `lsst.afw.image.ExposureF`
        Warped and PSF-matched template exposure. Attached PSF is that
        of the original science image.
    ``matchedScience`` : `lsst.afw.image.ExposureF`
        The science exposure after convolving with its own PSF.
        Attached PSF is that of the original science image.
    ``backgroundModel`` : `lsst.afw.math.Function2D`
        Background model that was fit while solving for the
        PSF-matching kernel
    ``psfMatchingKernel`` : `lsst.afw.math.Kernel`
        Final kernel used to PSF-match the template to the science
        image.

Reimplemented from lsst.ip.diffim.subtractImages.AlardLuptonSubtractTask.

Definition at line 1079 of file subtractImages.py.

1079 def run(self, template, science, sources, visitSummary=None):
1080 """Preconvolve the science image with its own PSF,
1081 convolve the template image with a PSF-matching kernel and subtract
1082 from the preconvolved science image.
1083
1084 Parameters
1085 ----------
1086 template : `lsst.afw.image.ExposureF`
1087 The template image, which has previously been warped to the science
1088 image. The template bbox will be padded by a few pixels compared to
1089 the science bbox.
1090 science : `lsst.afw.image.ExposureF`
1091 The science exposure.
1092 sources : `lsst.afw.table.SourceCatalog`
1093 Identified sources on the science exposure. This catalog is used to
1094 select sources in order to perform the AL PSF matching on stamp
1095 images around them.
1096 visitSummary : `lsst.afw.table.ExposureCatalog`, optional
1097 Exposure catalog with complete external calibrations. Catalog uses
1098 the detector id for the catalog id, sorted on id for fast lookup.
1099
1100 Returns
1101 -------
1102 results : `lsst.pipe.base.Struct`
1103 ``scoreExposure`` : `lsst.afw.image.ExposureF`
1104 Result of subtracting the convolved template and science
1105 images. Attached PSF is that of the original science image.
1106 ``matchedTemplate`` : `lsst.afw.image.ExposureF`
1107 Warped and PSF-matched template exposure. Attached PSF is that
1108 of the original science image.
1109 ``matchedScience`` : `lsst.afw.image.ExposureF`
1110 The science exposure after convolving with its own PSF.
1111 Attached PSF is that of the original science image.
1112 ``backgroundModel`` : `lsst.afw.math.Function2D`
1113 Background model that was fit while solving for the
1114 PSF-matching kernel
1115 ``psfMatchingKernel`` : `lsst.afw.math.Kernel`
1116 Final kernel used to PSF-match the template to the science
1117 image.
1118 """
1119 self._prepareInputs(template, science, visitSummary=visitSummary)
1120
1121 # TODO: DM-37212 we need to mirror the kernel in order to get correct cross correlation
1122 scienceKernel = science.psf.getKernel()
1123 matchedScience = self._convolveExposure(science, scienceKernel, self.convolutionControl,
1124 interpolateBadMaskPlanes=True)
1125 self.metadata["convolvedExposure"] = "Preconvolution"
1126 try:
1127 selectSources = self._sourceSelector(sources, matchedScience.mask)
1128 subtractResults = self.runPreconvolve(template, science, matchedScience,
1129 selectSources, scienceKernel)
1130
1131 except (RuntimeError, lsst.pex.exceptions.Exception) as e:
1132 self.log.warning("Failed to match template. Checking coverage")
1133 # Raise NoWorkFound if template fraction is insufficient
1134 checkTemplateIsSufficient(template[science.getBBox()], self.log,
1135 self.config.minTemplateFractionForExpectedSuccess,
1136 exceptionMessage="Template coverage lower than expected to succeed."
1137 f" Failure is tolerable: {e}")
1138 # checkTemplateIsSufficient did not raise NoWorkFound, so raise original exception
1139 raise e
1140
1141 return subtractResults
1142
Provides consistent interface for LSST exceptions.
Definition Exception.h:107

◆ runConvolveScience()

lsst.ip.diffim.subtractImages.AlardLuptonSubtractTask.runConvolveScience ( self,
template,
science,
selectSources )
inherited
Convolve the science image with a PSF-matching kernel and subtract
the template image.

Parameters
----------
template : `lsst.afw.image.ExposureF`
    Template exposure, warped to match the science exposure.
science : `lsst.afw.image.ExposureF`
    Science exposure to subtract from the template.
selectSources : `lsst.afw.table.SourceCatalog`
    Identified sources on the science exposure. This catalog is used to
    select sources in order to perform the AL PSF matching on stamp
    images around them.

Returns
-------
results : `lsst.pipe.base.Struct`

    ``difference`` : `lsst.afw.image.ExposureF`
        Result of subtracting template and science.
    ``matchedTemplate`` : `lsst.afw.image.ExposureF`
        Warped template exposure. Note that in this case, the template
        is not PSF-matched to the science image.
    ``backgroundModel`` : `lsst.afw.math.Function2D`
        Background model that was fit while solving for the PSF-matching kernel
    ``psfMatchingKernel`` : `lsst.afw.math.Kernel`
       Kernel used to PSF-match the science image to the template.

Definition at line 613 of file subtractImages.py.

613 def runConvolveScience(self, template, science, selectSources):
614 """Convolve the science image with a PSF-matching kernel and subtract
615 the template image.
616
617 Parameters
618 ----------
619 template : `lsst.afw.image.ExposureF`
620 Template exposure, warped to match the science exposure.
621 science : `lsst.afw.image.ExposureF`
622 Science exposure to subtract from the template.
623 selectSources : `lsst.afw.table.SourceCatalog`
624 Identified sources on the science exposure. This catalog is used to
625 select sources in order to perform the AL PSF matching on stamp
626 images around them.
627
628 Returns
629 -------
630 results : `lsst.pipe.base.Struct`
631
632 ``difference`` : `lsst.afw.image.ExposureF`
633 Result of subtracting template and science.
634 ``matchedTemplate`` : `lsst.afw.image.ExposureF`
635 Warped template exposure. Note that in this case, the template
636 is not PSF-matched to the science image.
637 ``backgroundModel`` : `lsst.afw.math.Function2D`
638 Background model that was fit while solving for the PSF-matching kernel
639 ``psfMatchingKernel`` : `lsst.afw.math.Kernel`
640 Kernel used to PSF-match the science image to the template.
641 """
642 bbox = science.getBBox()
643 kernelSources = self.makeKernel.selectKernelSources(science, template,
644 candidateList=selectSources,
645 preconvolved=False)
646 kernelResult = self.makeKernel.run(science, template, kernelSources,
647 preconvolved=False)
648 modelParams = kernelResult.backgroundModel.getParameters()
649 # We must invert the background model if the matching kernel is solved for the science image.
650 kernelResult.backgroundModel.setParameters([-p for p in modelParams])
651
652 kernelImage = lsst.afw.image.ImageD(kernelResult.psfMatchingKernel.getDimensions())
653 norm = kernelResult.psfMatchingKernel.computeImage(kernelImage, doNormalize=False)
654
655 matchedScience = self._convolveExposure(science, kernelResult.psfMatchingKernel,
656 self.convolutionControl,
657 psf=template.psf)
658
659 # Place back on native photometric scale
660 matchedScience.maskedImage /= norm
661 matchedTemplate = template.clone()[bbox]
662 matchedTemplate.maskedImage /= norm
663 matchedTemplate.setPhotoCalib(science.photoCalib)
664
665 difference = _subtractImages(matchedScience, matchedTemplate,
666 backgroundModel=(kernelResult.backgroundModel
667 if self.config.doSubtractBackground else None))
668
669 correctedExposure = self.finalize(template, science, difference,
670 kernelResult.psfMatchingKernel,
671 templateMatched=False)
672
673 return lsst.pipe.base.Struct(difference=correctedExposure,
674 matchedTemplate=matchedTemplate,
675 matchedScience=matchedScience,
676 backgroundModel=kernelResult.backgroundModel,
677 psfMatchingKernel=kernelResult.psfMatchingKernel,
678 kernelSources=kernelSources)
679

◆ runConvolveTemplate()

lsst.ip.diffim.subtractImages.AlardLuptonSubtractTask.runConvolveTemplate ( self,
template,
science,
selectSources )
inherited
Convolve the template image with a PSF-matching kernel and subtract
from the science image.

Parameters
----------
template : `lsst.afw.image.ExposureF`
    Template exposure, warped to match the science exposure.
science : `lsst.afw.image.ExposureF`
    Science exposure to subtract from the template.
selectSources : `lsst.afw.table.SourceCatalog`
    Identified sources on the science exposure. This catalog is used to
    select sources in order to perform the AL PSF matching on stamp
    images around them.

Returns
-------
results : `lsst.pipe.base.Struct`

    ``difference`` : `lsst.afw.image.ExposureF`
        Result of subtracting template and science.
    ``matchedTemplate`` : `lsst.afw.image.ExposureF`
        Warped and PSF-matched template exposure.
    ``backgroundModel`` : `lsst.afw.math.Function2D`
        Background model that was fit while solving for the PSF-matching kernel
    ``psfMatchingKernel`` : `lsst.afw.math.Kernel`
        Kernel used to PSF-match the template to the science image.

Definition at line 541 of file subtractImages.py.

541 def runConvolveTemplate(self, template, science, selectSources):
542 """Convolve the template image with a PSF-matching kernel and subtract
543 from the science image.
544
545 Parameters
546 ----------
547 template : `lsst.afw.image.ExposureF`
548 Template exposure, warped to match the science exposure.
549 science : `lsst.afw.image.ExposureF`
550 Science exposure to subtract from the template.
551 selectSources : `lsst.afw.table.SourceCatalog`
552 Identified sources on the science exposure. This catalog is used to
553 select sources in order to perform the AL PSF matching on stamp
554 images around them.
555
556 Returns
557 -------
558 results : `lsst.pipe.base.Struct`
559
560 ``difference`` : `lsst.afw.image.ExposureF`
561 Result of subtracting template and science.
562 ``matchedTemplate`` : `lsst.afw.image.ExposureF`
563 Warped and PSF-matched template exposure.
564 ``backgroundModel`` : `lsst.afw.math.Function2D`
565 Background model that was fit while solving for the PSF-matching kernel
566 ``psfMatchingKernel`` : `lsst.afw.math.Kernel`
567 Kernel used to PSF-match the template to the science image.
568 """
569 try:
570 kernelSources = self.makeKernel.selectKernelSources(template, science,
571 candidateList=selectSources,
572 preconvolved=False)
573 kernelResult = self.makeKernel.run(template, science, kernelSources,
574 preconvolved=False)
575 except Exception as e:
576 if self.config.allowKernelSourceDetection:
577 self.log.warning("Error encountered trying to construct the matching kernel"
578 f" Running source detection and retrying. {e}")
579 kernelSize = self.makeKernel.makeKernelBasisList(
580 self.templatePsfSize, self.sciencePsfSize)[0].getWidth()
581 sigmaToFwhm = 2*np.log(2*np.sqrt(2))
582 candidateList = self.makeKernel.makeCandidateList(template, science, kernelSize,
583 candidateList=None,
584 sigma=self.sciencePsfSize/sigmaToFwhm)
585 kernelSources = self.makeKernel.selectKernelSources(template, science,
586 candidateList=candidateList,
587 preconvolved=False)
588 kernelResult = self.makeKernel.run(template, science, kernelSources,
589 preconvolved=False)
590 else:
591 raise e
592
593 matchedTemplate = self._convolveExposure(template, kernelResult.psfMatchingKernel,
594 self.convolutionControl,
595 bbox=science.getBBox(),
596 psf=science.psf,
597 photoCalib=science.photoCalib)
598
599 difference = _subtractImages(science, matchedTemplate,
600 backgroundModel=(kernelResult.backgroundModel
601 if self.config.doSubtractBackground else None))
602 correctedExposure = self.finalize(template, science, difference,
603 kernelResult.psfMatchingKernel,
604 templateMatched=True)
605
606 return lsst.pipe.base.Struct(difference=correctedExposure,
607 matchedTemplate=matchedTemplate,
608 matchedScience=science,
609 backgroundModel=kernelResult.backgroundModel,
610 psfMatchingKernel=kernelResult.psfMatchingKernel,
611 kernelSources=kernelSources)
612

◆ runPreconvolve()

lsst.ip.diffim.subtractImages.AlardLuptonPreconvolveSubtractTask.runPreconvolve ( self,
template,
science,
matchedScience,
selectSources,
preConvKernel )
Convolve the science image with its own PSF, then convolve the
template with a matching kernel and subtract to form the Score
exposure.

Parameters
----------
template : `lsst.afw.image.ExposureF`
    Template exposure, warped to match the science exposure.
science : `lsst.afw.image.ExposureF`
    Science exposure to subtract from the template.
matchedScience : `lsst.afw.image.ExposureF`
    The science exposure, convolved with the reflection of its own PSF.
selectSources : `lsst.afw.table.SourceCatalog`
    Identified sources on the science exposure. This catalog is used to
    select sources in order to perform the AL PSF matching on stamp
    images around them.
preConvKernel : `lsst.afw.math.Kernel`
    The reflection of the kernel that was used to preconvolve the
    `science` exposure. Must be normalized to sum to 1.

Returns
-------
results : `lsst.pipe.base.Struct`

    ``scoreExposure`` : `lsst.afw.image.ExposureF`
        Result of subtracting the convolved template and science
        images. Attached PSF is that of the original science image.
    ``matchedTemplate`` : `lsst.afw.image.ExposureF`
        Warped and PSF-matched template exposure. Attached PSF is that
        of the original science image.
    ``matchedScience`` : `lsst.afw.image.ExposureF`
        The science exposure after convolving with its own PSF.
        Attached PSF is that of the original science image.
    ``backgroundModel`` : `lsst.afw.math.Function2D`
        Background model that was fit while solving for the
        PSF-matching kernel
    ``psfMatchingKernel`` : `lsst.afw.math.Kernel`
        Final kernel used to PSF-match the template to the science
        image.

Definition at line 1143 of file subtractImages.py.

1143 def runPreconvolve(self, template, science, matchedScience, selectSources, preConvKernel):
1144 """Convolve the science image with its own PSF, then convolve the
1145 template with a matching kernel and subtract to form the Score
1146 exposure.
1147
1148 Parameters
1149 ----------
1150 template : `lsst.afw.image.ExposureF`
1151 Template exposure, warped to match the science exposure.
1152 science : `lsst.afw.image.ExposureF`
1153 Science exposure to subtract from the template.
1154 matchedScience : `lsst.afw.image.ExposureF`
1155 The science exposure, convolved with the reflection of its own PSF.
1156 selectSources : `lsst.afw.table.SourceCatalog`
1157 Identified sources on the science exposure. This catalog is used to
1158 select sources in order to perform the AL PSF matching on stamp
1159 images around them.
1160 preConvKernel : `lsst.afw.math.Kernel`
1161 The reflection of the kernel that was used to preconvolve the
1162 `science` exposure. Must be normalized to sum to 1.
1163
1164 Returns
1165 -------
1166 results : `lsst.pipe.base.Struct`
1167
1168 ``scoreExposure`` : `lsst.afw.image.ExposureF`
1169 Result of subtracting the convolved template and science
1170 images. Attached PSF is that of the original science image.
1171 ``matchedTemplate`` : `lsst.afw.image.ExposureF`
1172 Warped and PSF-matched template exposure. Attached PSF is that
1173 of the original science image.
1174 ``matchedScience`` : `lsst.afw.image.ExposureF`
1175 The science exposure after convolving with its own PSF.
1176 Attached PSF is that of the original science image.
1177 ``backgroundModel`` : `lsst.afw.math.Function2D`
1178 Background model that was fit while solving for the
1179 PSF-matching kernel
1180 ``psfMatchingKernel`` : `lsst.afw.math.Kernel`
1181 Final kernel used to PSF-match the template to the science
1182 image.
1183 """
1184 bbox = science.getBBox()
1185 innerBBox = preConvKernel.shrinkBBox(bbox)
1186
1187 kernelSources = self.makeKernel.selectKernelSources(template[innerBBox], matchedScience[innerBBox],
1188 candidateList=selectSources,
1189 preconvolved=True)
1190 kernelResult = self.makeKernel.run(template[innerBBox], matchedScience[innerBBox], kernelSources,
1191 preconvolved=True)
1192
1193 matchedTemplate = self._convolveExposure(template, kernelResult.psfMatchingKernel,
1194 self.convolutionControl,
1195 bbox=bbox,
1196 psf=science.psf,
1197 interpolateBadMaskPlanes=True,
1198 photoCalib=science.photoCalib)
1199 score = _subtractImages(matchedScience, matchedTemplate,
1200 backgroundModel=(kernelResult.backgroundModel
1201 if self.config.doSubtractBackground else None))
1202 correctedScore = self.finalize(template[bbox], science, score,
1203 kernelResult.psfMatchingKernel,
1204 templateMatched=True, preConvMode=True,
1205 preConvKernel=preConvKernel)
1206
1207 return lsst.pipe.base.Struct(scoreExposure=correctedScore,
1208 matchedTemplate=matchedTemplate,
1209 matchedScience=matchedScience,
1210 backgroundModel=kernelResult.backgroundModel,
1211 psfMatchingKernel=kernelResult.psfMatchingKernel,
1212 kernelSources=kernelSources)
1213
1214

◆ updateMasks()

lsst.ip.diffim.subtractImages.AlardLuptonSubtractTask.updateMasks ( self,
template,
science )
inherited
Update the science and template mask planes before differencing.

Parameters
----------
template : `lsst.afw.image.Exposure`
    Template exposure, warped to match the science exposure.
    The template mask planes will be erased, except for a few specified
    in the task config.
science : `lsst.afw.image.Exposure`
    Science exposure to subtract from the template.
    The DETECTED and DETECTED_NEGATIVE mask planes of the science image
    will be erased.

Definition at line 980 of file subtractImages.py.

980 def updateMasks(self, template, science):
981 """Update the science and template mask planes before differencing.
982
983 Parameters
984 ----------
985 template : `lsst.afw.image.Exposure`
986 Template exposure, warped to match the science exposure.
987 The template mask planes will be erased, except for a few specified
988 in the task config.
989 science : `lsst.afw.image.Exposure`
990 Science exposure to subtract from the template.
991 The DETECTED and DETECTED_NEGATIVE mask planes of the science image
992 will be erased.
993 """
994 self._clearMask(science.mask, clearMaskPlanes=["DETECTED", "DETECTED_NEGATIVE"])
995
996 # We will clear ALL template mask planes, except for those specified
997 # via the `preserveTemplateMask` config. Mask planes specified via
998 # the `renameTemplateMask` config will be copied to new planes with
999 # "_TEMPLATE" appended to their names, and the original mask plane will
1000 # be cleared.
1001 clearMaskPlanes = [mp for mp in template.mask.getMaskPlaneDict().keys()
1002 if mp not in self.config.preserveTemplateMask]
1003 renameMaskPlanes = [mp for mp in self.config.renameTemplateMask
1004 if mp in template.mask.getMaskPlaneDict().keys()]
1005
1006 # propagate the mask plane related to Fake source injection
1007 # NOTE: the fake source injection sets FAKE plane, but it should be INJECTED
1008 # NOTE: This can be removed in DM-40796
1009 if "FAKE" in science.mask.getMaskPlaneDict().keys():
1010 self.log.info("Adding injected mask plane to science image")
1011 self._renameMaskPlanes(science.mask, "FAKE", "INJECTED")
1012 if "FAKE" in template.mask.getMaskPlaneDict().keys():
1013 self.log.info("Adding injected mask plane to template image")
1014 self._renameMaskPlanes(template.mask, "FAKE", "INJECTED_TEMPLATE")
1015 if "INJECTED" in renameMaskPlanes:
1016 renameMaskPlanes.remove("INJECTED")
1017 if "INJECTED_TEMPLATE" in clearMaskPlanes:
1018 clearMaskPlanes.remove("INJECTED_TEMPLATE")
1019
1020 for maskPlane in renameMaskPlanes:
1021 self._renameMaskPlanes(template.mask, maskPlane, maskPlane + "_TEMPLATE")
1022 self._clearMask(template.mask, clearMaskPlanes=clearMaskPlanes)
1023

Member Data Documentation

◆ _DefaultName

str lsst.ip.diffim.subtractImages.AlardLuptonSubtractTask._DefaultName = "alardLuptonSubtract"
staticprotectedinherited

Definition at line 274 of file subtractImages.py.

◆ ConfigClass

lsst.ip.diffim.subtractImages.AlardLuptonSubtractTask.ConfigClass = AlardLuptonSubtractConfig
staticinherited

Definition at line 273 of file subtractImages.py.

◆ convolutionControl

lsst.ip.diffim.subtractImages.AlardLuptonSubtractTask.convolutionControl = lsst.afw.math.ConvolutionControl()
inherited

Definition at line 284 of file subtractImages.py.

◆ log

lsst.ip.diffim.subtractImages.AlardLuptonSubtractTask.log = self.runConvolveScience(template, science, selectSources)
inherited

Definition at line 460 of file subtractImages.py.

◆ sciencePsfSize

lsst.ip.diffim.subtractImages.AlardLuptonSubtractTask.sciencePsfSize = getPsfFwhm(science.psf)
inherited

Definition at line 392 of file subtractImages.py.

◆ templatePsfSize

lsst.ip.diffim.subtractImages.AlardLuptonSubtractTask.templatePsfSize = getPsfFwhm(template.psf)
inherited

Definition at line 391 of file subtractImages.py.


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