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.

Sign up for GitHub

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

Jump to bottom

DM-14738: Support DCR-matched templates in image differencing #208

Merged
merged 5 commits into from
Aug 7, 2018
Merged

DM-14738: Support DCR-matched templates in image differencing #208

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
358f55e
Add check that `coaddName` is consistent.
isullivan Jul 25, 2018
1ac3bc1
Make `psf`, `filter`, and `bbox` properties of `DcrModel`
isullivan Jul 25, 2018
554068d
Allow default `warpCtrl` for `DcrModel`s
isullivan Jul 25, 2018
c618387
Make `DcrModel`s compatible with image differencing.
isullivan Jul 25, 2018
296a4d5
Move DcrModel to ip_diffim
isullivan Jul 26, 2018
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
33 changes: 19 additions & 14 deletions python/lsst/pipe/tasks/dcrAssembleCoadd.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -25,10 +25,10 @@
import lsst.afw.image as afwImage
import lsst.afw.math as afwMath
import lsst.coadd.utils as coaddUtils
from lsst.ip.diffim.dcrModel import applyDcr, calculateDcr, DcrModel
import lsst.pex.config as pexConfig
import lsst.pipe.base as pipeBase
from .assembleCoadd import AssembleCoaddTask, CompareWarpAssembleCoaddTask, CompareWarpAssembleCoaddConfig
from .dcrModel import applyDcr, calculateDcr, DcrModel

__all__ = ["DcrAssembleCoaddTask", "DcrAssembleCoaddConfig"]

Expand Down Expand Up @@ -232,8 +232,8 @@ def prepareDcrInputs(self, templateCoadd, tempExpRefList, weightList):
NotImplementedError
If ``lambdaMin`` is missing from the Mapper class of the obs package being used.
"""
self.filterInfo = templateCoadd.getFilter()
if np.isnan(self.filterInfo.getFilterProperty().getLambdaMin()):
filterInfo = templateCoadd.getFilter()
if np.isnan(filterInfo.getFilterProperty().getLambdaMin()):
raise NotImplementedError("No minimum/maximum wavelength information found"
" in the filter definition! Please add lambdaMin and lambdaMax"
" to the Mapper class in your obs package.")
Expand All @@ -245,7 +245,7 @@ def prepareDcrInputs(self, templateCoadd, tempExpRefList, weightList):
if self.config.doAirmassWeight:
weightList[visitNum] *= airmass
dcrShifts.append(np.max(np.abs(calculateDcr(visitInfo, templateCoadd.getWcs(),
self.filterInfo, self.config.dcrNumSubfilters))))
filterInfo, self.config.dcrNumSubfilters))))
self.bufferSize = int(np.ceil(np.max(dcrShifts)) + 1)
# Turn off the warping cache, since we set the linear interpolation length to the entire subregion
# This warper is only used for applying DCR shifts, which are assumed to be uniform across a patch
Expand All @@ -256,7 +256,8 @@ def prepareDcrInputs(self, templateCoadd, tempExpRefList, weightList):
cacheSize=warpCache, interpLength=warpInterpLength)
dcrModels = DcrModel.fromImage(templateCoadd.maskedImage,
self.config.dcrNumSubfilters,
self.filterInfo)
filterInfo=filterInfo,
psf=templateCoadd.getPsf())
return dcrModels

def run(self, skyInfo, tempExpRefList, imageScalerList, weightList,
Expand Down Expand Up @@ -414,7 +415,7 @@ def calculateNImage(self, dcrModels, bbox, tempExpRefList, spanSetMaskList, stat
mask = exposure.mask
if altMaskSpans is not None:
self.applyAltMaskPlanes(mask, altMaskSpans)
dcrShift = calculateDcr(visitInfo, wcs, self.filterInfo, self.config.dcrNumSubfilters)
dcrShift = calculateDcr(visitInfo, wcs, dcrModels.filter, self.config.dcrNumSubfilters)
for dcr, dcrNImage in zip(dcrShift, dcrNImages):
shiftedImage = applyDcr(exposure.maskedImage, dcr, self.warpCtrl, useInverse=True)
dcrNImage.array[shiftedImage.mask.array & statsCtrl.getAndMask() == 0] += 1
Expand Down Expand Up @@ -464,7 +465,7 @@ def dcrAssembleSubregion(self, dcrModels, bbox, tempExpRefList, imageScalerList,
"""
bboxGrow = afwGeom.Box2I(bbox)
bboxGrow.grow(self.bufferSize)
bboxGrow.clip(dcrModels.getBBox())
bboxGrow.clip(dcrModels.bbox)

tempExpName = self.getTempExpDatasetName(self.warpType)
residualGeneratorList = []
Expand All @@ -474,24 +475,25 @@ def dcrAssembleSubregion(self, dcrModels, bbox, tempExpRefList, imageScalerList,
visitInfo = exposure.getInfo().getVisitInfo()
wcs = exposure.getInfo().getWcs()
maskedImage = exposure.maskedImage
templateImage = dcrModels.buildMatchedTemplate(self.warpCtrl, visitInfo=visitInfo, bbox=bboxGrow,
wcs=wcs, mask=baseMask)
templateImage = dcrModels.buildMatchedTemplate(warpCtrl=self.warpCtrl, visitInfo=visitInfo,
bbox=bboxGrow, wcs=wcs, mask=baseMask)
imageScaler.scaleMaskedImage(maskedImage)
if altMaskSpans is not None:
self.applyAltMaskPlanes(maskedImage.mask, altMaskSpans)

if self.config.removeMaskPlanes:
self.removeMaskPlanes(maskedImage)
maskedImage -= templateImage
residualGeneratorList.append(self.dcrResiduals(maskedImage, visitInfo, bboxGrow, wcs))
residualGeneratorList.append(self.dcrResiduals(maskedImage, visitInfo, bboxGrow, wcs,
dcrModels.filter))

dcrSubModelOut = self.newModelFromResidual(dcrModels, residualGeneratorList, bboxGrow,
statsFlags, statsCtrl, weightList)
dcrSubModelOut.regularizeModel(bboxGrow, baseMask, statsCtrl, self.config.regularizeSigma,
self.config.clampFrequency, self.config.convergenceMaskPlanes)
dcrModels.assign(dcrSubModelOut, bbox)

def dcrResiduals(self, residual, visitInfo, bbox, wcs):
def dcrResiduals(self, residual, visitInfo, bbox, wcs, filterInfo):
"""Prepare a residual image for stacking in each subfilter by applying the reverse DCR shifts.

Parameters
Expand All @@ -505,13 +507,16 @@ def dcrResiduals(self, residual, visitInfo, bbox, wcs):
Sub-region of the coadd
wcs : `lsst.afw.geom.SkyWcs`
Coordinate system definition (wcs) for the exposure.
filterInfo : `lsst.afw.image.Filter`
The filter definition, set in the current instruments' obs package.
Required for any calculation of DCR, including making matched templates.

Yields
------
residualImage : `lsst.afw.image.maskedImageF`
The residual image for the next subfilter, shifted for DCR.
"""
dcrShift = calculateDcr(visitInfo, wcs, self.filterInfo, self.config.dcrNumSubfilters)
dcrShift = calculateDcr(visitInfo, wcs, filterInfo, self.config.dcrNumSubfilters)
for dcr in dcrShift:
yield applyDcr(residual, dcr, self.warpCtrl, bbox=bbox, useInverse=True)

Expand Down Expand Up @@ -553,7 +558,7 @@ def newModelFromResidual(self, dcrModels, residualGeneratorList, bbox, statsFlag
self.config.modelClampFactor, self.config.convergenceMaskPlanes)
dcrModels.conditionDcrModel(subfilter, newModel, bbox, gain=1.)
newModelImages.append(newModel)
return DcrModel(newModelImages, self.filterInfo)
return DcrModel(newModelImages, dcrModels.filter, dcrModels.psf)

def calculateConvergence(self, dcrModels, bbox, tempExpRefList, imageScalerList,
weightList, spanSetMaskList, statsCtrl):
Expand Down Expand Up @@ -622,7 +627,7 @@ def calculateSingleConvergence(self, dcrModels, exposure, significanceImage,
Quality of fit metric for one exposure, within the sub-region.
"""
convergeMask = exposure.mask.getPlaneBitMask(self.config.convergenceMaskPlanes)
templateImage = dcrModels.buildMatchedTemplate(self.warpCtrl,
templateImage = dcrModels.buildMatchedTemplate(warpCtrl=self.warpCtrl,
visitInfo=exposure.getInfo().getVisitInfo(),
bbox=exposure.getBBox(),
wcs=exposure.getInfo().getWcs())
Expand Down