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