LSST Applications g0b6bd0c080+a72a5dd7e6,g1182afd7b4+2a019aa3bb,g17e5ecfddb+2b8207f7de,g1d67935e3f+06cf436103,g38293774b4+ac198e9f13,g396055baef+6a2097e274,g3b44f30a73+6611e0205b,g480783c3b1+98f8679e14,g48ccf36440+89c08d0516,g4b93dc025c+98f8679e14,g5c4744a4d9+a302e8c7f0,g613e996a0d+e1c447f2e0,g6c8d09e9e7+25247a063c,g7271f0639c+98f8679e14,g7a9cd813b8+124095ede6,g9d27549199+a302e8c7f0,ga1cf026fa3+ac198e9f13,ga32aa97882+7403ac30ac,ga786bb30fb+7a139211af,gaa63f70f4e+9994eb9896,gabf319e997+ade567573c,gba47b54d5d+94dc90c3ea,gbec6a3398f+06cf436103,gc6308e37c7+07dd123edb,gc655b1545f+ade567573c,gcc9029db3c+ab229f5caf,gd01420fc67+06cf436103,gd877ba84e5+06cf436103,gdb4cecd868+6f279b5b48,ge2d134c3d5+cc4dbb2e3f,ge448b5faa6+86d1ceac1d,gecc7e12556+98f8679e14,gf3ee170dca+25247a063c,gf4ac96e456+ade567573c,gf9f5ea5b4d+ac198e9f13,gff490e6085+8c2580be5c,w.2022.27
LSST Data Management Base Package
Classes | Functions | Variables
lsst.pipe.tasks.dcrAssembleCoadd Namespace Reference

Classes

class  DcrAssembleCoaddConnections
 

Functions

def dcrAssembleSubregion (self, dcrModels, subExposures, bbox, dcrBBox, warpRefList, statsCtrl, convergenceMetric, gain, modelWeights, refImage, dcrWeights)
 
def dcrResiduals (self, residual, visitInfo, wcs, effectiveWavelength, bandwidth)
 
def newModelFromResidual (self, dcrModels, residualGeneratorList, dcrBBox, statsCtrl, gain, modelWeights, refImage, dcrWeights)
 
def calculateConvergence (self, dcrModels, subExposures, bbox, warpRefList, weightList, statsCtrl)
 
def calculateSingleConvergence (self, dcrModels, exposure, significanceImage, statsCtrl)
 
def stackCoadd (self, dcrCoadds)
 
def fillCoadd (self, dcrModels, skyInfo, warpRefList, weightList, calibration=None, coaddInputs=None, mask=None, variance=None)
 
def calculateGain (self, convergenceList, gainList)
 
def calculateModelWeights (self, dcrModels, dcrBBox)
 
def applyModelWeights (self, modelImages, refImage, modelWeights)
 
def loadSubExposures (self, bbox, statsCtrl, warpRefList, imageScalerList, spanSetMaskList)
 
def selectCoaddPsf (self, templateCoadd, warpRefList)
 

Variables

int weightsThreshold = 1.
 
int goodPix = dcrWeights[0].array > weightsThreshold
 

Function Documentation

◆ applyModelWeights()

def lsst.pipe.tasks.dcrAssembleCoadd.applyModelWeights (   self,
  modelImages,
  refImage,
  modelWeights 
)
Smoothly replace model pixel values with those from a
reference at locations away from detected sources.

Parameters
----------
modelImages : `list` of `lsst.afw.image.Image`
    The new DCR model images from the current iteration.
    The values will be modified in place.
refImage : `lsst.afw.image.MaskedImage`
    A reference image used to supply the default pixel values.
modelWeights : `numpy.ndarray` or `float`
    A 2D array of weight values that tapers smoothly to zero away from detected sources.
    Set to a placeholder value of 1.0 if ``self.config.useModelWeights`` is False.

Definition at line 1247 of file dcrAssembleCoadd.py.

1247 def applyModelWeights(self, modelImages, refImage, modelWeights):
1248 """Smoothly replace model pixel values with those from a
1249 reference at locations away from detected sources.
1250
1251 Parameters
1252 ----------
1253 modelImages : `list` of `lsst.afw.image.Image`
1254 The new DCR model images from the current iteration.
1255 The values will be modified in place.
1256 refImage : `lsst.afw.image.MaskedImage`
1257 A reference image used to supply the default pixel values.
1258 modelWeights : `numpy.ndarray` or `float`
1259 A 2D array of weight values that tapers smoothly to zero away from detected sources.
1260 Set to a placeholder value of 1.0 if ``self.config.useModelWeights`` is False.
1261 """
1262 if self.config.useModelWeights:
1263 for model in modelImages:
1264 model.array *= modelWeights
1265 model.array += refImage.array*(1. - modelWeights)/self.config.dcrNumSubfilters
1266
A class to represent a 2-dimensional array of pixels.
Definition: Image.h:51
A class to manipulate images, masks, and variance as a single object.
Definition: MaskedImage.h:73
def applyModelWeights(self, modelImages, refImage, modelWeights)

◆ calculateConvergence()

def lsst.pipe.tasks.dcrAssembleCoadd.calculateConvergence (   self,
  dcrModels,
  subExposures,
  bbox,
  warpRefList,
  weightList,
  statsCtrl 
)
Calculate a quality of fit metric for the matched templates.

Parameters
----------
dcrModels : `lsst.pipe.tasks.DcrModel`
    Best fit model of the true sky after correcting chromatic effects.
subExposures : `dict` of `lsst.afw.image.ExposureF`
    The pre-loaded exposures for the current subregion.
bbox : `lsst.geom.box.Box2I`
    Sub-region to coadd
warpRefList : `list` of `lsst.daf.butler.DeferredDatasetHandle` or
    `lsst.daf.persistence.ButlerDataRef`
    The data references to the input warped exposures.
weightList : `list` of `float`
    The weight to give each input exposure in the coadd
statsCtrl : `lsst.afw.math.StatisticsControl`
    Statistics control object for coadd

Returns
-------
convergenceMetric : `float`
    Quality of fit metric for all input exposures, within the sub-region

Definition at line 971 of file dcrAssembleCoadd.py.

971 def calculateConvergence(self, dcrModels, subExposures, bbox, warpRefList, weightList, statsCtrl):
972 """Calculate a quality of fit metric for the matched templates.
973
974 Parameters
975 ----------
976 dcrModels : `lsst.pipe.tasks.DcrModel`
977 Best fit model of the true sky after correcting chromatic effects.
978 subExposures : `dict` of `lsst.afw.image.ExposureF`
979 The pre-loaded exposures for the current subregion.
980 bbox : `lsst.geom.box.Box2I`
981 Sub-region to coadd
982 warpRefList : `list` of `lsst.daf.butler.DeferredDatasetHandle` or
984 The data references to the input warped exposures.
985 weightList : `list` of `float`
986 The weight to give each input exposure in the coadd
988 Statistics control object for coadd
989
990 Returns
991 -------
992 convergenceMetric : `float`
993 Quality of fit metric for all input exposures, within the sub-region
994 """
995 significanceImage = np.abs(dcrModels.getReferenceImage(bbox))
996 nSigma = 3.
997 significanceImage += nSigma*dcrModels.calculateNoiseCutoff(dcrModels[1], statsCtrl,
998 bufferSize=self.bufferSize)
999 if np.max(significanceImage) == 0:
1000 significanceImage += 1.
1001 weight = 0
1002 metric = 0.
1003 metricList = {}
1004 for warpExpRef, expWeight in zip(warpRefList, weightList):
1005 visit = warpExpRef.dataId["visit"]
1006 exposure = subExposures[visit][bbox]
1007 singleMetric = self.calculateSingleConvergence(dcrModels, exposure, significanceImage, statsCtrl)
1008 metric += singleMetric
1009 metricList[visit] = singleMetric
1010 weight += 1.
1011 self.log.info("Individual metrics:\n%s", metricList)
1012 return 1.0 if weight == 0.0 else metric/weight
1013
Pass parameters to a Statistics object.
Definition: Statistics.h:92
def calculateConvergence(self, dcrModels, subExposures, bbox, warpRefList, weightList, statsCtrl)

◆ calculateGain()

def lsst.pipe.tasks.dcrAssembleCoadd.calculateGain (   self,
  convergenceList,
  gainList 
)
Calculate the gain to use for the current iteration.

After calculating a new DcrModel, each value is averaged with the
value in the corresponding pixel from the previous iteration. This
reduces oscillating solutions that iterative techniques are plagued by,
and speeds convergence. By far the biggest changes to the model
happen in the first couple iterations, so we can also use a more
aggressive gain later when the model is changing slowly.

Parameters
----------
convergenceList : `list` of `float`
    The quality of fit metric from each previous iteration.
gainList : `list` of `float`
    The gains used in each previous iteration: appended with the new
    gain value.
    Gains are numbers between ``self.config.baseGain`` and 1.

Returns
-------
gain : `float`
    Relative weight to give the new solution when updating the model.
    A value of 1.0 gives equal weight to both solutions.

Raises
------
ValueError
    If ``len(convergenceList) != len(gainList)+1``.

Definition at line 1133 of file dcrAssembleCoadd.py.

1133 def calculateGain(self, convergenceList, gainList):
1134 """Calculate the gain to use for the current iteration.
1135
1136 After calculating a new DcrModel, each value is averaged with the
1137 value in the corresponding pixel from the previous iteration. This
1138 reduces oscillating solutions that iterative techniques are plagued by,
1139 and speeds convergence. By far the biggest changes to the model
1140 happen in the first couple iterations, so we can also use a more
1141 aggressive gain later when the model is changing slowly.
1142
1143 Parameters
1144 ----------
1145 convergenceList : `list` of `float`
1146 The quality of fit metric from each previous iteration.
1147 gainList : `list` of `float`
1148 The gains used in each previous iteration: appended with the new
1149 gain value.
1150 Gains are numbers between ``self.config.baseGain`` and 1.
1151
1152 Returns
1153 -------
1154 gain : `float`
1155 Relative weight to give the new solution when updating the model.
1156 A value of 1.0 gives equal weight to both solutions.
1157
1158 Raises
1159 ------
1160 ValueError
1161 If ``len(convergenceList) != len(gainList)+1``.
1162 """
1163 nIter = len(convergenceList)
1164 if nIter != len(gainList) + 1:
1165 raise ValueError("convergenceList (%d) must be one element longer than gainList (%d)."
1166 % (len(convergenceList), len(gainList)))
1167
1168 if self.config.baseGain is None:
1169 # If ``baseGain`` is not set, calculate it from the number of DCR subfilters
1170 # The more subfilters being modeled, the lower the gain should be.
1171 baseGain = 1./(self.config.dcrNumSubfilters - 1)
1172 else:
1173 baseGain = self.config.baseGain
1174
1175 if self.config.useProgressiveGain and nIter > 2:
1176 # To calculate the best gain to use, compare the past gains that have been used
1177 # with the resulting convergences to estimate the best gain to use.
1178 # Algorithmically, this is a Kalman filter.
1179 # If forward modeling proceeds perfectly, the convergence metric should
1180 # asymptotically approach a final value.
1181 # We can estimate that value from the measured changes in convergence
1182 # weighted by the gains used in each previous iteration.
1183 estFinalConv = [((1 + gainList[i])*convergenceList[i + 1] - convergenceList[i])/gainList[i]
1184 for i in range(nIter - 1)]
1185 # The convergence metric is strictly positive, so if the estimated final convergence is
1186 # less than zero, force it to zero.
1187 estFinalConv = np.array(estFinalConv)
1188 estFinalConv[estFinalConv < 0] = 0
1189 # Because the estimate may slowly change over time, only use the most recent measurements.
1190 estFinalConv = np.median(estFinalConv[max(nIter - 5, 0):])
1191 lastGain = gainList[-1]
1192 lastConv = convergenceList[-2]
1193 newConv = convergenceList[-1]
1194 # The predicted convergence is the value we would get if the new model calculated
1195 # in the previous iteration was perfect. Recall that the updated model that is
1196 # actually used is the gain-weighted average of the new and old model,
1197 # so the convergence would be similarly weighted.
1198 predictedConv = (estFinalConv*lastGain + lastConv)/(1. + lastGain)
1199 # If the measured and predicted convergence are very close, that indicates
1200 # that our forward model is accurate and we can use a more aggressive gain
1201 # If the measured convergence is significantly worse (or better!) than predicted,
1202 # that indicates that the model is not converging as expected and
1203 # we should use a more conservative gain.
1204 delta = (predictedConv - newConv)/((lastConv - estFinalConv)/(1 + lastGain))
1205 newGain = 1 - abs(delta)
1206 # Average the gains to prevent oscillating solutions.
1207 newGain = (newGain + lastGain)/2.
1208 gain = max(baseGain, newGain)
1209 else:
1210 gain = baseGain
1211 gainList.append(gain)
1212 return gain
1213
int max
def calculateGain(self, convergenceList, gainList)
Angle abs(Angle const &a)
Definition: Angle.h:106

◆ calculateModelWeights()

def lsst.pipe.tasks.dcrAssembleCoadd.calculateModelWeights (   self,
  dcrModels,
  dcrBBox 
)
Build an array that smoothly tapers to 0 away from detected sources.

Parameters
----------
dcrModels : `lsst.pipe.tasks.DcrModel`
    Best fit model of the true sky after correcting chromatic effects.
dcrBBox : `lsst.geom.box.Box2I`
    Sub-region of the coadd which includes a buffer to allow for DCR.

Returns
-------
weights : `numpy.ndarray` or `float`
    A 2D array of weight values that tapers smoothly to zero away from detected sources.
    Set to a placeholder value of 1.0 if ``self.config.useModelWeights`` is False.

Raises
------
ValueError
    If ``useModelWeights`` is set and ``modelWeightsWidth`` is negative.

Definition at line 1214 of file dcrAssembleCoadd.py.

1214 def calculateModelWeights(self, dcrModels, dcrBBox):
1215 """Build an array that smoothly tapers to 0 away from detected sources.
1216
1217 Parameters
1218 ----------
1219 dcrModels : `lsst.pipe.tasks.DcrModel`
1220 Best fit model of the true sky after correcting chromatic effects.
1221 dcrBBox : `lsst.geom.box.Box2I`
1222 Sub-region of the coadd which includes a buffer to allow for DCR.
1223
1224 Returns
1225 -------
1226 weights : `numpy.ndarray` or `float`
1227 A 2D array of weight values that tapers smoothly to zero away from detected sources.
1228 Set to a placeholder value of 1.0 if ``self.config.useModelWeights`` is False.
1229
1230 Raises
1231 ------
1232 ValueError
1233 If ``useModelWeights`` is set and ``modelWeightsWidth`` is negative.
1234 """
1235 if not self.config.useModelWeights:
1236 return 1.0
1237 if self.config.modelWeightsWidth < 0:
1238 raise ValueError("modelWeightsWidth must not be negative if useModelWeights is set")
1239 convergeMask = dcrModels.mask.getPlaneBitMask(self.config.convergenceMaskPlanes)
1240 convergeMaskPixels = dcrModels.mask[dcrBBox].array & convergeMask > 0
1241 weights = np.zeros_like(dcrModels[0][dcrBBox].array)
1242 weights[convergeMaskPixels] = 1.
1243 weights = ndimage.gaussian_filter(weights, self.config.modelWeightsWidth)
1244 weights /= np.max(weights)
1245 return weights
1246
def calculateModelWeights(self, dcrModels, dcrBBox)

◆ calculateSingleConvergence()

def lsst.pipe.tasks.dcrAssembleCoadd.calculateSingleConvergence (   self,
  dcrModels,
  exposure,
  significanceImage,
  statsCtrl 
)
Calculate a quality of fit metric for a single matched template.

Parameters
----------
dcrModels : `lsst.pipe.tasks.DcrModel`
    Best fit model of the true sky after correcting chromatic effects.
exposure : `lsst.afw.image.ExposureF`
    The input warped exposure to evaluate.
significanceImage : `numpy.ndarray`
    Array of weights for each pixel corresponding to its significance
    for the convergence calculation.
statsCtrl : `lsst.afw.math.StatisticsControl`
    Statistics control object for coadd

Returns
-------
convergenceMetric : `float`
    Quality of fit metric for one exposure, within the sub-region.

Definition at line 1014 of file dcrAssembleCoadd.py.

1014 def calculateSingleConvergence(self, dcrModels, exposure, significanceImage, statsCtrl):
1015 """Calculate a quality of fit metric for a single matched template.
1016
1017 Parameters
1018 ----------
1019 dcrModels : `lsst.pipe.tasks.DcrModel`
1020 Best fit model of the true sky after correcting chromatic effects.
1021 exposure : `lsst.afw.image.ExposureF`
1022 The input warped exposure to evaluate.
1023 significanceImage : `numpy.ndarray`
1024 Array of weights for each pixel corresponding to its significance
1025 for the convergence calculation.
1027 Statistics control object for coadd
1028
1029 Returns
1030 -------
1031 convergenceMetric : `float`
1032 Quality of fit metric for one exposure, within the sub-region.
1033 """
1034 convergeMask = exposure.mask.getPlaneBitMask(self.config.convergenceMaskPlanes)
1035 templateImage = dcrModels.buildMatchedTemplate(exposure=exposure,
1036 bbox=exposure.getBBox(),
1037 order=self.config.imageInterpOrder,
1038 splitSubfilters=self.config.splitSubfilters,
1039 splitThreshold=self.config.splitThreshold,
1040 amplifyModel=self.config.accelerateModel)
1041 diffVals = np.abs(exposure.image.array - templateImage.array)*significanceImage
1042 refVals = np.abs(exposure.image.array + templateImage.array)*significanceImage/2.
1043
1044 finitePixels = np.isfinite(diffVals)
1045 goodMaskPixels = (exposure.mask.array & statsCtrl.getAndMask()) == 0
1046 convergeMaskPixels = exposure.mask.array & convergeMask > 0
1047 usePixels = finitePixels & goodMaskPixels & convergeMaskPixels
1048 if np.sum(usePixels) == 0:
1049 metric = 0.
1050 else:
1051 diffUse = diffVals[usePixels]
1052 refUse = refVals[usePixels]
1053 metric = np.sum(diffUse/np.median(diffUse))/np.sum(refUse/np.median(diffUse))
1054 return metric
1055
def calculateSingleConvergence(self, dcrModels, exposure, significanceImage, statsCtrl)

◆ dcrAssembleSubregion()

def lsst.pipe.tasks.dcrAssembleCoadd.dcrAssembleSubregion (   self,
  dcrModels,
  subExposures,
  bbox,
  dcrBBox,
  warpRefList,
  statsCtrl,
  convergenceMetric,
  gain,
  modelWeights,
  refImage,
  dcrWeights 
)
Assemble the DCR coadd for a sub-region.

Build a DCR-matched template for each input exposure, then shift the
residuals according to the DCR in each subfilter.
Stack the shifted residuals and apply them as a correction to the
solution from the previous iteration.
Restrict the new model solutions from varying by more than a factor of
`modelClampFactor` from the last solution, and additionally restrict the
individual subfilter models from varying by more than a factor of
`frequencyClampFactor` from their average.
Finally, mitigate potentially oscillating solutions by averaging the new
solution with the solution from the previous iteration, weighted by
their convergence metric.

Parameters
----------
dcrModels : `lsst.pipe.tasks.DcrModel`
    Best fit model of the true sky after correcting chromatic effects.
subExposures : `dict` of `lsst.afw.image.ExposureF`
    The pre-loaded exposures for the current subregion.
bbox : `lsst.geom.box.Box2I`
    Bounding box of the subregion to coadd.
dcrBBox : `lsst.geom.box.Box2I`
    Sub-region of the coadd which includes a buffer to allow for DCR.
warpRefList : `list` of `lsst.daf.butler.DeferredDatasetHandle` or
    `lsst.daf.persistence.ButlerDataRef`
    The data references to the input warped exposures.
statsCtrl : `lsst.afw.math.StatisticsControl`
    Statistics control object for coadd
convergenceMetric : `float`
    Quality of fit metric for the matched templates of the input images.
gain : `float`, optional
    Relative weight to give the new solution when updating the model.
modelWeights : `numpy.ndarray` or `float`
    A 2D array of weight values that tapers smoothly to zero away from detected sources.
    Set to a placeholder value of 1.0 if ``self.config.useModelWeights`` is False.
refImage : `lsst.afw.image.Image`
    A reference image used to supply the default pixel values.
dcrWeights : `list` of `lsst.afw.image.Image`
    Per-pixel weights for each subfilter.
    Equal to 1/(number of unmasked images contributing to each pixel).

Definition at line 807 of file dcrAssembleCoadd.py.

809 gain, modelWeights, refImage, dcrWeights):
810 """Assemble the DCR coadd for a sub-region.
811
812 Build a DCR-matched template for each input exposure, then shift the
813 residuals according to the DCR in each subfilter.
814 Stack the shifted residuals and apply them as a correction to the
815 solution from the previous iteration.
816 Restrict the new model solutions from varying by more than a factor of
817 `modelClampFactor` from the last solution, and additionally restrict the
818 individual subfilter models from varying by more than a factor of
819 `frequencyClampFactor` from their average.
820 Finally, mitigate potentially oscillating solutions by averaging the new
821 solution with the solution from the previous iteration, weighted by
822 their convergence metric.
823
824 Parameters
825 ----------
826 dcrModels : `lsst.pipe.tasks.DcrModel`
827 Best fit model of the true sky after correcting chromatic effects.
828 subExposures : `dict` of `lsst.afw.image.ExposureF`
829 The pre-loaded exposures for the current subregion.
830 bbox : `lsst.geom.box.Box2I`
831 Bounding box of the subregion to coadd.
832 dcrBBox : `lsst.geom.box.Box2I`
833 Sub-region of the coadd which includes a buffer to allow for DCR.
834 warpRefList : `list` of `lsst.daf.butler.DeferredDatasetHandle` or
836 The data references to the input warped exposures.
838 Statistics control object for coadd
839 convergenceMetric : `float`
840 Quality of fit metric for the matched templates of the input images.
841 gain : `float`, optional
842 Relative weight to give the new solution when updating the model.
843 modelWeights : `numpy.ndarray` or `float`
844 A 2D array of weight values that tapers smoothly to zero away from detected sources.
845 Set to a placeholder value of 1.0 if ``self.config.useModelWeights`` is False.
846 refImage : `lsst.afw.image.Image`
847 A reference image used to supply the default pixel values.
848 dcrWeights : `list` of `lsst.afw.image.Image`
849 Per-pixel weights for each subfilter.
850 Equal to 1/(number of unmasked images contributing to each pixel).
851 """
852 residualGeneratorList = []
853
854 for warpExpRef in warpRefList:
855 visit = warpExpRef.dataId["visit"]
856 exposure = subExposures[visit]
857 visitInfo = exposure.getInfo().getVisitInfo()
858 wcs = exposure.getInfo().getWcs()
859 templateImage = dcrModels.buildMatchedTemplate(exposure=exposure,
860 bbox=exposure.getBBox(),
861 order=self.config.imageInterpOrder,
862 splitSubfilters=self.config.splitSubfilters,
863 splitThreshold=self.config.splitThreshold,
864 amplifyModel=self.config.accelerateModel)
865 residual = exposure.image.array - templateImage.array
866 # Note that the variance plane here is used to store weights, not the actual variance
867 residual *= exposure.variance.array
868 # The residuals are stored as a list of generators.
869 # This allows the residual for a given subfilter and exposure to be created
870 # on the fly, instead of needing to store them all in memory.
871 residualGeneratorList.append(self.dcrResiduals(residual, visitInfo, wcs,
872 dcrModels.effectiveWavelength,
873 dcrModels.bandwidth))
874
875 dcrSubModelOut = self.newModelFromResidual(dcrModels, residualGeneratorList, dcrBBox, statsCtrl,
876 gain=gain,
877 modelWeights=modelWeights,
878 refImage=refImage,
879 dcrWeights=dcrWeights)
880 dcrModels.assign(dcrSubModelOut, bbox)
881

◆ dcrResiduals()

def lsst.pipe.tasks.dcrAssembleCoadd.dcrResiduals (   self,
  residual,
  visitInfo,
  wcs,
  effectiveWavelength,
  bandwidth 
)
Prepare a residual image for stacking in each subfilter by applying the reverse DCR shifts.

Parameters
----------
residual : `numpy.ndarray`
    The residual masked image for one exposure,
    after subtracting the matched template
visitInfo : `lsst.afw.image.VisitInfo`
    Metadata for the exposure.
wcs : `lsst.afw.geom.SkyWcs`
    Coordinate system definition (wcs) for the exposure.

Yields
------
residualImage : `numpy.ndarray`
    The residual image for the next subfilter, shifted for DCR.

Definition at line 882 of file dcrAssembleCoadd.py.

882 def dcrResiduals(self, residual, visitInfo, wcs, effectiveWavelength, bandwidth):
883 """Prepare a residual image for stacking in each subfilter by applying the reverse DCR shifts.
884
885 Parameters
886 ----------
887 residual : `numpy.ndarray`
888 The residual masked image for one exposure,
889 after subtracting the matched template
890 visitInfo : `lsst.afw.image.VisitInfo`
891 Metadata for the exposure.
893 Coordinate system definition (wcs) for the exposure.
894
895 Yields
896 ------
897 residualImage : `numpy.ndarray`
898 The residual image for the next subfilter, shifted for DCR.
899 """
900 if self.config.imageInterpOrder > 1:
901 # Pre-calculate the spline-filtered residual image, so that step can be
902 # skipped in the shift calculation in `applyDcr`.
903 filteredResidual = ndimage.spline_filter(residual, order=self.config.imageInterpOrder)
904 else:
905 # No need to prefilter if order=1 (it will also raise an error)
906 filteredResidual = residual
907 # Note that `splitSubfilters` is always turned off in the reverse direction.
908 # This option introduces additional blurring if applied to the residuals.
909 dcrShift = calculateDcr(visitInfo, wcs, effectiveWavelength, bandwidth, self.config.dcrNumSubfilters,
910 splitSubfilters=False)
911 for dcr in dcrShift:
912 yield applyDcr(filteredResidual, dcr, useInverse=True, splitSubfilters=False,
913 doPrefilter=False, order=self.config.imageInterpOrder)
914
A 2-dimensional celestial WCS that transform pixels to ICRS RA/Dec, using the LSST standard for pixel...
Definition: SkyWcs.h:117
Information about a single exposure of an imaging camera.
Definition: VisitInfo.h:68
def applyDcr(image, dcr, useInverse=False, splitSubfilters=False, splitThreshold=0., doPrefilter=True, order=3)
Definition: dcrModel.py:667
def calculateDcr(visitInfo, wcs, effectiveWavelength, bandwidth, dcrNumSubfilters, splitSubfilters=False)
Definition: dcrModel.py:732
def dcrResiduals(self, residual, visitInfo, wcs, effectiveWavelength, bandwidth)

◆ fillCoadd()

def lsst.pipe.tasks.dcrAssembleCoadd.fillCoadd (   self,
  dcrModels,
  skyInfo,
  warpRefList,
  weightList,
  calibration = None,
  coaddInputs = None,
  mask = None,
  variance = None 
)
Create a list of coadd exposures from a list of masked images.

Parameters
----------
dcrModels : `lsst.pipe.tasks.DcrModel`
    Best fit model of the true sky after correcting chromatic effects.
skyInfo : `lsst.pipe.base.Struct`
    Patch geometry information, from getSkyInfo
warpRefList : `list` of `lsst.daf.butler.DeferredDatasetHandle` or
    `lsst.daf.persistence.ButlerDataRef`
    The data references to the input warped exposures.
weightList : `list` of `float`
    The weight to give each input exposure in the coadd
calibration : `lsst.afw.Image.PhotoCalib`, optional
    Scale factor to set the photometric calibration of an exposure.
coaddInputs : `lsst.afw.Image.CoaddInputs`, optional
    A record of the observations that are included in the coadd.
mask : `lsst.afw.image.Mask`, optional
    Optional mask to override the values in the final coadd.
variance : `lsst.afw.image.Image`, optional
    Optional variance plane to override the values in the final coadd.

Returns
-------
dcrCoadds : `list` of `lsst.afw.image.ExposureF`
    A list of coadd exposures, each exposure containing
    the model for one subfilter.

Definition at line 1075 of file dcrAssembleCoadd.py.

1076 mask=None, variance=None):
1077 """Create a list of coadd exposures from a list of masked images.
1078
1079 Parameters
1080 ----------
1081 dcrModels : `lsst.pipe.tasks.DcrModel`
1082 Best fit model of the true sky after correcting chromatic effects.
1083 skyInfo : `lsst.pipe.base.Struct`
1084 Patch geometry information, from getSkyInfo
1085 warpRefList : `list` of `lsst.daf.butler.DeferredDatasetHandle` or
1087 The data references to the input warped exposures.
1088 weightList : `list` of `float`
1089 The weight to give each input exposure in the coadd
1090 calibration : `lsst.afw.Image.PhotoCalib`, optional
1091 Scale factor to set the photometric calibration of an exposure.
1092 coaddInputs : `lsst.afw.Image.CoaddInputs`, optional
1093 A record of the observations that are included in the coadd.
1094 mask : `lsst.afw.image.Mask`, optional
1095 Optional mask to override the values in the final coadd.
1096 variance : `lsst.afw.image.Image`, optional
1097 Optional variance plane to override the values in the final coadd.
1098
1099 Returns
1100 -------
1101 dcrCoadds : `list` of `lsst.afw.image.ExposureF`
1102 A list of coadd exposures, each exposure containing
1103 the model for one subfilter.
1104 """
1105 dcrCoadds = []
1106 refModel = dcrModels.getReferenceImage()
1107 for model in dcrModels:
1108 if self.config.accelerateModel > 1:
1109 model.array = (model.array - refModel)*self.config.accelerateModel + refModel
1110 coaddExposure = afwImage.ExposureF(skyInfo.bbox, skyInfo.wcs)
1111 if calibration is not None:
1112 coaddExposure.setPhotoCalib(calibration)
1113 if coaddInputs is not None:
1114 coaddExposure.getInfo().setCoaddInputs(coaddInputs)
1115 # Set the metadata for the coadd, including PSF and aperture corrections.
1116 self.assembleMetadata(coaddExposure, warpRefList, weightList)
1117 # Overwrite the PSF
1118 coaddExposure.setPsf(dcrModels.psf)
1119 coaddUtils.setCoaddEdgeBits(dcrModels.mask[skyInfo.bbox], dcrModels.variance[skyInfo.bbox])
1120 maskedImage = afwImage.MaskedImageF(dcrModels.bbox)
1121 maskedImage.image = model
1122 maskedImage.mask = dcrModels.mask
1123 maskedImage.variance = dcrModels.variance
1124 coaddExposure.setMaskedImage(maskedImage[skyInfo.bbox])
1125 coaddExposure.setPhotoCalib(self.scaleZeroPoint.getPhotoCalib())
1126 if mask is not None:
1127 coaddExposure.setMask(mask)
1128 if variance is not None:
1129 coaddExposure.setVariance(variance)
1130 dcrCoadds.append(coaddExposure)
1131 return dcrCoadds
1132
Represent a 2-dimensional array of bitmask pixels.
Definition: Mask.h:77
void setCoaddEdgeBits(lsst::afw::image::Mask< lsst::afw::image::MaskPixel > &coaddMask, lsst::afw::image::Image< WeightPixelT > const &weightMap)
set edge bits of coadd mask based on weight map

◆ loadSubExposures()

def lsst.pipe.tasks.dcrAssembleCoadd.loadSubExposures (   self,
  bbox,
  statsCtrl,
  warpRefList,
  imageScalerList,
  spanSetMaskList 
)
Pre-load sub-regions of a list of exposures.

Parameters
----------
bbox : `lsst.geom.box.Box2I`
    Sub-region to coadd
statsCtrl : `lsst.afw.math.StatisticsControl`
    Statistics control object for coadd
warpRefList : `list` of `lsst.daf.butler.DeferredDatasetHandle` or
    `lsst.daf.persistence.ButlerDataRef`
    The data references to the input warped exposures.
imageScalerList : `list` of `lsst.pipe.task.ImageScaler`
    The image scalars correct for the zero point of the exposures.
spanSetMaskList : `list` of `dict` containing spanSet lists, or None
    Each element is dict with keys = mask plane name to add the spans to

Returns
-------
subExposures : `dict`
    The `dict` keys are the visit IDs,
    and the values are `lsst.afw.image.ExposureF`
    The pre-loaded exposures for the current subregion.
    The variance plane contains weights, and not the variance

Definition at line 1267 of file dcrAssembleCoadd.py.

1267 def loadSubExposures(self, bbox, statsCtrl, warpRefList, imageScalerList, spanSetMaskList):
1268 """Pre-load sub-regions of a list of exposures.
1269
1270 Parameters
1271 ----------
1272 bbox : `lsst.geom.box.Box2I`
1273 Sub-region to coadd
1275 Statistics control object for coadd
1276 warpRefList : `list` of `lsst.daf.butler.DeferredDatasetHandle` or
1278 The data references to the input warped exposures.
1279 imageScalerList : `list` of `lsst.pipe.task.ImageScaler`
1280 The image scalars correct for the zero point of the exposures.
1281 spanSetMaskList : `list` of `dict` containing spanSet lists, or None
1282 Each element is dict with keys = mask plane name to add the spans to
1283
1284 Returns
1285 -------
1286 subExposures : `dict`
1287 The `dict` keys are the visit IDs,
1288 and the values are `lsst.afw.image.ExposureF`
1289 The pre-loaded exposures for the current subregion.
1290 The variance plane contains weights, and not the variance
1291 """
1292 tempExpName = self.getTempExpDatasetName(self.warpType)
1293 zipIterables = zip(warpRefList, imageScalerList, spanSetMaskList)
1294 subExposures = {}
1295 for warpExpRef, imageScaler, altMaskSpans in zipIterables:
1296 if isinstance(warpExpRef, DeferredDatasetHandle):
1297 exposure = warpExpRef.get(parameters={'bbox': bbox})
1298 else:
1299 exposure = warpExpRef.get(tempExpName + "_sub", bbox=bbox)
1300 visit = warpExpRef.dataId["visit"]
1301 if altMaskSpans is not None:
1302 self.applyAltMaskPlanes(exposure.mask, altMaskSpans)
1303 imageScaler.scaleMaskedImage(exposure.maskedImage)
1304 # Note that the variance plane here is used to store weights, not the actual variance
1305 exposure.variance.array[:, :] = 0.
1306 # Set the weight of unmasked pixels to 1.
1307 exposure.variance.array[(exposure.mask.array & statsCtrl.getAndMask()) == 0] = 1.
1308 # Set the image value of masked pixels to zero.
1309 # This eliminates needing the mask plane when stacking images in ``newModelFromResidual``
1310 exposure.image.array[(exposure.mask.array & statsCtrl.getAndMask()) > 0] = 0.
1311 subExposures[visit] = exposure
1312 return subExposures
1313
def loadSubExposures(self, bbox, statsCtrl, warpRefList, imageScalerList, spanSetMaskList)

◆ newModelFromResidual()

def lsst.pipe.tasks.dcrAssembleCoadd.newModelFromResidual (   self,
  dcrModels,
  residualGeneratorList,
  dcrBBox,
  statsCtrl,
  gain,
  modelWeights,
  refImage,
  dcrWeights 
)
Calculate a new DcrModel from a set of image residuals.

Parameters
----------
dcrModels : `lsst.pipe.tasks.DcrModel`
    Current model of the true sky after correcting chromatic effects.
residualGeneratorList : `generator` of `numpy.ndarray`
    The residual image for the next subfilter, shifted for DCR.
dcrBBox : `lsst.geom.box.Box2I`
    Sub-region of the coadd which includes a buffer to allow for DCR.
statsCtrl : `lsst.afw.math.StatisticsControl`
    Statistics control object for coadd
gain : `float`
    Relative weight to give the new solution when updating the model.
modelWeights : `numpy.ndarray` or `float`
    A 2D array of weight values that tapers smoothly to zero away from detected sources.
    Set to a placeholder value of 1.0 if ``self.config.useModelWeights`` is False.
refImage : `lsst.afw.image.Image`
    A reference image used to supply the default pixel values.
dcrWeights : `list` of `lsst.afw.image.Image`
    Per-pixel weights for each subfilter.
    Equal to 1/(number of unmasked images contributing to each pixel).

Returns
-------
dcrModel : `lsst.pipe.tasks.DcrModel`
    New model of the true sky after correcting chromatic effects.

Definition at line 915 of file dcrAssembleCoadd.py.

916 gain, modelWeights, refImage, dcrWeights):
917 """Calculate a new DcrModel from a set of image residuals.
918
919 Parameters
920 ----------
921 dcrModels : `lsst.pipe.tasks.DcrModel`
922 Current model of the true sky after correcting chromatic effects.
923 residualGeneratorList : `generator` of `numpy.ndarray`
924 The residual image for the next subfilter, shifted for DCR.
925 dcrBBox : `lsst.geom.box.Box2I`
926 Sub-region of the coadd which includes a buffer to allow for DCR.
928 Statistics control object for coadd
929 gain : `float`
930 Relative weight to give the new solution when updating the model.
931 modelWeights : `numpy.ndarray` or `float`
932 A 2D array of weight values that tapers smoothly to zero away from detected sources.
933 Set to a placeholder value of 1.0 if ``self.config.useModelWeights`` is False.
934 refImage : `lsst.afw.image.Image`
935 A reference image used to supply the default pixel values.
936 dcrWeights : `list` of `lsst.afw.image.Image`
937 Per-pixel weights for each subfilter.
938 Equal to 1/(number of unmasked images contributing to each pixel).
939
940 Returns
941 -------
942 dcrModel : `lsst.pipe.tasks.DcrModel`
943 New model of the true sky after correcting chromatic effects.
944 """
945 newModelImages = []
946 for subfilter, model in enumerate(dcrModels):
947 residualsList = [next(residualGenerator) for residualGenerator in residualGeneratorList]
948 residual = np.sum(residualsList, axis=0)
949 residual *= dcrWeights[subfilter][dcrBBox].array
950 # `MaskedImage`s only support in-place addition, so rename for readability
951 newModel = model[dcrBBox].clone()
952 newModel.array += residual
953 # Catch any invalid values
954 badPixels = ~np.isfinite(newModel.array)
955 newModel.array[badPixels] = model[dcrBBox].array[badPixels]
956 if self.config.regularizeModelIterations > 0:
957 dcrModels.regularizeModelIter(subfilter, newModel, dcrBBox,
958 self.config.regularizeModelIterations,
959 self.config.regularizationWidth)
960 newModelImages.append(newModel)
961 if self.config.regularizeModelFrequency > 0:
962 dcrModels.regularizeModelFreq(newModelImages, dcrBBox, statsCtrl,
963 self.config.regularizeModelFrequency,
964 self.config.regularizationWidth)
965 dcrModels.conditionDcrModel(newModelImages, dcrBBox, gain=gain)
966 self.applyModelWeights(newModelImages, refImage[dcrBBox], modelWeights)
967 return DcrModel(newModelImages, dcrModels.filter, dcrModels.effectiveWavelength,
968 dcrModels.bandwidth, dcrModels.psf,
969 dcrModels.mask, dcrModels.variance)
970

◆ selectCoaddPsf()

def lsst.pipe.tasks.dcrAssembleCoadd.selectCoaddPsf (   self,
  templateCoadd,
  warpRefList 
)
Compute the PSF of the coadd from the exposures with the best seeing.

Parameters
----------
templateCoadd : `lsst.afw.image.ExposureF`
    The initial coadd exposure before accounting for DCR.
warpRefList : `list` of `lsst.daf.butler.DeferredDatasetHandle` or
    `lsst.daf.persistence.ButlerDataRef`
    The data references to the input warped exposures.

Returns
-------
psf : `lsst.meas.algorithms.CoaddPsf`
    The average PSF of the input exposures with the best seeing.

Definition at line 1314 of file dcrAssembleCoadd.py.

1314 def selectCoaddPsf(self, templateCoadd, warpRefList):
1315 """Compute the PSF of the coadd from the exposures with the best seeing.
1316
1317 Parameters
1318 ----------
1319 templateCoadd : `lsst.afw.image.ExposureF`
1320 The initial coadd exposure before accounting for DCR.
1321 warpRefList : `list` of `lsst.daf.butler.DeferredDatasetHandle` or
1323 The data references to the input warped exposures.
1324
1325 Returns
1326 -------
1328 The average PSF of the input exposures with the best seeing.
1329 """
1330 sigma2fwhm = 2.*np.sqrt(2.*np.log(2.))
1331 tempExpName = self.getTempExpDatasetName(self.warpType)
1332 # Note: ``ccds`` is a `lsst.afw.table.ExposureCatalog` with one entry per ccd and per visit
1333 # If there are multiple ccds, it will have that many times more elements than ``warpExpRef``
1334 ccds = templateCoadd.getInfo().getCoaddInputs().ccds
1335 templatePsf = templateCoadd.getPsf()
1336 # Just need a rough estimate; average positions are fine
1337 templateAvgPos = templatePsf.getAveragePosition()
1338 psfRefSize = templatePsf.computeShape(templateAvgPos).getDeterminantRadius()*sigma2fwhm
1339 psfSizes = np.zeros(len(ccds))
1340 ccdVisits = np.array(ccds["visit"])
1341 for warpExpRef in warpRefList:
1342 if isinstance(warpExpRef, DeferredDatasetHandle):
1343 # Gen 3 API
1344 psf = warpExpRef.get(component="psf")
1345 else:
1346 # Gen 2 API. Delete this when Gen 2 retired
1347 psf = warpExpRef.get(tempExpName).getPsf()
1348 visit = warpExpRef.dataId["visit"]
1349 psfAvgPos = psf.getAveragePosition()
1350 psfSize = psf.computeShape(psfAvgPos).getDeterminantRadius()*sigma2fwhm
1351 psfSizes[ccdVisits == visit] = psfSize
1352 # Note that the input PSFs include DCR, which should be absent from the DcrCoadd
1353 # The selected PSFs are those that have a FWHM less than or equal to the smaller
1354 # of the mean or median FWHM of the input exposures.
1355 sizeThreshold = min(np.median(psfSizes), psfRefSize)
1356 goodPsfs = psfSizes <= sizeThreshold
1357 psf = measAlg.CoaddPsf(ccds[goodPsfs], templateCoadd.getWcs(),
1358 self.config.coaddPsf.makeControl())
1359 return psf
int min
CoaddPsf is the Psf derived to be used for non-PSF-matched Coadd images.
Definition: CoaddPsf.h:58
def selectCoaddPsf(self, templateCoadd, warpRefList)

◆ stackCoadd()

def lsst.pipe.tasks.dcrAssembleCoadd.stackCoadd (   self,
  dcrCoadds 
)
Add a list of sub-band coadds together.

Parameters
----------
dcrCoadds : `list` of `lsst.afw.image.ExposureF`
    A list of coadd exposures, each exposure containing
    the model for one subfilter.

Returns
-------
coaddExposure : `lsst.afw.image.ExposureF`
    A single coadd exposure that is the sum of the sub-bands.

Definition at line 1056 of file dcrAssembleCoadd.py.

1056 def stackCoadd(self, dcrCoadds):
1057 """Add a list of sub-band coadds together.
1058
1059 Parameters
1060 ----------
1061 dcrCoadds : `list` of `lsst.afw.image.ExposureF`
1062 A list of coadd exposures, each exposure containing
1063 the model for one subfilter.
1064
1065 Returns
1066 -------
1067 coaddExposure : `lsst.afw.image.ExposureF`
1068 A single coadd exposure that is the sum of the sub-bands.
1069 """
1070 coaddExposure = dcrCoadds[0].clone()
1071 for coadd in dcrCoadds[1:]:
1072 coaddExposure.maskedImage += coadd.maskedImage
1073 return coaddExposure
1074

Variable Documentation

◆ goodPix

tuple lsst.pipe.tasks.dcrAssembleCoadd.goodPix = dcrWeights[0].array > weightsThreshold

Definition at line 798 of file dcrAssembleCoadd.py.

◆ weightsThreshold

int lsst.pipe.tasks.dcrAssembleCoadd.weightsThreshold = 1.

Definition at line 797 of file dcrAssembleCoadd.py.