LSSTApplications  18.1.0
LSSTDataManagementBasePackage
makeCoaddTempExp.py
Go to the documentation of this file.
1 #
2 # LSST Data Management System
3 # Copyright 2008, 2009, 2010, 2011, 2012 LSST Corporation.
4 #
5 # This product includes software developed by the
6 # LSST Project (http://www.lsst.org/).
7 #
8 # This program is free software: you can redistribute it and/or modify
9 # it under the terms of the GNU General Public License as published by
10 # the Free Software Foundation, either version 3 of the License, or
11 # (at your option) any later version.
12 #
13 # This program is distributed in the hope that it will be useful,
14 # but WITHOUT ANY WARRANTY; without even the implied warranty of
15 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16 # GNU General Public License for more details.
17 #
18 # You should have received a copy of the LSST License Statement and
19 # the GNU General Public License along with this program. If not,
20 # see <http://www.lsstcorp.org/LegalNotices/>.
21 #
22 import numpy
23 
24 import lsst.pex.config as pexConfig
25 import lsst.daf.persistence as dafPersist
26 import lsst.afw.image as afwImage
27 import lsst.coadd.utils as coaddUtils
28 import lsst.pipe.base as pipeBase
29 import lsst.log as log
30 from lsst.meas.algorithms import CoaddPsf, CoaddPsfConfig
31 from .coaddBase import CoaddBaseTask, makeSkyInfo
32 from .warpAndPsfMatch import WarpAndPsfMatchTask
33 from .coaddHelpers import groupPatchExposures, getGroupDataRef
34 
35 __all__ = ["MakeCoaddTempExpTask", "MakeWarpTask", "MakeWarpConfig"]
36 
37 
38 class MissingExposureError(Exception):
39  """Raised when data cannot be retrieved for an exposure.
40  When processing patches, sometimes one exposure is missing; this lets us
41  distinguish bewteen that case, and other errors.
42  """
43  pass
44 
45 
46 class MakeCoaddTempExpConfig(CoaddBaseTask.ConfigClass):
47  """Config for MakeCoaddTempExpTask
48  """
49  warpAndPsfMatch = pexConfig.ConfigurableField(
50  target=WarpAndPsfMatchTask,
51  doc="Task to warp and PSF-match calexp",
52  )
53  doWrite = pexConfig.Field(
54  doc="persist <coaddName>Coadd_<warpType>Warp",
55  dtype=bool,
56  default=True,
57  )
58  bgSubtracted = pexConfig.Field(
59  doc="Work with a background subtracted calexp?",
60  dtype=bool,
61  default=True,
62  )
63  coaddPsf = pexConfig.ConfigField(
64  doc="Configuration for CoaddPsf",
65  dtype=CoaddPsfConfig,
66  )
67  makeDirect = pexConfig.Field(
68  doc="Make direct Warp/Coadds",
69  dtype=bool,
70  default=True,
71  )
72  makePsfMatched = pexConfig.Field(
73  doc="Make Psf-Matched Warp/Coadd?",
74  dtype=bool,
75  default=False,
76  )
77 
78  doWriteEmptyWarps = pexConfig.Field(
79  dtype=bool,
80  default=False,
81  doc="Write out warps even if they are empty"
82  )
83 
84  hasFakes = pexConfig.Field(
85  doc="Should be set to True if fakes ources have been inserted into the input data.",
86  dtype=bool,
87  default=False,
88  )
89  doApplySkyCorr = pexConfig.Field(dtype=bool, default=False, doc="Apply sky correction?")
90 
91  def validate(self):
92  CoaddBaseTask.ConfigClass.validate(self)
93  if not self.makePsfMatched and not self.makeDirect:
94  raise RuntimeError("At least one of config.makePsfMatched and config.makeDirect must be True")
95  if self.doPsfMatch:
96  # Backwards compatibility.
97  log.warn("Config doPsfMatch deprecated. Setting makePsfMatched=True and makeDirect=False")
98  self.makePsfMatched = True
99  self.makeDirect = False
100 
101  def setDefaults(self):
102  CoaddBaseTask.ConfigClass.setDefaults(self)
103  self.warpAndPsfMatch.psfMatch.kernel.active.kernelSize = self.matchingKernelSize
104 
105 
111 
112 
114  r"""!Warp and optionally PSF-Match calexps onto an a common projection.
115 
116  @anchor MakeCoaddTempExpTask_
117 
118  @section pipe_tasks_makeCoaddTempExp_Contents Contents
119 
120  - @ref pipe_tasks_makeCoaddTempExp_Purpose
121  - @ref pipe_tasks_makeCoaddTempExp_Initialize
122  - @ref pipe_tasks_makeCoaddTempExp_IO
123  - @ref pipe_tasks_makeCoaddTempExp_Config
124  - @ref pipe_tasks_makeCoaddTempExp_Debug
125  - @ref pipe_tasks_makeCoaddTempExp_Example
126 
127  @section pipe_tasks_makeCoaddTempExp_Purpose Description
128 
129  Warp and optionally PSF-Match calexps onto a common projection, by
130  performing the following operations:
131  - Group calexps by visit/run
132  - For each visit, generate a Warp by calling method @ref makeTempExp.
133  makeTempExp loops over the visit's calexps calling @ref WarpAndPsfMatch
134  on each visit
135 
136  The result is a `directWarp` (and/or optionally a `psfMatchedWarp`).
137 
138  @section pipe_tasks_makeCoaddTempExp_Initialize Task Initialization
139 
140  @copydoc \_\_init\_\_
141 
142  This task has one special keyword argument: passing reuse=True will cause
143  the task to skip the creation of warps that are already present in the
144  output repositories.
145 
146  @section pipe_tasks_makeCoaddTempExp_IO Invoking the Task
147 
148  This task is primarily designed to be run from the command line.
149 
150  The main method is `runDataRef`, which takes a single butler data reference for the patch(es)
151  to process.
152 
153  @copydoc run
154 
155  WarpType identifies the types of convolutions applied to Warps (previously CoaddTempExps).
156  Only two types are available: direct (for regular Warps/Coadds) and psfMatched
157  (for Warps/Coadds with homogenized PSFs). We expect to add a third type, likelihood,
158  for generating likelihood Coadds with Warps that have been correlated with their own PSF.
159 
160  @section pipe_tasks_makeCoaddTempExp_Config Configuration parameters
161 
162  See @ref MakeCoaddTempExpConfig and parameters inherited from
163  @link lsst.pipe.tasks.coaddBase.CoaddBaseConfig CoaddBaseConfig @endlink
164 
165  @subsection pipe_tasks_MakeCoaddTempExp_psfMatching Guide to PSF-Matching Configs
166 
167  To make `psfMatchedWarps`, select `config.makePsfMatched=True`. The subtask
168  @link lsst.ip.diffim.modelPsfMatch.ModelPsfMatchTask ModelPsfMatchTask @endlink
169  is responsible for the PSF-Matching, and its config is accessed via `config.warpAndPsfMatch.psfMatch`.
170  The optimal configuration depends on aspects of dataset: the pixel scale, average PSF FWHM and
171  dimensions of the PSF kernel. These configs include the requested model PSF, the matching kernel size,
172  padding of the science PSF thumbnail and spatial sampling frequency of the PSF.
173 
174  *Config Guidelines*: The user must specify the size of the model PSF to which to match by setting
175  `config.modelPsf.defaultFwhm` in units of pixels. The appropriate values depends on science case.
176  In general, for a set of input images, this config should equal the FWHM of the visit
177  with the worst seeing. The smallest it should be set to is the median FWHM. The defaults
178  of the other config options offer a reasonable starting point.
179  The following list presents the most common problems that arise from a misconfigured
180  @link lsst.ip.diffim.modelPsfMatch.ModelPsfMatchTask ModelPsfMatchTask @endlink
181  and corresponding solutions. All assume the default Alard-Lupton kernel, with configs accessed via
182  ```config.warpAndPsfMatch.psfMatch.kernel['AL']```. Each item in the list is formatted as:
183  Problem: Explanation. *Solution*
184 
185  *Troublshooting PSF-Matching Configuration:*
186  - Matched PSFs look boxy: The matching kernel is too small. _Increase the matching kernel size.
187  For example:_
188 
189  config.warpAndPsfMatch.psfMatch.kernel['AL'].kernelSize=27 # default 21
190 
191  Note that increasing the kernel size also increases runtime.
192  - Matched PSFs look ugly (dipoles, quadropoles, donuts): unable to find good solution
193  for matching kernel. _Provide the matcher with more data by either increasing
194  the spatial sampling by decreasing the spatial cell size,_
195 
196  config.warpAndPsfMatch.psfMatch.kernel['AL'].sizeCellX = 64 # default 128
197  config.warpAndPsfMatch.psfMatch.kernel['AL'].sizeCellY = 64 # default 128
198 
199  _or increasing the padding around the Science PSF, for example:_
200 
201  config.warpAndPsfMatch.psfMatch.autoPadPsfTo=1.6 # default 1.4
202 
203  Increasing `autoPadPsfTo` increases the minimum ratio of input PSF dimensions to the
204  matching kernel dimensions, thus increasing the number of pixels available to fit
205  after convolving the PSF with the matching kernel.
206  Optionally, for debugging the effects of padding, the level of padding may be manually
207  controlled by setting turning off the automatic padding and setting the number
208  of pixels by which to pad the PSF:
209 
210  config.warpAndPsfMatch.psfMatch.doAutoPadPsf = False # default True
211  config.warpAndPsfMatch.psfMatch.padPsfBy = 6 # pixels. default 0
212 
213  - Deconvolution: Matching a large PSF to a smaller PSF produces
214  a telltale noise pattern which looks like ripples or a brain.
215  _Increase the size of the requested model PSF. For example:_
216 
217  config.modelPsf.defaultFwhm = 11 # Gaussian sigma in units of pixels.
218 
219  - High frequency (sometimes checkered) noise: The matching basis functions are too small.
220  _Increase the width of the Gaussian basis functions. For example:_
221 
222  config.warpAndPsfMatch.psfMatch.kernel['AL'].alardSigGauss=[1.5, 3.0, 6.0]
223  # from default [0.7, 1.5, 3.0]
224 
225 
226  @section pipe_tasks_makeCoaddTempExp_Debug Debug variables
227 
228  MakeCoaddTempExpTask has no debug output, but its subtasks do.
229 
230  @section pipe_tasks_makeCoaddTempExp_Example A complete example of using MakeCoaddTempExpTask
231 
232  This example uses the package ci_hsc to show how MakeCoaddTempExp fits
233  into the larger Data Release Processing.
234  Set up by running:
235 
236  setup ci_hsc
237  cd $CI_HSC_DIR
238  # if not built already:
239  python $(which scons) # this will take a while
240 
241  The following assumes that `processCcd.py` and `makeSkyMap.py` have previously been run
242  (e.g. by building `ci_hsc` above) to generate a repository of calexps and an
243  output respository with the desired SkyMap. The command,
244 
245  makeCoaddTempExp.py $CI_HSC_DIR/DATA --rerun ci_hsc \
246  --id patch=5,4 tract=0 filter=HSC-I \
247  --selectId visit=903988 ccd=16 --selectId visit=903988 ccd=17 \
248  --selectId visit=903988 ccd=23 --selectId visit=903988 ccd=24 \
249  --config doApplyUberCal=False makePsfMatched=True modelPsf.defaultFwhm=11
250 
251  writes a direct and PSF-Matched Warp to
252  - `$CI_HSC_DIR/DATA/rerun/ci_hsc/deepCoadd/HSC-I/0/5,4/warp-HSC-I-0-5,4-903988.fits` and
253  - `$CI_HSC_DIR/DATA/rerun/ci_hsc/deepCoadd/HSC-I/0/5,4/psfMatchedWarp-HSC-I-0-5,4-903988.fits`
254  respectively.
255 
256  @note PSF-Matching in this particular dataset would benefit from adding
257  `--configfile ./matchingConfig.py` to
258  the command line arguments where `matchingConfig.py` is defined by:
259 
260  echo "
261  config.warpAndPsfMatch.psfMatch.kernel['AL'].kernelSize=27
262  config.warpAndPsfMatch.psfMatch.kernel['AL'].alardSigGauss=[1.5, 3.0, 6.0]" > matchingConfig.py
263 
264 
265  Add the option `--help` to see more options.
266  """
267  ConfigClass = MakeCoaddTempExpConfig
268  _DefaultName = "makeCoaddTempExp"
269 
270  def __init__(self, reuse=False, **kwargs):
271  CoaddBaseTask.__init__(self, **kwargs)
272  self.reuse = reuse
273  self.makeSubtask("warpAndPsfMatch")
274  if self.config.hasFakes:
275  self.calexpType = "fakes_calexp"
276  else:
277  self.calexpType = "calexp"
278 
279  @pipeBase.timeMethod
280  def runDataRef(self, patchRef, selectDataList=[]):
281  """!Produce <coaddName>Coadd_<warpType>Warp images by warping and optionally PSF-matching.
282 
283  @param[in] patchRef: data reference for sky map patch. Must include keys "tract", "patch",
284  plus the camera-specific filter key (e.g. "filter" or "band")
285  @return: dataRefList: a list of data references for the new <coaddName>Coadd_directWarps
286  if direct or both warp types are requested and <coaddName>Coadd_psfMatchedWarps if only psfMatched
287  warps are requested.
288 
289  @warning: this task assumes that all exposures in a warp (coaddTempExp) have the same filter.
290 
291  @warning: this task sets the PhotoCalib of the coaddTempExp to the PhotoCalib of the first calexp
292  with any good pixels in the patch. For a mosaic camera the resulting PhotoCalib should be ignored
293  (assembleCoadd should determine zeropoint scaling without referring to it).
294  """
295  skyInfo = self.getSkyInfo(patchRef)
296 
297  # DataRefs to return are of type *_directWarp unless only *_psfMatchedWarp requested
298  if self.config.makePsfMatched and not self.config.makeDirect:
299  primaryWarpDataset = self.getTempExpDatasetName("psfMatched")
300  else:
301  primaryWarpDataset = self.getTempExpDatasetName("direct")
302 
303  calExpRefList = self.selectExposures(patchRef, skyInfo, selectDataList=selectDataList)
304 
305  if len(calExpRefList) == 0:
306  self.log.warn("No exposures to coadd for patch %s", patchRef.dataId)
307  return None
308  self.log.info("Selected %d calexps for patch %s", len(calExpRefList), patchRef.dataId)
309  calExpRefList = [calExpRef for calExpRef in calExpRefList if calExpRef.datasetExists(self.calexpType)]
310  self.log.info("Processing %d existing calexps for patch %s", len(calExpRefList), patchRef.dataId)
311 
312  groupData = groupPatchExposures(patchRef, calExpRefList, self.getCoaddDatasetName(),
313  primaryWarpDataset)
314  self.log.info("Processing %d warp exposures for patch %s", len(groupData.groups), patchRef.dataId)
315 
316  dataRefList = []
317  for i, (tempExpTuple, calexpRefList) in enumerate(groupData.groups.items()):
318  tempExpRef = getGroupDataRef(patchRef.getButler(), primaryWarpDataset,
319  tempExpTuple, groupData.keys)
320  if self.reuse and tempExpRef.datasetExists(datasetType=primaryWarpDataset, write=True):
321  self.log.info("Skipping makeCoaddTempExp for %s; output already exists.", tempExpRef.dataId)
322  dataRefList.append(tempExpRef)
323  continue
324  self.log.info("Processing Warp %d/%d: id=%s", i, len(groupData.groups), tempExpRef.dataId)
325 
326  # TODO: mappers should define a way to go from the "grouping keys" to a numeric ID (#2776).
327  # For now, we try to get a long integer "visit" key, and if we can't, we just use the index
328  # of the visit in the list.
329  try:
330  visitId = int(tempExpRef.dataId["visit"])
331  except (KeyError, ValueError):
332  visitId = i
333 
334  calExpList = []
335  ccdIdList = []
336  dataIdList = []
337 
338  for calExpInd, calExpRef in enumerate(calexpRefList):
339  self.log.info("Reading calexp %s of %s for Warp id=%s", calExpInd+1, len(calexpRefList),
340  calExpRef.dataId)
341  try:
342  ccdId = calExpRef.get("ccdExposureId", immediate=True)
343  except Exception:
344  ccdId = calExpInd
345  try:
346  # We augment the dataRef here with the tract, which is harmless for loading things
347  # like calexps that don't need the tract, and necessary for meas_mosaic outputs,
348  # which do.
349  calExpRef = calExpRef.butlerSubset.butler.dataRef(self.calexpType,
350  dataId=calExpRef.dataId,
351  tract=skyInfo.tractInfo.getId())
352  calExp = self.getCalibratedExposure(calExpRef, bgSubtracted=self.config.bgSubtracted)
353  except Exception as e:
354  self.log.warn("Calexp %s not found; skipping it: %s", calExpRef.dataId, e)
355  continue
356 
357  if self.config.doApplySkyCorr:
358  self.applySkyCorr(calExpRef, calExp)
359 
360  calExpList.append(calExp)
361  ccdIdList.append(ccdId)
362  dataIdList.append(calExpRef.dataId)
363 
364  exps = self.run(calExpList, ccdIdList, skyInfo, visitId, dataIdList).exposures
365 
366  if any(exps.values()):
367  dataRefList.append(tempExpRef)
368  else:
369  self.log.warn("Warp %s could not be created", tempExpRef.dataId)
370 
371  if self.config.doWrite:
372  for (warpType, exposure) in exps.items(): # compatible w/ Py3
373  if exposure is not None:
374  self.log.info("Persisting %s" % self.getTempExpDatasetName(warpType))
375  tempExpRef.put(exposure, self.getTempExpDatasetName(warpType))
376 
377  return dataRefList
378 
379  def run(self, calExpList, ccdIdList, skyInfo, visitId=0, dataIdList=None, **kwargs):
380  """Create a Warp from inputs
381 
382  We iterate over the multiple calexps in a single exposure to construct
383  the warp (previously called a coaddTempExp) of that exposure to the
384  supplied tract/patch.
385 
386  Pixels that receive no pixels are set to NAN; this is not correct
387  (violates LSST algorithms group policy), but will be fixed up by
388  interpolating after the coaddition.
389 
390  @param calexpRefList: List of data references for calexps that (may)
391  overlap the patch of interest
392  @param skyInfo: Struct from CoaddBaseTask.getSkyInfo() with geometric
393  information about the patch
394  @param visitId: integer identifier for visit, for the table that will
395  produce the CoaddPsf
396  @return a pipeBase Struct containing:
397  - exposures: a dictionary containing the warps requested:
398  "direct": direct warp if config.makeDirect
399  "psfMatched": PSF-matched warp if config.makePsfMatched
400  """
401  warpTypeList = self.getWarpTypeList()
402 
403  totGoodPix = {warpType: 0 for warpType in warpTypeList}
404  didSetMetadata = {warpType: False for warpType in warpTypeList}
405  coaddTempExps = {warpType: self._prepareEmptyExposure(skyInfo) for warpType in warpTypeList}
406  inputRecorder = {warpType: self.inputRecorder.makeCoaddTempExpRecorder(visitId, len(calExpList))
407  for warpType in warpTypeList}
408 
409  modelPsf = self.config.modelPsf.apply() if self.config.makePsfMatched else None
410  if dataIdList is None:
411  dataIdList = ccdIdList
412 
413  for calExpInd, (calExp, ccdId, dataId) in enumerate(zip(calExpList, ccdIdList, dataIdList)):
414  self.log.info("Processing calexp %d of %d for this Warp: id=%s",
415  calExpInd+1, len(calExpList), dataId)
416 
417  try:
418  warpedAndMatched = self.warpAndPsfMatch.run(calExp, modelPsf=modelPsf,
419  wcs=skyInfo.wcs, maxBBox=skyInfo.bbox,
420  makeDirect=self.config.makeDirect,
421  makePsfMatched=self.config.makePsfMatched)
422  except Exception as e:
423  self.log.warn("WarpAndPsfMatch failed for calexp %s; skipping it: %s", dataId, e)
424  continue
425  try:
426  numGoodPix = {warpType: 0 for warpType in warpTypeList}
427  for warpType in warpTypeList:
428  exposure = warpedAndMatched.getDict()[warpType]
429  if exposure is None:
430  continue
431  coaddTempExp = coaddTempExps[warpType]
432  if didSetMetadata[warpType]:
433  mimg = exposure.getMaskedImage()
434  mimg *= (coaddTempExp.getPhotoCalib().getInstFluxAtZeroMagnitude() /
435  exposure.getPhotoCalib().getInstFluxAtZeroMagnitude())
436  del mimg
437  numGoodPix[warpType] = coaddUtils.copyGoodPixels(
438  coaddTempExp.getMaskedImage(), exposure.getMaskedImage(), self.getBadPixelMask())
439  totGoodPix[warpType] += numGoodPix[warpType]
440  self.log.debug("Calexp %s has %d good pixels in this patch (%.1f%%) for %s",
441  dataId, numGoodPix[warpType],
442  100.0*numGoodPix[warpType]/skyInfo.bbox.getArea(), warpType)
443  if numGoodPix[warpType] > 0 and not didSetMetadata[warpType]:
444  coaddTempExp.setPhotoCalib(exposure.getPhotoCalib())
445  coaddTempExp.setFilter(exposure.getFilter())
446  coaddTempExp.getInfo().setVisitInfo(exposure.getInfo().getVisitInfo())
447  # PSF replaced with CoaddPsf after loop if and only if creating direct warp
448  coaddTempExp.setPsf(exposure.getPsf())
449  didSetMetadata[warpType] = True
450 
451  # Need inputRecorder for CoaddApCorrMap for both direct and PSF-matched
452  inputRecorder[warpType].addCalExp(calExp, ccdId, numGoodPix[warpType])
453 
454  except Exception as e:
455  self.log.warn("Error processing calexp %s; skipping it: %s", dataId, e)
456  continue
457 
458  for warpType in warpTypeList:
459  self.log.info("%sWarp has %d good pixels (%.1f%%)",
460  warpType, totGoodPix[warpType], 100.0*totGoodPix[warpType]/skyInfo.bbox.getArea())
461 
462  if totGoodPix[warpType] > 0 and didSetMetadata[warpType]:
463  inputRecorder[warpType].finish(coaddTempExps[warpType], totGoodPix[warpType])
464  if warpType == "direct":
465  coaddTempExps[warpType].setPsf(
466  CoaddPsf(inputRecorder[warpType].coaddInputs.ccds, skyInfo.wcs,
467  self.config.coaddPsf.makeControl()))
468  else:
469  if not self.config.doWriteEmptyWarps:
470  # No good pixels. Exposure still empty
471  coaddTempExps[warpType] = None
472 
473  result = pipeBase.Struct(exposures=coaddTempExps)
474  return result
475 
476  def getCalibratedExposure(self, dataRef, bgSubtracted):
477  """Return one calibrated Exposure, possibly with an updated SkyWcs.
478 
479  @param[in] dataRef a sensor-level data reference
480  @param[in] bgSubtracted return calexp with background subtracted? If False get the
481  calexp's background background model and add it to the calexp.
482  @return calibrated exposure
483 
484  @raises MissingExposureError If data for the exposure is not available.
485 
486  If config.doApplyUberCal, the exposure will be photometrically
487  calibrated via the `jointcal_photoCalib` dataset and have its SkyWcs
488  updated to the `jointcal_wcs`, otherwise it will be calibrated via the
489  Exposure's own PhotoCalib and have the original SkyWcs.
490  """
491  try:
492  exposure = dataRef.get(self.calexpType, immediate=True)
493  except dafPersist.NoResults as e:
494  raise MissingExposureError('Exposure not found: %s ' % str(e)) from e
495 
496  if not bgSubtracted:
497  background = dataRef.get("calexpBackground", immediate=True)
498  mi = exposure.getMaskedImage()
499  mi += background.getImage()
500  del mi
501 
502  if self.config.doApplyUberCal:
503  if self.config.useMeasMosaic:
504  from lsst.meas.mosaic import applyMosaicResultsExposure
505  # NOTE: this changes exposure in-place, updating its Calib and Wcs.
506  # Save the calibration error, as it gets overwritten with zero.
507  calibrationErr = exposure.getPhotoCalib().getCalibrationErr()
508  try:
509  applyMosaicResultsExposure(dataRef, calexp=exposure)
510  except dafPersist.NoResults as e:
511  raise MissingExposureError('Mosaic calibration not found: %s ' % str(e)) from e
512  photoCalib = afwImage.PhotoCalib(exposure.getPhotoCalib().getCalibrationMean(),
513  calibrationErr,
514  exposure.getBBox())
515  else:
516  photoCalib = dataRef.get("jointcal_photoCalib")
517  skyWcs = dataRef.get("jointcal_wcs")
518  exposure.setWcs(skyWcs)
519  else:
520  photoCalib = exposure.getPhotoCalib()
521 
522  exposure.maskedImage = photoCalib.calibrateImage(exposure.maskedImage,
523  includeScaleUncertainty=self.config.includeCalibVar)
524  exposure.maskedImage /= photoCalib.getCalibrationMean()
525  exposure.setPhotoCalib(photoCalib)
526  # TODO: The images will have a calibration of 1.0 everywhere once RFC-545 is implemented.
527  # exposure.setCalib(afwImage.Calib(1.0))
528  return exposure
529 
530  @staticmethod
531  def _prepareEmptyExposure(skyInfo):
532  """Produce an empty exposure for a given patch"""
533  exp = afwImage.ExposureF(skyInfo.bbox, skyInfo.wcs)
534  exp.getMaskedImage().set(numpy.nan, afwImage.Mask
535  .getPlaneBitMask("NO_DATA"), numpy.inf)
536  return exp
537 
538  def getWarpTypeList(self):
539  """Return list of requested warp types per the config.
540  """
541  warpTypeList = []
542  if self.config.makeDirect:
543  warpTypeList.append("direct")
544  if self.config.makePsfMatched:
545  warpTypeList.append("psfMatched")
546  return warpTypeList
547 
548  def applySkyCorr(self, dataRef, calexp):
549  """Apply correction to the sky background level
550 
551  Sky corrections can be generated with the 'skyCorrection.py'
552  executable in pipe_drivers. Because the sky model used by that
553  code extends over the entire focal plane, this can produce
554  better sky subtraction.
555 
556  The calexp is updated in-place.
557 
558  Parameters
559  ----------
560  dataRef : `lsst.daf.persistence.ButlerDataRef`
561  Data reference for calexp.
562  calexp : `lsst.afw.image.Exposure` or `lsst.afw.image.MaskedImage`
563  Calibrated exposure.
564  """
565  bg = dataRef.get("skyCorr")
566  if isinstance(calexp, afwImage.Exposure):
567  calexp = calexp.getMaskedImage()
568  calexp -= bg.getImage()
569 
570 
571 class MakeWarpConfig(pipeBase.PipelineTaskConfig, MakeCoaddTempExpConfig):
572  calExpList = pipeBase.InputDatasetField(
573  doc="Input exposures to be resampled and optionally PSF-matched onto a SkyMap projection/patch",
574  name="calexp",
575  storageClass="ExposureF",
576  dimensions=("instrument", "visit", "detector")
577  )
578  backgroundList = pipeBase.InputDatasetField(
579  doc="Input backgrounds to be added back into the calexp if bgSubtracted=False",
580  name="calexpBackground",
581  storageClass="Background",
582  dimensions=("instrument", "visit", "detector")
583  )
584  skyCorrList = pipeBase.InputDatasetField(
585  doc="Input Sky Correction to be subtracted from the calexp if doApplySkyCorr=True",
586  name="skyCorr",
587  storageClass="Background",
588  dimensions=("instrument", "visit", "detector")
589  )
590  skyMap = pipeBase.InputDatasetField(
591  doc="Input definition of geometry/bbox and projection/wcs for warped exposures",
592  nameTemplate="{coaddName}Coadd_skyMap",
593  storageClass="SkyMap",
594  dimensions=("skymap",),
595  scalar=True
596  )
597  direct = pipeBase.OutputDatasetField(
598  doc=("Output direct warped exposure (previously called CoaddTempExp), produced by resampling ",
599  "calexps onto the skyMap patch geometry."),
600  nameTemplate="{coaddName}Coadd_directWarp",
601  storageClass="ExposureF",
602  dimensions=("tract", "patch", "skymap", "visit", "instrument"),
603  scalar=True
604  )
605  psfMatched = pipeBase.OutputDatasetField(
606  doc=("Output PSF-Matched warped exposure (previously called CoaddTempExp), produced by resampling ",
607  "calexps onto the skyMap patch geometry and PSF-matching to a model PSF."),
608  nameTemplate="{coaddName}Coadd_psfMatchedWarp",
609  storageClass="ExposureF",
610  dimensions=("tract", "patch", "skymap", "visit", "instrument"),
611  scalar=True
612  )
613 
614  def setDefaults(self):
615  super().setDefaults()
616  self.formatTemplateNames({"coaddName": "deep"})
617  self.quantum.dimensions = ("tract", "patch", "skymap", "visit")
618 
619  def validate(self):
620  super().validate()
621  # TODO: Remove this constraint after DM-17062
622  if self.doApplyUberCal:
623  raise RuntimeError("Gen3 MakeWarpTask cannot apply meas_mosaic or jointcal results."
624  "Please set doApplyUbercal=False.")
625 
626 
627 class MakeWarpTask(MakeCoaddTempExpTask, pipeBase.PipelineTask):
628  """Warp and optionally PSF-Match calexps onto an a common projection
629 
630  First Draft of a Gen3 compatible MakeWarpTask which
631  currently does not handle doApplyUberCal=True.
632  """
633  ConfigClass = MakeWarpConfig
634  _DefaultName = "makeWarp"
635 
636  @classmethod
637  def getInputDatasetTypes(cls, config):
638  """Return input dataset type descriptors
639 
640  Remove input dataset types not used by the Task
641  """
642  inputTypeDict = super().getInputDatasetTypes(config)
643  if config.bgSubtracted:
644  inputTypeDict.pop("backgroundList", None)
645  if not config.doApplySkyCorr:
646  inputTypeDict.pop("skyCorrList", None)
647  return inputTypeDict
648 
649  @classmethod
650  def getOutputDatasetTypes(cls, config):
651  """Return output dataset type descriptors
652 
653  Remove output dataset types not produced by the Task
654  """
655  outputTypeDict = super().getOutputDatasetTypes(config)
656  if not config.makeDirect:
657  outputTypeDict.pop("direct", None)
658  if not config.makePsfMatched:
659  outputTypeDict.pop("psfMatched", None)
660  return outputTypeDict
661 
662  def adaptArgsAndRun(self, inputData, inputDataIds, outputDataIds, butler):
663  """Construct warps for requested warp type for single epoch
664 
665  PipelineTask (Gen3) entry point to warp and optionally PSF-match
666  calexps. This method is analogous to `runDataRef`, it prepares all
667  the data products to be passed to `run`.
668  Return a Struct with only requested warpTypes controlled by the configs
669  makePsfMatched and makeDirect.
670 
671  Parameters
672  ----------
673  inputData : `dict`
674  Keys are the names of the configs describing input dataset types.
675  Values are input Python-domain data objects (or lists of objects)
676  retrieved from data butler.
677  inputDataIds : `dict`
678  Keys are the names of the configs describing input dataset types.
679  Values are DataIds (or lists of DataIds) that task consumes for
680  corresponding dataset type.
681  outputDataIds : `dict`
682  Keys are the names of the configs describing input dataset types.
683  Values are DataIds (or lists of DataIds) that task is to produce
684  for corresponding dataset type.
685  butler : `lsst.daf.butler.Butler`
686  Gen3 Butler object for fetching additional data products before
687  running the Task
688 
689  Returns
690  -------
691  result : `lsst.pipe.base.Struct`
692  Result struct with components:
693 
694  - ``direct`` : (optional) direct Warp Exposure
695  (``lsst.afw.image.Exposure``)
696  - ``psfMatched``: (optional) PSF-Matched Warp Exposure
697  (``lsst.afw.image.Exposure``)
698  """
699  # Construct skyInfo expected by `run`
700  skyMap = inputData["skyMap"]
701  outputDataId = next(iter(outputDataIds.values()))
702  inputData['skyInfo'] = makeSkyInfo(skyMap,
703  tractId=outputDataId['tract'],
704  patchId=outputDataId['patch'])
705 
706  # Construct list of DataIds expected by `run`
707  dataIdList = inputDataIds['calExpList']
708  inputData['dataIdList'] = dataIdList
709 
710  # Construct list of ccdExposureIds expected by `run`
711  inputData['ccdIdList'] = [butler.registry.packDataId("visit_detector", dataId)
712  for dataId in dataIdList]
713 
714  # Extract integer visitId requested by `run`
715  visits = [dataId['visit'] for dataId in dataIdList]
716  assert(all(visits[0] == visit for visit in visits))
717  inputData["visitId"] = visits[0]
718 
719  self.prepareCalibratedExposures(**inputData)
720  results = self.run(**inputData)
721  return pipeBase.Struct(**results.exposures)
722 
723  def prepareCalibratedExposures(self, calExpList, backgroundList=None, skyCorrList=None, **kwargs):
724  """Calibrate and add backgrounds to input calExpList in place
725 
726  TODO DM-17062: apply jointcal/meas_mosaic here
727 
728  Parameters
729  ----------
730  calExpList : `list` of `lsst.afw.image.Exposure`
731  Sequence of calexps to be modified in place
732  backgroundList : `list` of `lsst.afw.math.backgroundList`
733  Sequence of backgrounds to be added back in if bgSubtracted=False
734  skyCorrList : `list` of `lsst.afw.math.backgroundList`
735  Sequence of background corrections to be subtracted if doApplySkyCorr=True
736  """
737  backgroundList = len(calExpList)*[None] if backgroundList is None else backgroundList
738  skyCorrList = len(calExpList)*[None] if skyCorrList is None else skyCorrList
739  for calexp, background, skyCorr in zip(calExpList, backgroundList, skyCorrList):
740  mi = calexp.maskedImage
741  if not self.config.bgSubtracted:
742  mi += background.getImage()
743  if self.config.doApplySkyCorr:
744  mi -= skyCorr.getImage()
def getCoaddDatasetName(self, warpType="direct")
Definition: coaddBase.py:149
def getGroupDataRef(butler, datasetType, groupTuple, keys)
Definition: coaddHelpers.py:99
Base class for coaddition.
Definition: coaddBase.py:100
The photometric calibration of an exposure.
Definition: PhotoCalib.h:116
def prepareCalibratedExposures(self, calExpList, backgroundList=None, skyCorrList=None, kwargs)
A class to contain the data, WCS, and other information needed to describe an image of the sky...
Definition: Exposure.h:72
def makeSkyInfo(skyMap, tractId, patchId)
Definition: coaddBase.py:249
Fit spatial kernel using approximate fluxes for candidates, and solving a linear system of equations...
daf::base::PropertySet * set
Definition: fits.cc:884
Warp and optionally PSF-Match calexps onto an a common projection.
bool any(CoordinateExpr< N > const &expr) noexcept
Return true if any elements are true.
Definition: Log.h:691
def adaptArgsAndRun(self, inputData, inputDataIds, outputDataIds, butler)
def getSkyInfo(self, patchRef)
Use getSkyinfo to return the skyMap, tract and patch information, wcs and the outer bbox of the patch...
Definition: coaddBase.py:133
def getTempExpDatasetName(self, warpType="direct")
Definition: coaddBase.py:164
bool all(CoordinateExpr< N > const &expr) noexcept
Return true if all elements are true.
def run(self, calExpList, ccdIdList, skyInfo, visitId=0, dataIdList=None, kwargs)
def getBadPixelMask(self)
Convenience method to provide the bitmask from the mask plane names.
Definition: coaddBase.py:199
Represent a 2-dimensional array of bitmask pixels.
Definition: Mask.h:78
CoaddPsf is the Psf derived to be used for non-PSF-matched Coadd images.
Definition: CoaddPsf.h:58
def getCalibratedExposure(self, dataRef, bgSubtracted)
def selectExposures(self, patchRef, skyInfo=None, selectDataList=[])
Select exposures to coadd.
Definition: coaddBase.py:113
def runDataRef(self, patchRef, selectDataList=[])
Produce <coaddName>Coadd_<warpType>Warp images by warping and optionally PSF-matching.
Backwards-compatibility support for depersisting the old Calib (FluxMag0/FluxMag0Err) objects...
int copyGoodPixels(lsst::afw::image::MaskedImage< ImagePixelT, lsst::afw::image::MaskPixel, lsst::afw::image::VariancePixel > &destImage, lsst::afw::image::MaskedImage< ImagePixelT, lsst::afw::image::MaskPixel, lsst::afw::image::VariancePixel > const &srcImage, lsst::afw::image::MaskPixel const badPixelMask)
copy good pixels from one masked image to another
def groupPatchExposures(patchDataRef, calexpDataRefList, coaddDatasetType="deepCoadd", tempExpDatasetType="deepCoadd_directWarp")
Definition: coaddHelpers.py:60