Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

DM-35175: Fix the variance plane calculation when the science image is convolved #222

Merged
merged 2 commits into from
Jun 25, 2022
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
50 changes: 39 additions & 11 deletions python/lsst/ip/diffim/subtractImages.py
Original file line number Diff line number Diff line change
Expand Up @@ -163,6 +163,7 @@ def __init__(self, **kwargs):
# Normalization is an extra, unnecessary, calculation and will result
# in mis-subtraction of the images if there are calibration errors.
self.convolutionControl.setDoNormalize(False)
self.convolutionControl.setDoCopyEdge(True)

def run(self, template, science, sources):
"""PSF match, subtract, and decorrelate two images.
Expand Down Expand Up @@ -204,10 +205,8 @@ def run(self, template, science, sources):
if self.config.forceCompatibility:
# Compatibility option to maintain old functionality
# This should be removed in the future!
self.log.warning("Running with `config.forceCompatibility=True`")
sources = None
kernelSources = self.makeKernel.selectKernelSources(template, science,
candidateList=sources,
preconvolved=False)
sciencePsfSize = getPsfFwhm(science.psf)
templatePsfSize = getPsfFwhm(template.psf)
self.log.info("Science PSF size: %f", sciencePsfSize)
Expand All @@ -228,12 +227,27 @@ def run(self, template, science, sources):
else:
raise RuntimeError("Cannot handle AlardLuptonSubtract mode: %s", self.config.mode)

if self.config.doScaleVariance and ~self.config.forceCompatibility:
# Scale the variance of the template and science images before
# convolution, subtraction, or decorrelation so that they have the
# correct ratio.
templateVarFactor = self.scaleVariance.run(template.maskedImage)
sciVarFactor = self.scaleVariance.run(science.maskedImage)
self.log.info("Template variance scaling factor: %.2f", templateVarFactor)
self.metadata.add("scaleTemplateVarianceFactor", templateVarFactor)
self.log.info("Science variance scaling factor: %.2f", sciVarFactor)
self.metadata.add("scaleScienceVarianceFactor", sciVarFactor)

kernelSources = self.makeKernel.selectKernelSources(template, science,
candidateList=sources,
preconvolved=False)
if convolveTemplate:
subtractResults = self.runConvolveTemplate(template, science, kernelSources)
else:
subtractResults = self.runConvolveScience(template, science, kernelSources)

if self.config.doScaleVariance:
if self.config.doScaleVariance and self.config.forceCompatibility:
# The old behavior scaled the variance of the final image difference.
diffimVarFactor = self.scaleVariance.run(subtractResults.difference.maskedImage)
self.log.info("Diffim variance scaling factor: %.2f", diffimVarFactor)
self.metadata.add("scaleDiffimVarianceFactor", diffimVarFactor)
Expand Down Expand Up @@ -277,7 +291,8 @@ def runConvolveTemplate(self, template, science, sources):
matchedTemplate = self._convolveExposure(template, kernelResult.psfMatchingKernel,
self.convolutionControl,
bbox=science.getBBox(),
psf=science.psf)
psf=science.psf,
photoCalib=science.getPhotoCalib())
difference = _subtractImages(science, matchedTemplate,
backgroundModel=(kernelResult.backgroundModel
if self.config.doSubtractBackground else None))
Expand All @@ -286,6 +301,7 @@ def runConvolveTemplate(self, template, science, sources):

return lsst.pipe.base.Struct(difference=correctedExposure,
matchedTemplate=matchedTemplate,
matchedScience=science,
backgroundModel=kernelResult.backgroundModel,
psfMatchingKernel=kernelResult.psfMatchingKernel)

Expand Down Expand Up @@ -326,22 +342,29 @@ def runConvolveScience(self, template, science, sources):
# We must invert the background model if the matching kernel is solved for the science image.
kernelResult.backgroundModel.setParameters([-p for p in modelParams])

kernelImage = lsst.afw.image.ImageD(kernelResult.psfMatchingKernel.getDimensions())
norm = kernelResult.psfMatchingKernel.computeImage(kernelImage, doNormalize=False)
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is there any reason one might want to allow doNormalize to be configurable here, as opposed to hard-coded false? I think not but doesn't hurt to think about it a bit.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It's actually important that it be hard-coded as False here, so that we get the flux scale right.


matchedScience = self._convolveExposure(science, kernelResult.psfMatchingKernel,
self.convolutionControl,
psf=template.psf)

difference = _subtractImages(matchedScience, template[science.getBBox()],
# Place back on native photometric scale
matchedScience.maskedImage /= norm
matchedTemplate = template.clone()[science.getBBox()]
matchedTemplate.maskedImage /= norm
matchedTemplate.setPhotoCalib(science.getPhotoCalib())

difference = _subtractImages(matchedScience, matchedTemplate,
backgroundModel=(kernelResult.backgroundModel
if self.config.doSubtractBackground else None))

# Place back on native photometric scale
difference.maskedImage /= kernelResult.psfMatchingKernel.computeImage(
lsst.afw.image.ImageD(kernelResult.psfMatchingKernel.getDimensions()), False)
correctedExposure = self.finalize(template, science, difference, kernelResult.psfMatchingKernel,
templateMatched=False)

return lsst.pipe.base.Struct(difference=correctedExposure,
matchedTemplate=template,
matchedTemplate=matchedTemplate,
matchedScience=matchedScience,
backgroundModel=kernelResult.backgroundModel,
psfMatchingKernel=kernelResult.psfMatchingKernel,)

Expand Down Expand Up @@ -426,7 +449,8 @@ def _validateExposures(template, science):
@staticmethod
def _convolveExposure(exposure, kernel, convolutionControl,
bbox=None,
psf=None):
psf=None,
photoCalib=None):
"""Convolve an exposure with the given kernel.

Parameters
Expand All @@ -441,6 +465,8 @@ def _convolveExposure(exposure, kernel, convolutionControl,
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
-------
Expand All @@ -450,6 +476,8 @@ def _convolveExposure(exposure, kernel, convolutionControl,
convolvedExposure = exposure.clone()
if psf is not None:
convolvedExposure.setPsf(psf)
if photoCalib is not None:
convolvedExposure.setPhotoCalib(photoCalib)
convolvedImage = lsst.afw.image.MaskedImageF(exposure.getBBox())
lsst.afw.math.convolve(convolvedImage, exposure.maskedImage, kernel, convolutionControl)
convolvedExposure.setMaskedImage(convolvedImage)
Expand Down