Skip to content
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.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
92 changes: 48 additions & 44 deletions python/lsst/ip/diffim/makeKernel.py
Original file line number Diff line number Diff line change
Expand Up @@ -105,7 +105,8 @@ def __init__(self, *args, **kwargs):
self.makeSubtask("selectDetection", schema=self.selectSchema)
self.makeSubtask("selectMeasurement", schema=self.selectSchema, algMetadata=self.selectAlgMetadata)

def run(self, template, science, kernelSources, preconvolved=False):
def run(self, template, science, kernelSources, preconvolved=False,
templateFwhmPix=None, scienceFwhmPix=None):
"""Solve for the kernel and background model that best match two
Exposures evaluated at the given source locations.

Expand All @@ -131,27 +132,28 @@ def run(self, template, science, kernelSources, preconvolved=False):
Spatially varying background-matching function.
"""
kernelCellSet = self._buildCellSet(template.maskedImage, science.maskedImage, kernelSources)
# Calling getPsfFwhm on template.psf fails on some rare occasions when
# the template has no input exposures at the average position of the
# stars. So we try getPsfFwhm first on template, and if that fails we
# evaluate the PSF on a grid specified by fwhmExposure* fields.
# To keep consistent definitions for PSF size on the template and
# science images, we use the same method for both.
try:
templateFwhmPix = getPsfFwhm(template.psf)
scienceFwhmPix = getPsfFwhm(science.psf)
except (InvalidParameterError, RangeError):
self.log.debug("Unable to evaluate PSF at the average position. "
"Evaluting PSF on a grid of points."
)
templateFwhmPix = evaluateMeanPsfFwhm(template,
fwhmExposureBuffer=self.config.fwhmExposureBuffer,
fwhmExposureGrid=self.config.fwhmExposureGrid
)
scienceFwhmPix = evaluateMeanPsfFwhm(science,
fwhmExposureBuffer=self.config.fwhmExposureBuffer,
fwhmExposureGrid=self.config.fwhmExposureGrid
)
if (scienceFwhmPix is None) or (templateFwhmPix is None):
# Calling getPsfFwhm on template.psf fails on some rare occasions when
# the template has no input exposures at the average position of the
# stars. So we try getPsfFwhm first on template, and if that fails we
# evaluate the PSF on a grid specified by fwhmExposure* fields.
# To keep consistent definitions for PSF size on the template and
# science images, we use the same method for both.
try:
templateFwhmPix = getPsfFwhm(template.psf)
scienceFwhmPix = getPsfFwhm(science.psf)
except (InvalidParameterError, RangeError):
self.log.debug("Unable to evaluate PSF at the average position. "
"Evaluting PSF on a grid of points."
)
templateFwhmPix = evaluateMeanPsfFwhm(template,
fwhmExposureBuffer=self.config.fwhmExposureBuffer,
fwhmExposureGrid=self.config.fwhmExposureGrid
)
scienceFwhmPix = evaluateMeanPsfFwhm(science,
fwhmExposureBuffer=self.config.fwhmExposureBuffer,
fwhmExposureGrid=self.config.fwhmExposureGrid
)

if preconvolved:
scienceFwhmPix *= np.sqrt(2)
Expand All @@ -163,7 +165,8 @@ def run(self, template, science, kernelSources, preconvolved=False):
backgroundModel=backgroundModel,
)

def selectKernelSources(self, template, science, candidateList=None, preconvolved=False):
def selectKernelSources(self, template, science, candidateList=None, preconvolved=False,
templateFwhmPix=None, scienceFwhmPix=None):
"""Select sources from a list of candidates, and extract footprints.

Parameters
Expand All @@ -182,27 +185,28 @@ def selectKernelSources(self, template, science, candidateList=None, preconvolve
kernelSources : `lsst.afw.table.SourceCatalog`
Kernel candidates with appropriate sized footprints.
"""
# Calling getPsfFwhm on template.psf fails on some rare occasions when
# the template has no input exposures at the average position of the
# stars. So we try getPsfFwhm first on template, and if that fails we
# evaluate the PSF on a grid specified by fwhmExposure* fields.
# To keep consistent definitions for PSF size on the template and
# science images, we use the same method for both.
try:
templateFwhmPix = getPsfFwhm(template.psf)
scienceFwhmPix = getPsfFwhm(science.psf)
except (InvalidParameterError, RangeError):
self.log.debug("Unable to evaluate PSF at the average position. "
"Evaluting PSF on a grid of points."
)
templateFwhmPix = evaluateMeanPsfFwhm(template,
fwhmExposureBuffer=self.config.fwhmExposureBuffer,
fwhmExposureGrid=self.config.fwhmExposureGrid
)
scienceFwhmPix = evaluateMeanPsfFwhm(science,
fwhmExposureBuffer=self.config.fwhmExposureBuffer,
fwhmExposureGrid=self.config.fwhmExposureGrid
)
if (scienceFwhmPix is None) or (templateFwhmPix is None):
# Calling getPsfFwhm on template.psf fails on some rare occasions when
# the template has no input exposures at the average position of the
# stars. So we try getPsfFwhm first on template, and if that fails we
# evaluate the PSF on a grid specified by fwhmExposure* fields.
# To keep consistent definitions for PSF size on the template and
# science images, we use the same method for both.
try:
templateFwhmPix = getPsfFwhm(template.psf)
scienceFwhmPix = getPsfFwhm(science.psf)
except (InvalidParameterError, RangeError):
self.log.debug("Unable to evaluate PSF at the average position. "
"Evaluting PSF on a grid of points."
)
templateFwhmPix = evaluateMeanPsfFwhm(template,
fwhmExposureBuffer=self.config.fwhmExposureBuffer,
fwhmExposureGrid=self.config.fwhmExposureGrid
)
scienceFwhmPix = evaluateMeanPsfFwhm(science,
fwhmExposureBuffer=self.config.fwhmExposureBuffer,
fwhmExposureGrid=self.config.fwhmExposureGrid
)
if preconvolved:
scienceFwhmPix *= np.sqrt(2)
kernelSize = self.makeKernelBasisList(templateFwhmPix, scienceFwhmPix)[0].getWidth()
Expand Down
112 changes: 65 additions & 47 deletions python/lsst/ip/diffim/subtractImages.py
Original file line number Diff line number Diff line change
Expand Up @@ -374,43 +374,10 @@ def run(self, template, science, sources, visitSummary=None):
"""
self._prepareInputs(template, science, visitSummary=visitSummary)

# In the event that getPsfFwhm fails, evaluate the PSF on a grid.
fwhmExposureBuffer = self.config.makeKernel.fwhmExposureBuffer
fwhmExposureGrid = self.config.makeKernel.fwhmExposureGrid

# Calling getPsfFwhm on template.psf fails on some rare occasions when
# the template has no input exposures at the average position of the
# stars. So we try getPsfFwhm first on template, and if that fails we
# evaluate the PSF on a grid specified by fwhmExposure* fields.
# To keep consistent definitions for PSF size on the template and
# science images, we use the same method for both.
# In the try block below, we catch two exceptions:
# 1. InvalidParameterError, in case the point where we are evaluating
# the PSF lands in a gap in the template.
# 2. RangeError, in case the template coverage is so poor that we end
# up near a region with no data.
try:
self.templatePsfSize = getPsfFwhm(template.psf)
self.sciencePsfSize = getPsfFwhm(science.psf)
except (lsst.pex.exceptions.InvalidParameterError, lsst.pex.exceptions.RangeError):
self.log.info("Unable to evaluate PSF at the average position. "
"Evaluting PSF on a grid of points."
)
self.templatePsfSize = evaluateMeanPsfFwhm(template,
fwhmExposureBuffer=fwhmExposureBuffer,
fwhmExposureGrid=fwhmExposureGrid
)
self.sciencePsfSize = evaluateMeanPsfFwhm(science,
fwhmExposureBuffer=fwhmExposureBuffer,
fwhmExposureGrid=fwhmExposureGrid
)
self.log.info("Science PSF FWHM: %f pixels", self.sciencePsfSize)
self.log.info("Template PSF FWHM: %f pixels", self.templatePsfSize)
self.metadata["sciencePsfSize"] = self.sciencePsfSize
self.metadata["templatePsfSize"] = self.templatePsfSize

# Calculate estimated image depths, i.e., limiting magnitudes
maglim_science = self._calculateMagLim(science, fallbackPsfSize=self.sciencePsfSize)
if np.isnan(maglim_science):
self.log.warning("Limiting magnitude of the science image is NaN!")
fluxlim_science = (maglim_science*u.ABmag).to_value(u.nJy)
maglim_template = self._calculateMagLim(template, fallbackPsfSize=self.templatePsfSize)
if np.isnan(maglim_template):
Expand All @@ -426,8 +393,8 @@ def run(self, template, science, sources, visitSummary=None):
if self.config.mode == "auto":
convolveTemplate = _shapeTest(template,
science,
fwhmExposureBuffer=fwhmExposureBuffer,
fwhmExposureGrid=fwhmExposureGrid)
fwhmExposureBuffer=self.config.makeKernel.fwhmExposureBuffer,
fwhmExposureGrid=self.config.makeKernel.fwhmExposureGrid)
if convolveTemplate:
if self.sciencePsfSize < self.templatePsfSize:
self.log.info("Average template PSF size is greater, "
Expand Down Expand Up @@ -579,9 +546,13 @@ def runConvolveTemplate(self, template, science, selectSources):
try:
kernelSources = self.makeKernel.selectKernelSources(template, science,
candidateList=selectSources,
preconvolved=False)
preconvolved=False,
templateFwhmPix=self.templatePsfSize,
scienceFwhmPix=self.sciencePsfSize)
kernelResult = self.makeKernel.run(template, science, kernelSources,
preconvolved=False)
preconvolved=False,
templateFwhmPix=self.templatePsfSize,
scienceFwhmPix=self.sciencePsfSize)
except Exception as e:
if self.config.allowKernelSourceDetection:
self.log.warning("Error encountered trying to construct the matching kernel"
Expand All @@ -594,9 +565,13 @@ def runConvolveTemplate(self, template, science, selectSources):
sigma=self.sciencePsfSize/sigmaToFwhm)
kernelSources = self.makeKernel.selectKernelSources(template, science,
candidateList=candidateList,
preconvolved=False)
preconvolved=False,
templateFwhmPix=self.templatePsfSize,
scienceFwhmPix=self.sciencePsfSize)
kernelResult = self.makeKernel.run(template, science, kernelSources,
preconvolved=False)
preconvolved=False,
templateFwhmPix=self.templatePsfSize,
scienceFwhmPix=self.sciencePsfSize)
else:
raise e

Expand Down Expand Up @@ -652,9 +627,13 @@ def runConvolveScience(self, template, science, selectSources):
bbox = science.getBBox()
kernelSources = self.makeKernel.selectKernelSources(science, template,
candidateList=selectSources,
preconvolved=False)
preconvolved=False,
templateFwhmPix=self.templatePsfSize,
scienceFwhmPix=self.sciencePsfSize)
kernelResult = self.makeKernel.run(science, template, kernelSources,
preconvolved=False)
preconvolved=False,
templateFwhmPix=self.templatePsfSize,
scienceFwhmPix=self.sciencePsfSize)
modelParams = kernelResult.backgroundModel.getParameters()
# We must invert the background model if the matching kernel is solved for the science image.
kernelResult.backgroundModel.setParameters([-p for p in modelParams])
Expand Down Expand Up @@ -756,22 +735,24 @@ def _calculateMagLim(self, exposure, nsigma=5.0, fallbackPsfSize=None):
maglim : `astropy.units.Quantity`
The limiting magnitude of the exposure, or np.nan.
"""
if exposure.photoCalib is None:
return np.nan
try:
psf = exposure.getPsf()
psf_shape = psf.computeShape(psf.getAveragePosition())
except (lsst.pex.exceptions.InvalidParameterError, afwDetection.InvalidPsfError):
if fallbackPsfSize is not None:
self.log.info("Unable to evaluate PSF, using fallback FWHM %f", fallbackPsfSize)
psf_area = np.pi*(fallbackPsfSize/2)**2
zeropoint = exposure.getPhotoCalib().instFluxToMagnitude(1)
zeropoint = exposure.photoCalib.instFluxToMagnitude(1)
maglim = zeropoint - 2.5*np.log10(nsigma*np.sqrt(psf_area))
else:
self.log.info("Unable to evaluate PSF, setting maglim to nan")
maglim = np.nan
else:
# Get a more accurate area than `psf_shape.getArea()` via moments
psf_area = np.pi*np.sqrt(psf_shape.getIxx()*psf_shape.getIyy())
zeropoint = exposure.getPhotoCalib().instFluxToMagnitude(1)
zeropoint = exposure.photoCalib.instFluxToMagnitude(1)
maglim = zeropoint - 2.5*np.log10(nsigma*np.sqrt(psf_area))
finally:
return maglim
Expand Down Expand Up @@ -987,6 +968,39 @@ def _prepareInputs(self, template, science, visitSummary=None):
# We don't want the detection mask from the science image
self.updateMasks(template, science)

# Calling getPsfFwhm on template.psf fails on some rare occasions when
# the template has no input exposures at the average position of the
# stars. So we try getPsfFwhm first on template, and if that fails we
# evaluate the PSF on a grid specified by fwhmExposure* fields.
# To keep consistent definitions for PSF size on the template and
# science images, we use the same method for both.
# In the try block below, we catch two exceptions:
# 1. InvalidParameterError, in case the point where we are evaluating
# the PSF lands in a gap in the template.
# 2. RangeError, in case the template coverage is so poor that we end
# up near a region with no data.
try:
self.templatePsfSize = getPsfFwhm(template.psf)
self.sciencePsfSize = getPsfFwhm(science.psf)
except (lsst.pex.exceptions.InvalidParameterError, lsst.pex.exceptions.RangeError):
self.log.info("Unable to evaluate PSF at the average position. "
"Evaluting PSF on a grid of points."
)
self.templatePsfSize = evaluateMeanPsfFwhm(
template,
fwhmExposureBuffer=self.config.makeKernel.fwhmExposureBuffer,
fwhmExposureGrid=self.config.makeKernel.fwhmExposureGrid
)
self.sciencePsfSize = evaluateMeanPsfFwhm(
science,
fwhmExposureBuffer=self.config.makeKernel.fwhmExposureBuffer,
fwhmExposureGrid=self.config.makeKernel.fwhmExposureGrid
)
self.log.info("Science PSF FWHM: %f pixels", self.sciencePsfSize)
self.log.info("Template PSF FWHM: %f pixels", self.templatePsfSize)
self.metadata["sciencePsfSize"] = self.sciencePsfSize
self.metadata["templatePsfSize"] = self.templatePsfSize

def updateMasks(self, template, science):
"""Update the science and template mask planes before differencing.

Expand Down Expand Up @@ -1196,9 +1210,13 @@ def runPreconvolve(self, template, science, matchedScience, selectSources, preCo

kernelSources = self.makeKernel.selectKernelSources(template[innerBBox], matchedScience[innerBBox],
candidateList=selectSources,
preconvolved=True)
preconvolved=True,
templateFwhmPix=self.templatePsfSize,
scienceFwhmPix=self.sciencePsfSize)
kernelResult = self.makeKernel.run(template[innerBBox], matchedScience[innerBBox], kernelSources,
preconvolved=True)
preconvolved=True,
templateFwhmPix=self.templatePsfSize,
scienceFwhmPix=self.sciencePsfSize)

matchedTemplate = self._convolveExposure(template, kernelResult.psfMatchingKernel,
self.convolutionControl,
Expand Down