Hide keyboard shortcuts

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. 

24 

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 

32 

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 

45 

46from astro_metadata_translator import makeObservationInfo 

47 

48__all__ = ["MockWarpReference", "makeMockSkyInfo", "MockCoaddTestData", 

49 "MockGen2WarpReference"] 

50 

51 

52class MockGen2WarpReference: 

53 """Very simple object that looks like a Gen2 data reference to a warped 

54 exposure. 

55 

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." 

80 

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() 

88 

89 self.datasetTypes = (f"{coaddName}Coadd_directWarp", f"{coaddName}Coadd_psfMatchedWarp", 

90 f"{coaddName}Coadd_directWarp_sub", f"{coaddName}Coadd_psfMatchedWarp_sub") 

91 

92 self.metadataTypes = (f"{coaddName}Coadd_directWarp_visitInfo", 

93 f"{coaddName}Coadd_psfMatchedWarp_visitInfo") 

94 

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} 

101 

102 def get(self, datasetType, bbox=None, **kwargs): 

103 """Retrieve the specified dataset using the API of the Gen 2 Butler. 

104 

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. 

115 

116 Returns 

117 ------- 

118 `lsst.afw.image.Exposure` or `lsst.afw.image.VisitInfo` 

119 Either the exposure or its metadata, depending on the datasetType. 

120 

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}.") 

144 

145 @property 

146 def dataId(self): 

147 """Generate a valid data identifier. 

148 

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) 

155 

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 

163 

164 

165class MockWarpReference(lsst.daf.butler.DeferredDatasetHandle): 

166 """Very simple object that looks like a Gen 3 data reference to a warped 

167 exposure. 

168 

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 

189 

190 def get(self, bbox=None, component=None, parameters=None, **kwargs): 

191 """Retrieve the specified dataset using the API of the Gen 3 Butler. 

192 

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. 

205 

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 

224 

225 @property 

226 def dataId(self): 

227 """Generate a valid data identifier. 

228 

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) 

235 

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`") 

239 

240 

241def makeMockSkyInfo(bbox, wcs, patch): 

242 """Construct a `Struct` containing the geometry of the patch to be coadded. 

243 

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. 

250 

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 

261 

262 

263class MockCoaddTestData: 

264 """Generate repeatable simulated exposures with consistent metadata that 

265 are realistic enough to test the image coaddition algorithms. 

266 

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'. 

274 

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. 

314 

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`)." 

350 

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() 

370 

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 

374 

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 

382 

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 

388 

389 def setDummyCoaddInputs(self, exposure, expId): 

390 """Generate an `ExposureCatalog` as though the exposures had been 

391 processed using `warpAndPsfMatch`. 

392 

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) 

402 

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) 

408 

409 def makeCoaddTempExp(self, rawExposure, visitInfo, expId): 

410 """Add the metadata required by `AssembleCoaddTask` to an exposure. 

411 

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. 

420 

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) 

428 

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 

435 

436 def makeDummyWcs(self, rotAngle=None, pixelScale=None, crval=None, flipX=True): 

437 """Make a World Coordinate System object for testing. 

438 

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. 

449 

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 

465 

466 def makeDummyVisitInfo(self, exposureId, randomizeTime=False): 

467 """Make a self-consistent visitInfo object for testing. 

468 

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. 

475 

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) 

490 

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 

514 

515 def defineFilter(self, lambdaEff=476.31, lambdaMin=405., lambdaMax=552., filterName="gTest"): 

516 """Construct a `Filter` with sufficient information to calculate DCR. 

517 

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. 

530 

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) 

539 

540 def makeTestImage(self, expId, noiseLevel=None, psfSize=None, backgroundLevel=None, 

541 detectionSigma=5.): 

542 """Make a reproduceable PSF-convolved masked image for testing. 

543 

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) 

562 

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 

584 

585 exposure = self.makeCoaddTempExp(model, visitInfo, expId) 

586 matchedExposure = self.makeCoaddTempExp(modelPsfMatched, visitInfo, expId) 

587 return exposure, matchedExposure 

588 

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. 

593 

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'. 

602 

603 Returns 

604 ------- 

605 dataRefList : `list` of `MockGen2WarpReference` 

606 The data references. 

607 

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 

617 

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. 

622 

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'. 

633 

634 Returns 

635 ------- 

636 dataRefList : `list` of `MockWarpReference` 

637 The data references. 

638 

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