LSSTApplications  10.0-2-g4f67435,11.0.rc2+1,11.0.rc2+12,11.0.rc2+3,11.0.rc2+4,11.0.rc2+5,11.0.rc2+6,11.0.rc2+7,11.0.rc2+8
LSSTDataManagementBasePackage
assembleCoadd.py
Go to the documentation of this file.
1 from __future__ import division, absolute_import
2 #
3 # LSST Data Management System
4 # Copyright 2008, 2009, 2010, 2011, 2012 LSST Corporation.
5 #
6 # This product includes software developed by the
7 # LSST Project (http://www.lsst.org/).
8 #
9 # This program is free software: you can redistribute it and/or modify
10 # it under the terms of the GNU General Public License as published by
11 # the Free Software Foundation, either version 3 of the License, or
12 # (at your option) any later version.
13 #
14 # This program is distributed in the hope that it will be useful,
15 # but WITHOUT ANY WARRANTY; without even the implied warranty of
16 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17 # GNU General Public License for more details.
18 #
19 # You should have received a copy of the LSST License Statement and
20 # the GNU General Public License along with this program. If not,
21 # see <http://www.lsstcorp.org/LegalNotices/>.
22 #
23 import numpy
24 import lsst.pex.config as pexConfig
25 import lsst.afw.geom as afwGeom
26 import lsst.afw.image as afwImage
27 import lsst.afw.math as afwMath
28 import lsst.coadd.utils as coaddUtils
29 import lsst.pipe.base as pipeBase
30 import lsst.meas.algorithms as measAlg
31 from .coaddBase import CoaddBaseTask, SelectDataIdContainer
32 from .interpImage import InterpImageTask
33 from .matchBackgrounds import MatchBackgroundsTask
34 from .scaleZeroPoint import ScaleZeroPointTask
35 from .coaddHelpers import groupPatchExposures, getGroupDataRef
36 
37 __all__ = ["AssembleCoaddTask"]
38 
39 class AssembleCoaddConfig(CoaddBaseTask.ConfigClass):
40  subregionSize = pexConfig.ListField(
41  dtype = int,
42  doc = "Width, height of stack subregion size; " \
43  "make small enough that a full stack of images will fit into memory at once.",
44  length = 2,
45  default = (2000, 2000),
46  )
47  doSigmaClip = pexConfig.Field(
48  dtype = bool,
49  doc = "Perform sigma clipped outlier rejection? If False then compute a simple mean.",
50  default = True,
51  )
52  sigmaClip = pexConfig.Field(
53  dtype = float,
54  doc = "Sigma for outlier rejection; ignored if doSigmaClip false.",
55  default = 3.0,
56  )
57  clipIter = pexConfig.Field(
58  dtype = int,
59  doc = "Number of iterations of outlier rejection; ignored if doSigmaClip false.",
60  default = 2,
61  )
62  scaleZeroPoint = pexConfig.ConfigurableField(
63  target = ScaleZeroPointTask,
64  doc = "Task to adjust the photometric zero point of the coadd temp exposures",
65  )
66  doInterp = pexConfig.Field(
67  doc = "Interpolate over NaN pixels? Also extrapolate, if necessary, but the results are ugly.",
68  dtype = bool,
69  default = True,
70  )
71  interpImage = pexConfig.ConfigurableField(
72  target = InterpImageTask,
73  doc = "Task to interpolate (and extrapolate) over NaN pixels",
74  )
75  matchBackgrounds = pexConfig.ConfigurableField(
76  target = MatchBackgroundsTask,
77  doc = "Task to match backgrounds",
78  )
79  maxMatchResidualRatio = pexConfig.Field(
80  doc = "Maximum ratio of the mean squared error of the background matching model to the variance " \
81  "of the difference in backgrounds",
82  dtype = float,
83  default = 1.1
84  )
85  maxMatchResidualRMS = pexConfig.Field(
86  doc = "Maximum RMS of residuals of the background offset fit in matchBackgrounds.",
87  dtype = float,
88  default = 1.0
89  )
90  doWrite = pexConfig.Field(
91  doc = "Persist coadd?",
92  dtype = bool,
93  default = True,
94  )
95  doMatchBackgrounds = pexConfig.Field(
96  doc = "Match backgrounds of coadd temp exposures before coadding them? " \
97  "If False, the coadd temp expsosures must already have been background subtracted or matched",
98  dtype = bool,
99  default = True,
100  )
101  autoReference = pexConfig.Field(
102  doc = "Automatically select the coadd temp exposure to use as a reference for background matching? " \
103  "Ignored if doMatchBackgrounds false. " \
104  "If False you must specify the reference temp exposure as the data Id",
105  dtype = bool,
106  default = True,
107  )
108  maskPropagationThresholds = pexConfig.DictField(
109  keytype = str,
110  itemtype = float,
111  doc = ("Threshold (in fractional weight) of rejection at which we propagate a mask plane to "
112  "the coadd; that is, we set the mask bit on the coadd if the fraction the rejected frames "
113  "would have contributed exceeds this value."),
114  default = {"SAT": 0.1},
115  )
116 
117  def setDefaults(self):
118  CoaddBaseTask.ConfigClass.setDefaults(self)
119  self.badMaskPlanes = ["NO_DATA", "BAD", "CR",]
120 
121 
122 class AssembleCoaddTask(CoaddBaseTask):
123  """Assemble a coadd from a set of coaddTempExp
124  """
125  ConfigClass = AssembleCoaddConfig
126  _DefaultName = "assembleCoadd"
127 
128  def __init__(self, *args, **kwargs):
129  CoaddBaseTask.__init__(self, *args, **kwargs)
130  self.makeSubtask("interpImage")
131  self.makeSubtask("matchBackgrounds")
132  self.makeSubtask("scaleZeroPoint")
133 
134  @pipeBase.timeMethod
135  def run(self, dataRef, selectDataList=[]):
136  """Assemble a coadd from a set of coaddTempExp
137 
138  The coadd is computed as a mean with optional outlier rejection.
139 
140  AssumeCoaddTask performs coaddition of "coadd temporary exposures" ("coaddTempExps"). Each
141  coaddTempExp is the size of a patch and contains data for one run, visit or (for a non-mosaic camera)
142  exposure. The coaddTempExps to coadd are selected from the provided selectDataList based on their
143  overlap with the patch specified by dataRef.
144 
145  By default, coaddTempExps contain backgrounds and hence require config.doMatchBackgrounds=True.
146  Background-subtracted coaddTempExps can be coadded by setting config.doMatchBackgrounds=False.
147 
148  When background matching is enabled, the task may be configured to automatically select a reference
149  exposure (config.autoReference=True). If this is not done, then we require that the input dataRef
150  provides access to a coaddTempExp (dataset type coaddName + 'Coadd_tempExp') which is used as the
151  reference exposure.
152 
153  @param dataRef: Data reference defining the patch for coaddition and the reference coaddTempExp
154  (if config.autoReference=False). Used to access the following data products:
155  - [in] self.config.coaddName + "Coadd_skyMap"
156  - [in] self.config.coaddName + "Coadd_tempExp" (optionally)
157  - [out] self.config.coaddName + "Coadd"
158  @param selectDataList[in]: List of data references to coaddTempExp. Data to be coadded will be
159  selected from this list based on overlap with the patch defined by dataRef.
160 
161  @return: a pipeBase.Struct with fields:
162  - coaddExposure: coadded exposure
163  """
164  skyInfo = self.getSkyInfo(dataRef)
165  calExpRefList = self.selectExposures(dataRef, skyInfo, selectDataList=selectDataList)
166  if len(calExpRefList) == 0:
167  self.log.warn("No exposures to coadd")
168  return
169  self.log.info("Coadding %d exposures" % len(calExpRefList))
170 
171  tempExpRefList = self.getTempExpRefList(dataRef, calExpRefList)
172  inputData = self.prepareInputs(tempExpRefList)
173  self.log.info("Found %d %s" % (len(inputData.tempExpRefList), self.getTempExpDatasetName()))
174  if len(inputData.tempExpRefList) == 0:
175  self.log.warn("No coadd temporary exposures found")
176  return
177  if self.config.doMatchBackgrounds:
178  refImageScaler = self.getBackgroundReferenceScaler(dataRef)
179  inputData = self.backgroundMatching(inputData, dataRef, refImageScaler)
180  if len(inputData.tempExpRefList) == 0:
181  self.log.warn("No valid background models")
182  return
183 
184  coaddExp = self.assemble(skyInfo, inputData.tempExpRefList, inputData.imageScalerList,
185  inputData.weightList,
186  inputData.backgroundInfoList if self.config.doMatchBackgrounds else None)
187  if self.config.doMatchBackgrounds:
188  self.addBackgroundMatchingMetadata(coaddExp, inputData.tempExpRefList,
189  inputData.backgroundInfoList)
190 
191  if self.config.doInterp:
192  self.interpImage.interpolateOnePlane(maskedImage=coaddExp.getMaskedImage(), planeName="NO_DATA")
193  # The variance must be positive; work around for DM-3201.
194  varArray = coaddExp.getMaskedImage().getVariance().getArray()
195  varArray[:] = numpy.where(varArray > 0, varArray, numpy.inf)
196 
197  if self.config.doWrite:
198  self.writeCoaddOutput(dataRef, coaddExp)
199 
200  return pipeBase.Struct(coaddExposure=coaddExp)
201 
202  def getTempExpRefList(self, patchRef, calExpRefList):
203  """Generate list of coaddTempExp data references
204 
205  @param patchRef: Data reference for patch
206  @param calExpRefList: List of data references for input calexps
207  @return List of coaddTempExp data references
208  """
209  butler = patchRef.getButler()
210  groupData = groupPatchExposures(patchRef, calExpRefList, self.getCoaddDatasetName(),
211  self.getTempExpDatasetName())
212  tempExpRefList = [getGroupDataRef(butler, self.getTempExpDatasetName(), g, groupData.keys) for
213  g in groupData.groups.keys()]
214  return tempExpRefList
215 
216  def getBackgroundReferenceScaler(self, dataRef):
217  """Construct an image scaler for the background reference frame
218 
219  If there is no reference frame ('autoReference'), then this is a no-op
220  and None is returned
221 
222  @param dataRef: Data reference for the background reference frame, or None
223  @return image scaler, or None
224  """
225  if self.config.autoReference:
226  return None
227 
228  # We've been given the data reference
229  dataset = self.getTempExpDatasetName()
230  if not dataRef.datasetExists(dataset):
231  raise RuntimeError("Could not find reference exposure %s %s." % (dataset, dataRef.dataId))
232 
233  refExposure = dataRef.get(self.getTempExpDatasetName(), immediate=True)
234  refImageScaler = self.scaleZeroPoint.computeImageScaler(
235  exposure = refExposure,
236  dataRef = dataRef,
237  )
238  return refImageScaler
239 
240 
241  def prepareInputs(self, refList):
242  """Prepare the input warps for coaddition
243 
244  This involves measuring weights and constructing image scalers
245  for each of the inputs.
246 
247  @param refList: List of data references to tempExp
248  @return Struct:
249  - tempExprefList: List of data references to tempExp
250  - weightList: List of weightings
251  - imageScalerList: List of image scalers
252  """
253  statsCtrl = afwMath.StatisticsControl()
254  statsCtrl.setNumSigmaClip(self.config.sigmaClip)
255  statsCtrl.setNumIter(self.config.clipIter)
256  statsCtrl.setAndMask(self.getBadPixelMask())
257  statsCtrl.setNanSafe(True)
258 
259  # compute tempExpRefList: a list of tempExpRef that actually exist
260  # and weightList: a list of the weight of the associated coadd tempExp
261  # and imageScalerList: a list of scale factors for the associated coadd tempExp
262  tempExpRefList = []
263  weightList = []
264  imageScalerList = []
265  tempExpName = self.getTempExpDatasetName()
266  for tempExpRef in refList:
267  if not tempExpRef.datasetExists(tempExpName):
268  self.log.warn("Could not find %s %s; skipping it" % (tempExpName, tempExpRef.dataId))
269  continue
270 
271  tempExp = tempExpRef.get(tempExpName, immediate=True)
272  maskedImage = tempExp.getMaskedImage()
273  imageScaler = self.scaleZeroPoint.computeImageScaler(
274  exposure = tempExp,
275  dataRef = tempExpRef,
276  )
277  try:
278  imageScaler.scaleMaskedImage(maskedImage)
279  except Exception, e:
280  self.log.warn("Scaling failed for %s (skipping it): %s" % (tempExpRef.dataId, e))
281  continue
282  statObj = afwMath.makeStatistics(maskedImage.getVariance(), maskedImage.getMask(),
283  afwMath.MEANCLIP, statsCtrl)
284  meanVar, meanVarErr = statObj.getResult(afwMath.MEANCLIP);
285  weight = 1.0 / float(meanVar)
286  if not numpy.isfinite(weight):
287  self.log.warn("Non-finite weight for %s: skipping" % (tempExpRef.dataId,))
288  continue
289  self.log.info("Weight of %s %s = %0.3f" % (tempExpName, tempExpRef.dataId, weight))
290 
291  del maskedImage
292  del tempExp
293 
294  tempExpRefList.append(tempExpRef)
295  weightList.append(weight)
296  imageScalerList.append(imageScaler)
297 
298  return pipeBase.Struct(tempExpRefList=tempExpRefList, weightList=weightList,
299  imageScalerList=imageScalerList)
300 
301 
302  def backgroundMatching(self, inputData, refExpDataRef=None, refImageScaler=None):
303  """Perform background matching on the prepared inputs
304 
305  If no reference is provided, the background matcher will select one.
306 
307  This method returns a new inputData Struct that can replace the original.
308 
309  @param inputData: Struct from prepareInputs() with tempExpRefList, weightList, imageScalerList
310  @param refExpDataRef: Data reference for background reference tempExp, or None
311  @param refImageScaler: Image scaler for background reference tempExp, or None
312  @return Struct:
313  - tempExprefList: List of data references to tempExp
314  - weightList: List of weightings
315  - imageScalerList: List of image scalers
316  - backgroundInfoList: result from background matching
317  """
318  try:
319  backgroundInfoList = self.matchBackgrounds.run(
320  expRefList = inputData.tempExpRefList,
321  imageScalerList = inputData.imageScalerList,
322  refExpDataRef = refExpDataRef if not self.config.autoReference else None,
323  refImageScaler = refImageScaler,
324  expDatasetType = self.getTempExpDatasetName(),
325  ).backgroundInfoList
326  except Exception, e:
327  self.log.fatal("Cannot match backgrounds: %s" % (e))
328  raise pipeBase.TaskError("Background matching failed.")
329 
330  newWeightList = []
331  newTempExpRefList = []
332  newBackgroundStructList = []
333  newScaleList = []
334  # the number of good backgrounds may be < than len(tempExpList)
335  # sync these up and correct the weights
336  for tempExpRef, bgInfo, scaler, weight in zip(inputData.tempExpRefList, backgroundInfoList,
337  inputData.imageScalerList, inputData.weightList):
338  if not bgInfo.isReference:
339  # skip exposure if it has no backgroundModel
340  # or if fit was bad
341  if (bgInfo.backgroundModel is None):
342  self.log.info("No background offset model available for %s: skipping"%(
343  tempExpRef.dataId))
344  continue
345  try:
346  varianceRatio = bgInfo.matchedMSE / bgInfo.diffImVar
347  except Exception, e:
348  self.log.info("MSE/Var ratio not calculable (%s) for %s: skipping" %
349  (e, tempExpRef.dataId,))
350  continue
351  if not numpy.isfinite(varianceRatio):
352  self.log.info("MSE/Var ratio not finite (%.2f / %.2f) for %s: skipping" %
353  (bgInfo.matchedMSE, bgInfo.diffImVar,
354  tempExpRef.dataId,))
355  continue
356  elif (varianceRatio > self.config.maxMatchResidualRatio):
357  self.log.info("Bad fit. MSE/Var ratio %.2f > %.2f for %s: skipping" % (
358  varianceRatio, self.config.maxMatchResidualRatio, tempExpRef.dataId,))
359  continue
360  elif ( bgInfo.fitRMS > self.config.maxMatchResidualRMS):
361  self.log.info("Bad fit. RMS %.2f > %.2f for %s: skipping" % (
362  bgInfo.fitRMS, self.config.maxMatchResidualRMS, tempExpRef.dataId,))
363  continue
364  newWeightList.append(1 / (1 / weight + bgInfo.fitRMS**2))
365  newTempExpRefList.append(tempExpRef)
366  newBackgroundStructList.append(bgInfo)
367  newScaleList.append(scaler)
368 
369  return pipeBase.Struct(tempExpRefList=newTempExpRefList, weightList=newWeightList,
370  imageScalerList=newScaleList, backgroundInfoList=newBackgroundStructList)
371 
372  def assemble(self, skyInfo, tempExpRefList, imageScalerList, weightList, bgInfoList=None):
373  """Assemble a coadd from input warps
374 
375  The assembly is performed over small areas on the image at a time, to
376  conserve memory usage.
377 
378  @param skyInfo: Patch geometry information, from getSkyInfo
379  @param tempExpRefList: List of data references to tempExp
380  @param imageScalerList: List of image scalers
381  @param weightList: List of weights
382  @param bgInfoList: List of background data from background matching
383  @return coadded exposure
384  """
385  tempExpName = self.getTempExpDatasetName()
386  self.log.info("Assembling %s %s" % (len(tempExpRefList), tempExpName))
387 
388  statsCtrl = afwMath.StatisticsControl()
389  statsCtrl.setNumSigmaClip(self.config.sigmaClip)
390  statsCtrl.setNumIter(self.config.clipIter)
391  statsCtrl.setAndMask(self.getBadPixelMask())
392  statsCtrl.setNanSafe(True)
393  statsCtrl.setWeighted(True)
394  statsCtrl.setCalcErrorFromInputVariance(True)
395  for plane, threshold in self.config.maskPropagationThresholds.items():
396  bit = afwImage.MaskU.getMaskPlane(plane)
397  statsCtrl.setMaskPropagationThreshold(bit, threshold)
398 
399  if self.config.doSigmaClip:
400  statsFlags = afwMath.MEANCLIP
401  else:
402  statsFlags = afwMath.MEAN
403 
404  if bgInfoList is None:
405  bgInfoList = [None]*len(tempExpRefList)
406 
407  coaddExposure = afwImage.ExposureF(skyInfo.bbox, skyInfo.wcs)
408  coaddExposure.setCalib(self.scaleZeroPoint.getCalib())
409  coaddExposure.getInfo().setCoaddInputs(self.inputRecorder.makeCoaddInputs())
410  self.assembleMetadata(coaddExposure, tempExpRefList, weightList)
411  coaddMaskedImage = coaddExposure.getMaskedImage()
412  subregionSizeArr = self.config.subregionSize
413  subregionSize = afwGeom.Extent2I(subregionSizeArr[0], subregionSizeArr[1])
414  for subBBox in _subBBoxIter(skyInfo.bbox, subregionSize):
415  try:
416  self.assembleSubregion(coaddExposure, subBBox, tempExpRefList, imageScalerList,
417  weightList, bgInfoList, statsFlags, statsCtrl)
418  except Exception, e:
419  self.log.fatal("Cannot compute coadd %s: %s" % (subBBox, e,))
420 
421  coaddUtils.setCoaddEdgeBits(coaddMaskedImage.getMask(), coaddMaskedImage.getVariance())
422 
423  return coaddExposure
424 
425  def assembleMetadata(self, coaddExposure, tempExpRefList, weightList):
426  """Set the metadata for the coadd
427 
428  This basic implementation simply sets the filter from the
429  first input.
430 
431  @param coaddExposure: The target image for the coadd
432  @param tempExpRefList: List of data references to tempExp
433  @param weightList: List of weights
434  """
435  tempExpName = self.getTempExpDatasetName()
436  # We load a single pixel of each coaddTempExp, because we just want to get at the metadata
437  # (and we need more than just the PropertySet that contains the header). See #2777.
439  first = True
440  coaddInputs = coaddExposure.getInfo().getCoaddInputs()
441  for tempExpRef, weight in zip(tempExpRefList, weightList):
442  tempExp = tempExpRef.get(tempExpName + "_sub", bbox=bbox, imageOrigin="LOCAL", immediate=True)
443  if first:
444  coaddExposure.setFilter(tempExp.getFilter())
445  first = False
446  self.inputRecorder.addVisitToCoadd(coaddInputs, tempExp, weight)
447  coaddInputs.visits.sort()
448  if self.config.doPsfMatch:
449  psf = self.config.modelPsf.apply()
450  else:
451  psf = measAlg.CoaddPsf(coaddInputs.ccds, coaddExposure.getWcs())
452  coaddExposure.setPsf(psf)
453  apCorrMap = measAlg.makeCoaddApCorrMap(coaddInputs.ccds, coaddExposure.getBBox(afwImage.PARENT),
454  coaddExposure.getWcs())
455  coaddExposure.getInfo().setApCorrMap(apCorrMap)
456 
457  def assembleSubregion(self, coaddExposure, bbox, tempExpRefList, imageScalerList, weightList,
458  bgInfoList, statsFlags, statsCtrl):
459  """Assemble the coadd for a sub-region
460 
461  @param coaddExposure: The target image for the coadd
462  @param bbox: Sub-region to coadd
463  @param tempExpRefList: List of data reference to tempExp
464  @param imageScalerList: List of image scalers
465  @param weightList: List of weights
466  @param bgInfoList: List of background data from background matching
467  @param statsFlags: Statistic for coadd
468  @param statsCtrl: Statistics control object for coadd
469  """
470  self.log.logdebug("Computing coadd over %s" % bbox)
471  tempExpName = self.getTempExpDatasetName()
472  coaddMaskedImage = coaddExposure.getMaskedImage()
473  coaddView = afwImage.MaskedImageF(coaddMaskedImage, bbox, afwImage.PARENT, False)
474  maskedImageList = afwImage.vectorMaskedImageF() # [] is rejected by afwMath.statisticsStack
475  for tempExpRef, imageScaler, bgInfo in zip(tempExpRefList, imageScalerList, bgInfoList):
476  exposure = tempExpRef.get(tempExpName + "_sub", bbox=bbox)
477  maskedImage = exposure.getMaskedImage()
478  imageScaler.scaleMaskedImage(maskedImage)
479 
480  if self.config.doMatchBackgrounds and not bgInfo.isReference:
481  backgroundModel = bgInfo.backgroundModel
482  backgroundImage = backgroundModel.getImage() if \
483  self.matchBackgrounds.config.usePolynomial else \
484  backgroundModel.getImageF()
485  backgroundImage.setXY0(coaddMaskedImage.getXY0())
486  maskedImage += backgroundImage.Factory(backgroundImage, bbox, afwImage.PARENT, False)
487  var = maskedImage.getVariance()
488  var += (bgInfo.fitRMS)**2
489 
490  maskedImageList.append(maskedImage)
491 
492  with self.timer("stack"):
493  coaddSubregion = afwMath.statisticsStack(
494  maskedImageList, statsFlags, statsCtrl, weightList)
495 
496  coaddView <<= coaddSubregion
497 
498 
499  def addBackgroundMatchingMetadata(self, coaddExposure, tempExpRefList, backgroundInfoList):
500  """Add metadata from the background matching to the coadd
501 
502  @param coaddExposure: Coadd
503  @param backgroundInfoList: List of background info, results from background matching
504  """
505  self.log.info("Adding exposure information to metadata")
506  metadata = coaddExposure.getMetadata()
507  metadata.addString("CTExp_SDQA1_DESCRIPTION",
508  "Background matching: Ratio of matchedMSE / diffImVar")
509  for ind, (tempExpRef, backgroundInfo) in enumerate(zip(tempExpRefList, backgroundInfoList)):
510  tempExpStr = '&'.join('%s=%s' % (k,v) for k,v in tempExpRef.dataId.items())
511  if backgroundInfo.isReference:
512  metadata.addString("ReferenceExp_ID", tempExpStr)
513  else:
514  metadata.addString("CTExp_ID_%d" % (ind), tempExpStr)
515  metadata.addDouble("CTExp_SDQA1_%d" % (ind),
516  backgroundInfo.matchedMSE/backgroundInfo.diffImVar)
517  metadata.addDouble("CTExp_SDQA2_%d" % (ind),
518  backgroundInfo.fitRMS)
519  @classmethod
521  """Create an argument parser
522  """
523  parser = pipeBase.ArgumentParser(name=cls._DefaultName)
524  parser.add_id_argument("--id", cls.ConfigClass().coaddName + "Coadd_tempExp",
525  help="data ID, e.g. --id tract=12345 patch=1,2",
526  ContainerClass=AssembleCoaddDataIdContainer)
527  parser.add_id_argument("--selectId", "calexp", help="data ID, e.g. --selectId visit=6789 ccd=0..9",
528  ContainerClass=SelectDataIdContainer)
529  return parser
530 
531 
532 
533 def _subBBoxIter(bbox, subregionSize):
534  """Iterate over subregions of a bbox
535 
536  @param[in] bbox: bounding box over which to iterate: afwGeom.Box2I
537  @param[in] subregionSize: size of sub-bboxes
538 
539  @return subBBox: next sub-bounding box of size subregionSize or smaller;
540  each subBBox is contained within bbox, so it may be smaller than subregionSize at the edges of bbox,
541  but it will never be empty
542  """
543  if bbox.isEmpty():
544  raise RuntimeError("bbox %s is empty" % (bbox,))
545  if subregionSize[0] < 1 or subregionSize[1] < 1:
546  raise RuntimeError("subregionSize %s must be nonzero" % (subregionSize,))
547 
548  for rowShift in range(0, bbox.getHeight(), subregionSize[1]):
549  for colShift in range(0, bbox.getWidth(), subregionSize[0]):
550  subBBox = afwGeom.Box2I(bbox.getMin() + afwGeom.Extent2I(colShift, rowShift), subregionSize)
551  subBBox.clip(bbox)
552  if subBBox.isEmpty():
553  raise RuntimeError("Bug: empty bbox! bbox=%s, subregionSize=%s, colShift=%s, rowShift=%s" % \
554  (bbox, subregionSize, colShift, rowShift))
555  yield subBBox
556 
557 
558 
559 class AssembleCoaddDataIdContainer(pipeBase.DataIdContainer):
560  """A version of lsst.pipe.base.DataIdContainer specialized for assembleCoadd.
561  """
562  def makeDataRefList(self, namespace):
563  """Make self.refList from self.idList.
564 
565  Interpret the config.doMatchBackgrounds, config.autoReference,
566  and whether a visit/run supplied.
567  If a visit/run is supplied, config.autoReference is automatically set to False.
568  if config.doMatchBackgrounds == false, then a visit/run will be ignored if accidentally supplied.
569 
570  """
571  keysCoadd = namespace.butler.getKeys(datasetType=namespace.config.coaddName + "Coadd",
572  level=self.level)
573  keysCoaddTempExp = namespace.butler.getKeys(datasetType=namespace.config.coaddName + "Coadd_tempExp",
574  level=self.level)
575 
576  if namespace.config.doMatchBackgrounds:
577  if namespace.config.autoReference: #matcher will pick it's own reference image
578  datasetType = namespace.config.coaddName + "Coadd"
579  validKeys = keysCoadd
580  else:
581  datasetType = namespace.config.coaddName + "Coadd_tempExp"
582  validKeys = keysCoaddTempExp
583  else: #bkg subtracted coadd
584  datasetType = namespace.config.coaddName + "Coadd"
585  validKeys = keysCoadd
586 
587  for dataId in self.idList:
588  # tract and patch are required
589  for key in validKeys:
590  if key not in dataId:
591  raise RuntimeError("--id must include " + key)
592 
593  for key in dataId: # check if users supplied visit/run
594  if (key not in keysCoadd) and (key in keysCoaddTempExp): #user supplied a visit/run
595  if namespace.config.autoReference:
596  # user probably meant: autoReference = False
597  namespace.config.autoReference = False
598  datasetType = namespace.config.coaddName + "Coadd_tempExp"
599  print "Switching config.autoReference to False; applies only to background Matching."
600  break
601 
602  dataRef = namespace.butler.dataRef(
603  datasetType = datasetType,
604  dataId = dataId,
605  )
606  self.refList.append(dataRef)
607 
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
An integer coordinate rectangle.
Definition: Box.h:53
Pass parameters to a Statistics objectA class to pass parameters which control how the stats are calc...
Definition: Statistics.h:92
CoaddPsf is the Psf derived to be used for non-PSF-matched Coadd images.
Definition: CoaddPsf.h:45
Statistics makeStatistics(afwImage::Mask< afwImage::MaskPixel > const &msk, int const flags, StatisticsControl const &sctrl)
Specialization to handle Masks.
Definition: Statistics.cc:1082
lsst::afw::image::MaskedImage< PixelT >::Ptr statisticsStack(lsst::afw::image::MaskedImage< PixelT > const &image, Property flags, char dimension, StatisticsControl const &sctrl=StatisticsControl())
A function to compute statistics on the rows or columns of an image.