Coverage for tests/assembleCoaddTestUtils.py : 21%

Hot-keys on this page
r m x p toggle line displays
j k next/prev highlighted chunk
0 (zero) top of page
1 (one) first highlighted chunk
1# This file is part of pipe_tasks.
2#
3# LSST Data Management System
4# This product includes software developed by the
5# LSST Project (http://www.lsst.org/).
6# See COPYRIGHT file at the top of the source tree.
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 <https://www.lsstcorp.org/LegalNotices/>.
21#
22"""Set up simulated test data and simplified APIs for AssembleCoaddTask
23and its derived classes.
25This is not intended to test accessing data with the Butler and instead uses
26mock Butler data references to pass in the simulated data.
27"""
28from astropy.time import Time
29from astropy import units as u
30from astropy.coordinates import SkyCoord, EarthLocation, Angle
31import numpy as np
33from lsst.afw.cameraGeom.testUtils import DetectorWrapper
34import lsst.afw.geom as afwGeom
35import lsst.afw.image as afwImage
36import lsst.afw.image.utils as afwImageUtils
37import lsst.daf.butler
38import lsst.daf.persistence.dataId
39import lsst.geom as geom
40from lsst.geom import arcseconds, degrees
41from lsst.meas.algorithms.testUtils import plantSources
42from lsst.obs.base import MakeRawVisitInfoViaObsInfo
43import lsst.pipe.base as pipeBase
44from lsst.pipe.tasks.coaddInputRecorder import CoaddInputRecorderTask, CoaddInputRecorderConfig
46from astro_metadata_translator import makeObservationInfo
48__all__ = ["MockWarpReference", "makeMockSkyInfo", "MockCoaddTestData",
49 "MockGen2WarpReference"]
52class MockGen2WarpReference:
53 """Very simple object that looks like a Gen2 data reference to a warped
54 exposure.
56 Parameters
57 ----------
58 exposure : `lsst.afw.image.Exposure`
59 The exposure to be retrieved by the data reference.
60 exposurePsfMatched : `lsst.afw.image.Exposure`, optional
61 The exposure to be retrieved by the data reference, with a degraded PSF.
62 coaddName : `str`
63 The type of coadd being produced. Typically 'deep'.
64 patch : `str`
65 Unique identifier for a subdivision of a tract.
66 In Gen 2 the patch identifier consists
67 of two integers separated by a comma.
68 tract : `int`
69 Unique identifier for a tract of a skyMap
70 visit : `int`
71 Unique identifier for an observation,
72 potentially consisting of multiple ccds.
73 """
74 datasetTypes = None
75 "List of the names of exposures that can be retrieved."
76 metadataTypes = None
77 "List of the names of metadata objects that can be retrieved."
78 dataLookup = None
79 "Stores the data and metadata that can be retrieved."
81 def __init__(self, exposure, exposurePsfMatched=None, coaddName='deep',
82 patch="2,3", tract=0, visit=100):
83 self.coaddName = coaddName
84 self.tract = tract
85 self.patch = patch
86 self.visit = visit
87 visitInfo = exposure.getInfo().getVisitInfo()
89 self.datasetTypes = (f"{coaddName}Coadd_directWarp", f"{coaddName}Coadd_psfMatchedWarp",
90 f"{coaddName}Coadd_directWarp_sub", f"{coaddName}Coadd_psfMatchedWarp_sub")
92 self.metadataTypes = (f"{coaddName}Coadd_directWarp_visitInfo",
93 f"{coaddName}Coadd_psfMatchedWarp_visitInfo")
95 self.dataLookup = {f"{coaddName}Coadd_directWarp": exposure,
96 f"{coaddName}Coadd_psfMatchedWarp": exposurePsfMatched,
97 f"{coaddName}Coadd_directWarp_sub": exposure,
98 f"{coaddName}Coadd_psfMatchedWarp_sub": exposurePsfMatched,
99 f"{coaddName}Coadd_directWarp_visitInfo": visitInfo,
100 f"{coaddName}Coadd_psfMatchedWarp_visitInfo": visitInfo}
102 def get(self, datasetType, bbox=None, **kwargs):
103 """Retrieve the specified dataset using the API of the Gen 2 Butler.
105 Parameters
106 ----------
107 datasetType : `str`
108 Name of the type of exposure to retrieve.
109 bbox : `lsst.geom.box.Box2I`, optional
110 If supplied and the `datasetType ends in "_sub",
111 then retrieve only a subregion of the exposure.
112 **kwargs
113 Additional keyword arguments such as `immediate=True` that would
114 control internal butler behavior.
116 Returns
117 -------
118 `lsst.afw.image.Exposure` or `lsst.afw.image.VisitInfo`
119 Either the exposure or its metadata, depending on the datasetType.
121 Raises
122 ------
123 KeyError
124 If a bounding box is specified incorrectly for the datasetType.
125 ValueError
126 If an unknown datasetType is supplied.
127 """
128 if "_sub" in datasetType:
129 if bbox is None:
130 raise KeyError(f"A bbox must be supplied for dataset {datasetType}")
131 else:
132 if bbox is not None:
133 raise KeyError(f"A bbox cannot be supplied for dataset {datasetType}")
134 if datasetType in self.datasetTypes:
135 exp = self.dataLookup[datasetType].clone()
136 if "_sub" in datasetType:
137 return exp[bbox]
138 else:
139 return exp
140 elif datasetType in self.metadataTypes:
141 return self.dataLookup[datasetType]
142 else:
143 raise ValueError(f"Unknown datasetType {datasetType}. Must be one of {self.datasetTypes}.")
145 @property
146 def dataId(self):
147 """Generate a valid data identifier.
149 Returns
150 -------
151 dataId : `lsst.daf.persistence.dataId`
152 Data identifier dict for patch.
153 """
154 return lsst.daf.persistence.dataId.DataId(tract=self.tract, patch=self.patch, visit=self.visit)
156 def datasetExists(self, tempExpName):
157 """Mimic the Gen2 Butler API for determining whether a dataset exists.
158 """
159 if tempExpName in self.datasetTypes:
160 return True
161 else:
162 return False
165class MockWarpReference(lsst.daf.butler.DeferredDatasetHandle):
166 """Very simple object that looks like a Gen 3 data reference to a warped
167 exposure.
169 Parameters
170 ----------
171 exposure : `lsst.afw.image.Exposure`
172 The exposure to be retrieved by the data reference.
173 coaddName : `str`
174 The type of coadd being produced. Typically 'deep'.
175 patch : `int`
176 Unique identifier for a subdivision of a tract.
177 tract : `int`
178 Unique identifier for a tract of a skyMap
179 visit : `int`
180 Unique identifier for an observation,
181 potentially consisting of multiple ccds.
182 """
183 def __init__(self, exposure, coaddName='deep', patch=42, tract=0, visit=100):
184 self.coaddName = coaddName
185 self.exposure = exposure
186 self.tract = tract
187 self.patch = patch
188 self.visit = visit
190 def get(self, bbox=None, component=None, parameters=None, **kwargs):
191 """Retrieve the specified dataset using the API of the Gen 3 Butler.
193 Parameters
194 ----------
195 bbox : `lsst.geom.box.Box2I`, optional
196 If supplied, retrieve only a subregion of the exposure.
197 component : `str`, optional
198 If supplied, return the named metadata of the exposure.
199 parameters : `dict`, optional
200 If supplied, use the parameters to modify the exposure,
201 typically by taking a subset.
202 **kwargs
203 Additional keyword arguments such as `immediate=True` that would
204 control internal butler behavior.
206 Returns
207 -------
208 `lsst.afw.image.Exposure` or `lsst.afw.image.VisitInfo`
209 or `lsst.meas.algorithms.SingleGaussianPsf`
210 Either the exposure or its metadata, depending on the datasetType.
211 """
212 if component == 'psf':
213 return self.exposure.getPsf()
214 elif component == 'visitInfo':
215 return self.exposure.getInfo().getVisitInfo()
216 if parameters is not None:
217 if "bbox" in parameters:
218 bbox = parameters["bbox"]
219 exp = self.exposure.clone()
220 if bbox is not None:
221 return exp[bbox]
222 else:
223 return exp
225 @property
226 def dataId(self):
227 """Generate a valid data identifier.
229 Returns
230 -------
231 dataId : `lsst.daf.persistence.dataId`
232 Data identifier dict for the patch.
233 """
234 return lsst.daf.persistence.dataId.DataId(tract=self.tract, patch=self.patch, visit=self.visit)
236 def datasetExists(self, tempExpName):
237 """Raise a more informative error if this Gen 2 method is called."""
238 raise NotImplementedError("Gen3 butler data references don't support `datasetExists`")
241def makeMockSkyInfo(bbox, wcs, patch):
242 """Construct a `Struct` containing the geometry of the patch to be coadded.
244 Parameters
245 ----------
246 bbox : `lsst.geom.Box`
247 Bounding box of the patch to be coadded.
248 wcs : `lsst.afw.geom.SkyWcs`
249 Coordinate system definition (wcs) for the exposure.
251 Returns
252 -------
253 skyInfo : `lsst.pipe.base.Struct`
254 Patch geometry information.
255 """
256 def getIndex():
257 return patch
258 patchInfo = pipeBase.Struct(getIndex=getIndex)
259 skyInfo = pipeBase.Struct(bbox=bbox, wcs=wcs, patchInfo=patchInfo)
260 return skyInfo
263class MockCoaddTestData:
264 """Generate repeatable simulated exposures with consistent metadata that
265 are realistic enough to test the image coaddition algorithms.
267 Notes
268 -----
269 The simple GaussianPsf used by lsst.meas.algorithms.testUtils.plantSources
270 will always return an average position of (0, 0).
271 The bounding box of the exposures MUST include (0, 0), or else the PSF will
272 not be valid and `AssembleCoaddTask` will fail with the error
273 'Could not find a valid average position for CoaddPsf'.
275 Parameters
276 ----------
277 shape : `lsst.geom.Extent2I`, optional
278 Size of the bounding box of the exposures to be simulated, in pixels.
279 offset : `lsst.geom.Point2I`, optional
280 Pixel coordinate of the lower left corner of the bounding box.
281 backgroundLevel : `float`, optional
282 Background value added to all pixels in the simulated images.
283 seed : `int`, optional
284 Seed value to initialize the random number generator.
285 nSrc : `int`, optional
286 Number of sources to simulate.
287 fluxRange : `float`, optional
288 Range in flux amplitude of the simulated sources.
289 noiseLevel : `float`, optional
290 Standard deviation of the noise to add to each pixel.
291 sourceSigma : `float`, optional
292 Average amplitude of the simulated sources,
293 relative to ``noiseLevel``
294 minPsfSize : `float`, optional
295 The smallest PSF width (sigma) to use, in pixels.
296 maxPsfSize : `float`, optional
297 The largest PSF width (sigma) to use, in pixels.
298 pixelScale : `lsst.geom.Angle`, optional
299 The plate scale of the simulated images.
300 ra : `lsst.geom.Angle`, optional
301 Right Ascension of the boresight of the camera for the observation.
302 dec : `lsst.geom.Angle`, optional
303 Declination of the boresight of the camera for the observation.
304 ccd : `int`, optional
305 CCD number to put in the metadata of the exposure.
306 patch : `int`, optional
307 Unique identifier for a subdivision of a tract.
308 patchGen2 : `str`, optional
309 Unique identifier for a subdivision of a tract.
310 In Gen 2 the patch identifier consists
311 of two integers separated by a comma.
312 tract : `int`, optional
313 Unique identifier for a tract of a skyMap.
315 Raises
316 ------
317 ValueError
318 If the bounding box does not contain the pixel coordinate (0, 0).
319 This is due to `GaussianPsf` that is used by `lsst.meas.algorithms.testUtils.plantSources`
320 lacking the option to specify the pixel origin.
321 """
322 rotAngle = 0.*degrees
323 "Rotation of the pixel grid on the sky, East from North (`lsst.geom.Angle`)."
324 filterInfo = None
325 """The filter definition, usually set in the current instruments' obs package.
326 For these tests, a simple filter is defined without using an obs package (`lsst.afw.image.Filter`).
327 """
328 rngData = None
329 """Pre-initialized random number generator for constructing the test images
330 repeatably (`numpy.random.Generator`).
331 """
332 rngMods = None
333 """Pre-initialized random number generator for applying modifications to
334 the test images for only some test cases (`numpy.random.Generator`).
335 """
336 kernelSize = None
337 "Width of the kernel used for simulating sources, in pixels."
338 exposures = {}
339 "The simulated test data, with variable PSF sizes (`dict` of `lsst.afw.image.Exposure`)"
340 matchedExposures = {}
341 """The simulated exposures, all with PSF width set to `maxPsfSize`
342 (`dict` of `lsst.afw.image.Exposure`).
343 """
344 photoCalib = afwImage.makePhotoCalibFromCalibZeroPoint(27, 10)
345 "The photometric zero point to use for converting counts to flux units (`lsst.afw.image.PhotoCalib`)."
346 badMaskPlanes = ["NO_DATA", "BAD"]
347 "Mask planes that, if set, the associated pixel should not be included in the coaddTempExp."
348 detector = None
349 "Properties of the CCD for the exposure (`lsst.afw.cameraGeom.Detector`)."
351 def __init__(self, shape=geom.Extent2I(201, 301), offset=geom.Point2I(-123, -45),
352 backgroundLevel=314.592, seed=42, nSrc=37,
353 fluxRange=2., noiseLevel=5, sourceSigma=200.,
354 minPsfSize=1.5, maxPsfSize=3.,
355 pixelScale=0.2*arcseconds, ra=209.*degrees, dec=-20.25*degrees,
356 ccd=37, patch=42, patchGen2="2,3", tract=0):
357 self.ra = ra
358 self.dec = dec
359 self.pixelScale = pixelScale
360 self.patch = patch
361 self.patchGen2 = patchGen2
362 self.tract = tract
363 self.filterInfo = self.defineFilter()
364 self.rngData = np.random.default_rng(seed)
365 self.rngMods = np.random.default_rng(seed + 1)
366 self.bbox = geom.Box2I(offset, shape)
367 if not self.bbox.contains(0, 0):
368 raise ValueError(f"The bounding box must contain the coordinate (0, 0). {repr(self.bbox)}")
369 self.wcs = self.makeDummyWcs()
371 # Set up properties of the simulations
372 nSigmaForKernel = 5
373 self.kernelSize = (int(maxPsfSize*nSigmaForKernel + 0.5)//2)*2 + 1 # make sure it is odd
375 bufferSize = self.kernelSize//2
376 x0, y0 = self.bbox.getBegin()
377 xSize, ySize = self.bbox.getDimensions()
378 # Set the pixel coordinates and fluxes of the simulated sources.
379 self.xLoc = self.rngData.random(nSrc)*(xSize - 2*bufferSize) + bufferSize + x0
380 self.yLoc = self.rngData.random(nSrc)*(ySize - 2*bufferSize) + bufferSize + y0
381 self.flux = (self.rngData.random(nSrc)*(fluxRange - 1.) + 1.)*sourceSigma*noiseLevel
383 self.backgroundLevel = backgroundLevel
384 self.noiseLevel = noiseLevel
385 self.minPsfSize = minPsfSize
386 self.maxPsfSize = maxPsfSize
387 self.detector = DetectorWrapper(name=f"detector {ccd}", id=ccd).detector
389 def setDummyCoaddInputs(self, exposure, expId):
390 """Generate an `ExposureCatalog` as though the exposures had been
391 processed using `warpAndPsfMatch`.
393 Parameters
394 ----------
395 exposure : `lsst.afw.image.Exposure`
396 The exposure to construct a `CoaddInputs` `ExposureCatalog` for.
397 expId : `int`
398 A unique identifier for the visit.
399 """
400 badPixelMask = afwImage.Mask.getPlaneBitMask(self.badMaskPlanes)
401 nGoodPix = np.sum(exposure.getMask().getArray() & badPixelMask == 0)
403 config = CoaddInputRecorderConfig()
404 inputRecorder = CoaddInputRecorderTask(config=config, name="inputRecorder")
405 tempExpInputRecorder = inputRecorder.makeCoaddTempExpRecorder(expId, num=1)
406 tempExpInputRecorder.addCalExp(exposure, expId, nGoodPix)
407 tempExpInputRecorder.finish(exposure, nGoodPix=nGoodPix)
409 def makeCoaddTempExp(self, rawExposure, visitInfo, expId):
410 """Add the metadata required by `AssembleCoaddTask` to an exposure.
412 Parameters
413 ----------
414 rawExposure : `lsst.afw.image.Exposure`
415 The simulated exposure.
416 visitInfo : `lsst.afw.image.VisitInfo`
417 VisitInfo containing metadata for the exposure.
418 expId : `int`
419 A unique identifier for the visit.
421 Returns
422 -------
423 tempExp : `lsst.afw.image.Exposure`
424 The exposure, with all of the metadata needed for coaddition.
425 """
426 tempExp = rawExposure.clone()
427 tempExp.setWcs(self.wcs)
429 tempExp.setFilter(self.filterInfo)
430 tempExp.setPhotoCalib(self.photoCalib)
431 tempExp.getInfo().setVisitInfo(visitInfo)
432 tempExp.getInfo().setDetector(self.detector)
433 self.setDummyCoaddInputs(tempExp, expId)
434 return tempExp
436 def makeDummyWcs(self, rotAngle=None, pixelScale=None, crval=None, flipX=True):
437 """Make a World Coordinate System object for testing.
439 Parameters
440 ----------
441 rotAngle : `lsst.geom.Angle`
442 Rotation of the CD matrix, East from North
443 pixelScale : `lsst.geom.Angle`
444 Pixel scale of the projection.
445 crval : `lsst.afw.geom.SpherePoint`
446 Coordinates of the reference pixel of the wcs.
447 flipX : `bool`, optional
448 Flip the direction of increasing Right Ascension.
450 Returns
451 -------
452 wcs : `lsst.afw.geom.skyWcs.SkyWcs`
453 A wcs that matches the inputs.
454 """
455 if rotAngle is None:
456 rotAngle = self.rotAngle
457 if pixelScale is None:
458 pixelScale = self.pixelScale
459 if crval is None:
460 crval = geom.SpherePoint(self.ra, self.dec)
461 crpix = geom.Box2D(self.bbox).getCenter()
462 cdMatrix = afwGeom.makeCdMatrix(scale=pixelScale, orientation=rotAngle, flipX=flipX)
463 wcs = afwGeom.makeSkyWcs(crpix=crpix, crval=crval, cdMatrix=cdMatrix)
464 return wcs
466 def makeDummyVisitInfo(self, exposureId, randomizeTime=False):
467 """Make a self-consistent visitInfo object for testing.
469 Parameters
470 ----------
471 exposureId : `int`, optional
472 Unique integer identifier for this observation.
473 randomizeTime : `bool`, optional
474 Add a random offset within a 6 hour window to the observation time.
476 Returns
477 -------
478 visitInfo : `lsst.afw.image.VisitInfo`
479 VisitInfo for the exposure.
480 """
481 lsstLat = -30.244639*u.degree
482 lsstLon = -70.749417*u.degree
483 lsstAlt = 2663.*u.m
484 lsstTemperature = 20.*u.Celsius
485 lsstHumidity = 40. # in percent
486 lsstPressure = 73892.*u.pascal
487 loc = EarthLocation(lat=lsstLat,
488 lon=lsstLon,
489 height=lsstAlt)
491 time = Time(2000.0, format="jyear", scale="tt")
492 if randomizeTime:
493 # Pick a random time within a 6 hour window
494 time += 6*u.hour*(self.rngMods.random() - 0.5)
495 radec = SkyCoord(dec=self.dec.asDegrees(), ra=self.ra.asDegrees(),
496 unit='deg', obstime=time, frame='icrs', location=loc)
497 airmass = float(1.0/np.sin(radec.altaz.alt))
498 obsInfo = makeObservationInfo(location=loc,
499 detector_exposure_id=exposureId,
500 datetime_begin=time,
501 datetime_end=time,
502 boresight_airmass=airmass,
503 boresight_rotation_angle=Angle(0.*u.degree),
504 boresight_rotation_coord='sky',
505 temperature=lsstTemperature,
506 pressure=lsstPressure,
507 relative_humidity=lsstHumidity,
508 tracking_radec=radec,
509 altaz_begin=radec.altaz,
510 observation_type='science',
511 )
512 visitInfo = MakeRawVisitInfoViaObsInfo.observationInfo2visitInfo(obsInfo)
513 return visitInfo
515 def defineFilter(self, lambdaEff=476.31, lambdaMin=405., lambdaMax=552., filterName="gTest"):
516 """Construct a `Filter` with sufficient information to calculate DCR.
518 Parameters
519 ----------
520 lambdaEff : `float`, optional
521 The effective wavelength of the filter, defaults to LSST g-band value.
522 lambdaMin : float, optional
523 The minimum wavelength of the filter with greater than 1% transmission,
524 defaults to LSST g-band value.
525 lambdaMax : `float`, optional
526 The maximum wavelength of the filter with greater than 1% transmission,
527 defaults to LSST g-band value.
528 filterName : `str`, optional
529 The simplified name of the filter.
531 Returns
532 -------
533 filterInfo : `lsst.afw.image.Filter`
534 The filter definition.
535 """
536 afwImageUtils.defineFilter(filterName, lambdaEff,
537 lambdaMin=lambdaMin, lambdaMax=lambdaMax)
538 return afwImage.Filter(filterName)
540 def makeTestImage(self, expId, noiseLevel=None, psfSize=None, backgroundLevel=None,
541 detectionSigma=5.):
542 """Make a reproduceable PSF-convolved masked image for testing.
544 Parameters
545 ----------
546 expId : `int`
547 A unique identifier to use to refer to the visit.
548 noiseLevel : `float`, optional
549 Standard deviation of the noise to add to each pixel.
550 psfSize : `float`, optional
551 Width of the PSF of the simulated sources, in pixels.
552 backgroundLevel : `float`, optional
553 Background value added to all pixels in the simulated images.
554 detectionSigma : `float`, optional
555 Threshold amplitude of the image to set the "DETECTED" mask.
556 """
557 if backgroundLevel is None:
558 backgroundLevel = self.backgroundLevel
559 if noiseLevel is None:
560 noiseLevel = 5.
561 visitInfo = self.makeDummyVisitInfo(expId, randomizeTime=True)
563 if psfSize is None:
564 psfSize = self.rngMods.random()*(self.maxPsfSize - self.minPsfSize) + self.minPsfSize
565 nSrc = len(self.flux)
566 sigmas = [psfSize for src in range(nSrc)]
567 sigmasPsfMatched = [self.maxPsfSize for src in range(nSrc)]
568 coordList = list(zip(self.xLoc, self.yLoc, self.flux, sigmas))
569 coordListPsfMatched = list(zip(self.xLoc, self.yLoc, self.flux, sigmasPsfMatched))
570 xSize, ySize = self.bbox.getDimensions()
571 model = plantSources(self.bbox, self.kernelSize, self.backgroundLevel,
572 coordList, addPoissonNoise=False)
573 modelPsfMatched = plantSources(self.bbox, self.kernelSize, self.backgroundLevel,
574 coordListPsfMatched, addPoissonNoise=False)
575 model.variance.array = np.abs(model.image.array) + noiseLevel
576 modelPsfMatched.variance.array = np.abs(modelPsfMatched.image.array) + noiseLevel
577 noise = self.rngData.random((ySize, xSize))*noiseLevel
578 noise -= np.median(noise)
579 model.image.array += noise
580 modelPsfMatched.image.array += noise
581 detectedMask = afwImage.Mask.getPlaneBitMask("DETECTED")
582 detectionThreshold = self.backgroundLevel + detectionSigma*noiseLevel
583 model.mask.array[model.image.array > detectionThreshold] += detectedMask
585 exposure = self.makeCoaddTempExp(model, visitInfo, expId)
586 matchedExposure = self.makeCoaddTempExp(modelPsfMatched, visitInfo, expId)
587 return exposure, matchedExposure
589 @staticmethod
590 def makeGen2DataRefList(exposures, matchedExposures, tract=0, patch="2,3", coaddName="deep"):
591 """Make data references from the simulated exposures that can be
592 retrieved using the Gen 2 Butler API.
594 Parameters
595 ----------
596 tract : `int`
597 Unique identifier for a tract of a skyMap.
598 patch : `str`
599 Unique identifier for a subdivision of a tract.
600 coaddName : `str`
601 The type of coadd being produced. Typically 'deep'.
603 Returns
604 -------
605 dataRefList : `list` of `MockGen2WarpReference`
606 The data references.
608 """
609 dataRefList = []
610 for expId in exposures:
611 exposure = exposures[expId]
612 exposurePsfMatched = matchedExposures[expId]
613 dataRef = MockGen2WarpReference(exposure, exposurePsfMatched=exposurePsfMatched,
614 coaddName=coaddName, tract=tract, patch=patch, visit=expId)
615 dataRefList.append(dataRef)
616 return dataRefList
618 @staticmethod
619 def makeDataRefList(exposures, matchedExposures, warpType, tract=0, patch=42, coaddName="deep"):
620 """Make data references from the simulated exposures that can be
621 retrieved using the Gen 3 Butler API.
623 Parameters
624 ----------
625 warpType : `str`
626 Either 'direct' or 'psfMatched'.
627 tract : `int`, optional
628 Unique identifier for a tract of a skyMap.
629 patch : `int`, optional
630 Unique identifier for a subdivision of a tract.
631 coaddName : `str`, optional
632 The type of coadd being produced. Typically 'deep'.
634 Returns
635 -------
636 dataRefList : `list` of `MockWarpReference`
637 The data references.
639 Raises
640 ------
641 ValueError
642 If an unknown `warpType` is supplied.
643 """
644 dataRefList = []
645 for expId in exposures:
646 if warpType == 'direct':
647 exposure = exposures[expId]
648 elif warpType == 'psfMatched':
649 exposure = matchedExposures[expId]
650 else:
651 raise ValueError("warpType must be one of 'direct' or 'psfMatched'")
652 dataRef = MockWarpReference(exposure, coaddName=coaddName,
653 tract=tract, patch=patch, visit=expId)
654 dataRefList.append(dataRef)
655 return dataRefList