Coverage for python/lsst/pipe/tasks/insertFakes.py: 17%

Shortcuts 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

404 statements  

1# This file is part of pipe tasks 

2# 

3# Developed for the LSST Data Management System. 

4# This product includes software developed by the LSST Project 

5# (http://www.lsst.org). 

6# See the COPYRIGHT file at the top-level directory of this distribution 

7# for details of code ownership. 

8# 

9# This program is free software: you can redistribute it and/or modify 

10# it under the terms of the GNU General Public License as published by 

11# the Free Software Foundation, either version 3 of the License, or 

12# (at your option) any later version. 

13# 

14# This program is distributed in the hope that it will be useful, 

15# but WITHOUT ANY WARRANTY; without even the implied warranty of 

16# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 

17# GNU General Public License for more details. 

18# 

19# You should have received a copy of the GNU General Public License 

20# along with this program. If not, see <http://www.gnu.org/licenses/>. 

21 

22""" 

23Insert fakes into deepCoadds 

24""" 

25import galsim 

26from astropy.table import Table 

27import numpy as np 

28 

29import lsst.geom as geom 

30import lsst.afw.image as afwImage 

31import lsst.afw.math as afwMath 

32import lsst.pex.config as pexConfig 

33import lsst.pipe.base as pipeBase 

34 

35from lsst.pipe.base import CmdLineTask, PipelineTask, PipelineTaskConfig, PipelineTaskConnections 

36import lsst.pipe.base.connectionTypes as cT 

37from lsst.pex.exceptions import LogicError, InvalidParameterError 

38from lsst.coadd.utils.coaddDataIdContainer import ExistingCoaddDataIdContainer 

39from lsst.geom import SpherePoint, radians, Box2D, Point2D 

40 

41__all__ = ["InsertFakesConfig", "InsertFakesTask"] 

42 

43 

44def _add_fake_sources(exposure, objects, calibFluxRadius=12.0, logger=None): 

45 """Add fake sources to the given exposure 

46 

47 Parameters 

48 ---------- 

49 exposure : `lsst.afw.image.exposure.exposure.ExposureF` 

50 The exposure into which the fake sources should be added 

51 objects : `typing.Iterator` [`tuple` ['lsst.geom.SpherePoint`, `galsim.GSObject`]] 

52 An iterator of tuples that contains (or generates) locations and object 

53 surface brightness profiles to inject. 

54 calibFluxRadius : `float`, optional 

55 Aperture radius (in pixels) used to define the calibration for this 

56 exposure+catalog. This is used to produce the correct instrumental fluxes 

57 within the radius. The value should match that of the field defined in 

58 slot_CalibFlux_instFlux. 

59 logger : `lsst.log.log.log.Log` or `logging.Logger`, optional 

60 Logger. 

61 """ 

62 exposure.mask.addMaskPlane("FAKE") 

63 bitmask = exposure.mask.getPlaneBitMask("FAKE") 

64 if logger: 

65 logger.info(f"Adding mask plane with bitmask {bitmask}") 

66 

67 wcs = exposure.getWcs() 

68 psf = exposure.getPsf() 

69 

70 bbox = exposure.getBBox() 

71 fullBounds = galsim.BoundsI(bbox.minX, bbox.maxX, bbox.minY, bbox.maxY) 

72 gsImg = galsim.Image(exposure.image.array, bounds=fullBounds) 

73 

74 pixScale = wcs.getPixelScale().asArcseconds() 

75 

76 for spt, gsObj in objects: 

77 pt = wcs.skyToPixel(spt) 

78 posd = galsim.PositionD(pt.x, pt.y) 

79 posi = galsim.PositionI(pt.x//1, pt.y//1) 

80 if logger: 

81 logger.debug(f"Adding fake source at {pt}") 

82 

83 mat = wcs.linearizePixelToSky(spt, geom.arcseconds).getMatrix() 

84 gsWCS = galsim.JacobianWCS(mat[0, 0], mat[0, 1], mat[1, 0], mat[1, 1]) 

85 

86 # This check is here because sometimes the WCS 

87 # is multivalued and objects that should not be 

88 # were being included. 

89 gsPixScale = np.sqrt(gsWCS.pixelArea()) 

90 if gsPixScale < pixScale/2 or gsPixScale > pixScale*2: 

91 continue 

92 

93 try: 

94 psfArr = psf.computeKernelImage(pt).array 

95 except InvalidParameterError: 

96 # Try mapping to nearest point contained in bbox. 

97 contained_pt = Point2D( 

98 np.clip(pt.x, bbox.minX, bbox.maxX), 

99 np.clip(pt.y, bbox.minY, bbox.maxY) 

100 ) 

101 if pt == contained_pt: # no difference, so skip immediately 

102 if logger: 

103 logger.infof( 

104 "Cannot compute Psf for object at {}; skipping", 

105 pt 

106 ) 

107 continue 

108 # otherwise, try again with new point 

109 try: 

110 psfArr = psf.computeKernelImage(contained_pt).array 

111 except InvalidParameterError: 

112 if logger: 

113 logger.infof( 

114 "Cannot compute Psf for object at {}; skipping", 

115 pt 

116 ) 

117 continue 

118 

119 apCorr = psf.computeApertureFlux(calibFluxRadius) 

120 psfArr /= apCorr 

121 gsPSF = galsim.InterpolatedImage(galsim.Image(psfArr), wcs=gsWCS) 

122 

123 conv = galsim.Convolve(gsObj, gsPSF) 

124 stampSize = conv.getGoodImageSize(gsWCS.minLinearScale()) 

125 subBounds = galsim.BoundsI(posi).withBorder(stampSize//2) 

126 subBounds &= fullBounds 

127 

128 if subBounds.area() > 0: 

129 subImg = gsImg[subBounds] 

130 offset = posd - subBounds.true_center 

131 # Note, for calexp injection, pixel is already part of the PSF and 

132 # for coadd injection, it's incorrect to include the output pixel. 

133 # So for both cases, we draw using method='no_pixel'. 

134 

135 conv.drawImage( 

136 subImg, 

137 add_to_image=True, 

138 offset=offset, 

139 wcs=gsWCS, 

140 method='no_pixel' 

141 ) 

142 

143 subBox = geom.Box2I( 

144 geom.Point2I(subBounds.xmin, subBounds.ymin), 

145 geom.Point2I(subBounds.xmax, subBounds.ymax) 

146 ) 

147 exposure[subBox].mask.array |= bitmask 

148 

149 

150def _isWCSGalsimDefault(wcs, hdr): 

151 """Decide if wcs = galsim.PixelScale(1.0) is explicitly present in header, 

152 or if it's just the galsim default. 

153 

154 Parameters 

155 ---------- 

156 wcs : galsim.BaseWCS 

157 Potentially default WCS. 

158 hdr : galsim.fits.FitsHeader 

159 Header as read in by galsim. 

160 

161 Returns 

162 ------- 

163 isDefault : bool 

164 True if default, False if explicitly set in header. 

165 """ 

166 if wcs != galsim.PixelScale(1.0): 

167 return False 

168 if hdr.get('GS_WCS') is not None: 

169 return False 

170 if hdr.get('CTYPE1', 'LINEAR') == 'LINEAR': 

171 return not any(k in hdr for k in ['CD1_1', 'CDELT1']) 

172 for wcs_type in galsim.fitswcs.fits_wcs_types: 

173 # If one of these succeeds, then assume result is explicit 

174 try: 

175 wcs_type._readHeader(hdr) 

176 return False 

177 except Exception: 

178 pass 

179 else: 

180 return not any(k in hdr for k in ['CD1_1', 'CDELT1']) 

181 

182 

183class InsertFakesConnections(PipelineTaskConnections, 

184 defaultTemplates={"coaddName": "deep", 

185 "fakesType": "fakes_"}, 

186 dimensions=("tract", "patch", "band", "skymap")): 

187 

188 image = cT.Input( 

189 doc="Image into which fakes are to be added.", 

190 name="{coaddName}Coadd", 

191 storageClass="ExposureF", 

192 dimensions=("tract", "patch", "band", "skymap") 

193 ) 

194 

195 fakeCat = cT.Input( 

196 doc="Catalog of fake sources to draw inputs from.", 

197 name="{fakesType}fakeSourceCat", 

198 storageClass="DataFrame", 

199 dimensions=("tract", "skymap") 

200 ) 

201 

202 imageWithFakes = cT.Output( 

203 doc="Image with fake sources added.", 

204 name="{fakesType}{coaddName}Coadd", 

205 storageClass="ExposureF", 

206 dimensions=("tract", "patch", "band", "skymap") 

207 ) 

208 

209 

210class InsertFakesConfig(PipelineTaskConfig, 

211 pipelineConnections=InsertFakesConnections): 

212 """Config for inserting fake sources 

213 """ 

214 

215 # Unchanged 

216 

217 doCleanCat = pexConfig.Field( 

218 doc="If true removes bad sources from the catalog.", 

219 dtype=bool, 

220 default=True, 

221 ) 

222 

223 fakeType = pexConfig.Field( 

224 doc="What type of fake catalog to use, snapshot (includes variability in the magnitudes calculated " 

225 "from the MJD of the image), static (no variability) or filename for a user defined fits" 

226 "catalog.", 

227 dtype=str, 

228 default="static", 

229 ) 

230 

231 calibFluxRadius = pexConfig.Field( 

232 doc="Aperture radius (in pixels) that was used to define the calibration for this image+catalog. " 

233 "This will be used to produce the correct instrumental fluxes within the radius. " 

234 "This value should match that of the field defined in slot_CalibFlux_instFlux.", 

235 dtype=float, 

236 default=12.0, 

237 ) 

238 

239 coaddName = pexConfig.Field( 

240 doc="The name of the type of coadd used", 

241 dtype=str, 

242 default="deep", 

243 ) 

244 

245 doSubSelectSources = pexConfig.Field( 

246 doc="Set to True if you wish to sub select sources to be input based on the value in the column" 

247 "set in the sourceSelectionColName config option.", 

248 dtype=bool, 

249 default=False 

250 ) 

251 

252 insertImages = pexConfig.Field( 

253 doc="Insert images directly? True or False.", 

254 dtype=bool, 

255 default=False, 

256 ) 

257 

258 doProcessAllDataIds = pexConfig.Field( 

259 doc="If True, all input data IDs will be processed, even those containing no fake sources.", 

260 dtype=bool, 

261 default=False, 

262 ) 

263 

264 trimBuffer = pexConfig.Field( 

265 doc="Size of the pixel buffer surrounding the image. Only those fake sources with a centroid" 

266 "falling within the image+buffer region will be considered for fake source injection.", 

267 dtype=int, 

268 default=100, 

269 ) 

270 

271 sourceType = pexConfig.Field( 

272 doc="The column name for the source type used in the fake source catalog.", 

273 dtype=str, 

274 default="sourceType", 

275 ) 

276 

277 fits_alignment = pexConfig.ChoiceField( 

278 doc="How should injections from FITS files be aligned?", 

279 dtype=str, 

280 allowed={ 

281 "wcs": ( 

282 "Input image will be transformed such that the local WCS in " 

283 "the FITS header matches the local WCS in the target image. " 

284 "I.e., North, East, and angular distances in the input image " 

285 "will match North, East, and angular distances in the target " 

286 "image." 

287 ), 

288 "pixel": ( 

289 "Input image will _not_ be transformed. Up, right, and pixel " 

290 "distances in the input image will match up, right and pixel " 

291 "distances in the target image." 

292 ) 

293 }, 

294 default="pixel" 

295 ) 

296 

297 # New source catalog config variables 

298 

299 ra_col = pexConfig.Field( 

300 doc="Source catalog column name for RA (in radians).", 

301 dtype=str, 

302 default="ra", 

303 ) 

304 

305 dec_col = pexConfig.Field( 

306 doc="Source catalog column name for dec (in radians).", 

307 dtype=str, 

308 default="dec", 

309 ) 

310 

311 bulge_semimajor_col = pexConfig.Field( 

312 doc="Source catalog column name for the semimajor axis (in arcseconds) " 

313 "of the bulge half-light ellipse.", 

314 dtype=str, 

315 default="bulge_semimajor", 

316 ) 

317 

318 bulge_axis_ratio_col = pexConfig.Field( 

319 doc="Source catalog column name for the axis ratio of the bulge " 

320 "half-light ellipse.", 

321 dtype=str, 

322 default="bulge_axis_ratio", 

323 ) 

324 

325 bulge_pa_col = pexConfig.Field( 

326 doc="Source catalog column name for the position angle (measured from " 

327 "North through East in degrees) of the semimajor axis of the bulge " 

328 "half-light ellipse.", 

329 dtype=str, 

330 default="bulge_pa", 

331 ) 

332 

333 bulge_n_col = pexConfig.Field( 

334 doc="Source catalog column name for the Sersic index of the bulge.", 

335 dtype=str, 

336 default="bulge_n", 

337 ) 

338 

339 disk_semimajor_col = pexConfig.Field( 

340 doc="Source catalog column name for the semimajor axis (in arcseconds) " 

341 "of the disk half-light ellipse.", 

342 dtype=str, 

343 default="disk_semimajor", 

344 ) 

345 

346 disk_axis_ratio_col = pexConfig.Field( 

347 doc="Source catalog column name for the axis ratio of the disk " 

348 "half-light ellipse.", 

349 dtype=str, 

350 default="disk_axis_ratio", 

351 ) 

352 

353 disk_pa_col = pexConfig.Field( 

354 doc="Source catalog column name for the position angle (measured from " 

355 "North through East in degrees) of the semimajor axis of the disk " 

356 "half-light ellipse.", 

357 dtype=str, 

358 default="disk_pa", 

359 ) 

360 

361 disk_n_col = pexConfig.Field( 

362 doc="Source catalog column name for the Sersic index of the disk.", 

363 dtype=str, 

364 default="disk_n", 

365 ) 

366 

367 bulge_disk_flux_ratio_col = pexConfig.Field( 

368 doc="Source catalog column name for the bulge/disk flux ratio.", 

369 dtype=str, 

370 default="bulge_disk_flux_ratio", 

371 ) 

372 

373 mag_col = pexConfig.Field( 

374 doc="Source catalog column name template for magnitudes, in the format " 

375 "``filter name``_mag_col. E.g., if this config variable is set to " 

376 "``%s_mag``, then the i-band magnitude will be searched for in the " 

377 "``i_mag`` column of the source catalog.", 

378 dtype=str, 

379 default="%s_mag" 

380 ) 

381 

382 select_col = pexConfig.Field( 

383 doc="Source catalog column name to be used to select which sources to " 

384 "add.", 

385 dtype=str, 

386 default="select", 

387 ) 

388 

389 # Deprecated config variables 

390 

391 raColName = pexConfig.Field( 

392 doc="RA column name used in the fake source catalog.", 

393 dtype=str, 

394 default="raJ2000", 

395 deprecated="Use `ra_col` instead." 

396 ) 

397 

398 decColName = pexConfig.Field( 

399 doc="Dec. column name used in the fake source catalog.", 

400 dtype=str, 

401 default="decJ2000", 

402 deprecated="Use `dec_col` instead." 

403 ) 

404 

405 diskHLR = pexConfig.Field( 

406 doc="Column name for the disk half light radius used in the fake source catalog.", 

407 dtype=str, 

408 default="DiskHalfLightRadius", 

409 deprecated=( 

410 "Use `disk_semimajor_col`, `disk_axis_ratio_col`, and `disk_pa_col`" 

411 " to specify disk half-light ellipse." 

412 ) 

413 ) 

414 

415 aDisk = pexConfig.Field( 

416 doc="The column name for the semi major axis length of the disk component used in the fake source" 

417 "catalog.", 

418 dtype=str, 

419 default="a_d", 

420 deprecated=( 

421 "Use `disk_semimajor_col`, `disk_axis_ratio_col`, and `disk_pa_col`" 

422 " to specify disk half-light ellipse." 

423 ) 

424 ) 

425 

426 bDisk = pexConfig.Field( 

427 doc="The column name for the semi minor axis length of the disk component.", 

428 dtype=str, 

429 default="b_d", 

430 deprecated=( 

431 "Use `disk_semimajor_col`, `disk_axis_ratio_col`, and `disk_pa_col`" 

432 " to specify disk half-light ellipse." 

433 ) 

434 ) 

435 

436 paDisk = pexConfig.Field( 

437 doc="The column name for the PA of the disk component used in the fake source catalog.", 

438 dtype=str, 

439 default="pa_disk", 

440 deprecated=( 

441 "Use `disk_semimajor_col`, `disk_axis_ratio_col`, and `disk_pa_col`" 

442 " to specify disk half-light ellipse." 

443 ) 

444 ) 

445 

446 nDisk = pexConfig.Field( 

447 doc="The column name for the sersic index of the disk component used in the fake source catalog.", 

448 dtype=str, 

449 default="disk_n", 

450 deprecated="Use `disk_n_col` instead." 

451 ) 

452 

453 bulgeHLR = pexConfig.Field( 

454 doc="Column name for the bulge half light radius used in the fake source catalog.", 

455 dtype=str, 

456 default="BulgeHalfLightRadius", 

457 deprecated=( 

458 "Use `bulge_semimajor_col`, `bulge_axis_ratio_col`, and " 

459 "`bulge_pa_col` to specify disk half-light ellipse." 

460 ) 

461 ) 

462 

463 aBulge = pexConfig.Field( 

464 doc="The column name for the semi major axis length of the bulge component.", 

465 dtype=str, 

466 default="a_b", 

467 deprecated=( 

468 "Use `bulge_semimajor_col`, `bulge_axis_ratio_col`, and " 

469 "`bulge_pa_col` to specify disk half-light ellipse." 

470 ) 

471 ) 

472 

473 bBulge = pexConfig.Field( 

474 doc="The column name for the semi minor axis length of the bulge component used in the fake source " 

475 "catalog.", 

476 dtype=str, 

477 default="b_b", 

478 deprecated=( 

479 "Use `bulge_semimajor_col`, `bulge_axis_ratio_col`, and " 

480 "`bulge_pa_col` to specify disk half-light ellipse." 

481 ) 

482 ) 

483 

484 paBulge = pexConfig.Field( 

485 doc="The column name for the PA of the bulge component used in the fake source catalog.", 

486 dtype=str, 

487 default="pa_bulge", 

488 deprecated=( 

489 "Use `bulge_semimajor_col`, `bulge_axis_ratio_col`, and " 

490 "`bulge_pa_col` to specify disk half-light ellipse." 

491 ) 

492 ) 

493 

494 nBulge = pexConfig.Field( 

495 doc="The column name for the sersic index of the bulge component used in the fake source catalog.", 

496 dtype=str, 

497 default="bulge_n", 

498 deprecated="Use `bulge_n_col` instead." 

499 ) 

500 

501 magVar = pexConfig.Field( 

502 doc="The column name for the magnitude calculated taking variability into account. In the format " 

503 "``filter name``magVar, e.g. imagVar for the magnitude in the i band.", 

504 dtype=str, 

505 default="%smagVar", 

506 deprecated="Use `mag_col` instead." 

507 ) 

508 

509 sourceSelectionColName = pexConfig.Field( 

510 doc="The name of the column in the input fakes catalogue to be used to determine which sources to" 

511 "add, default is none and when this is used all sources are added.", 

512 dtype=str, 

513 default="templateSource", 

514 deprecated="Use `select_col` instead." 

515 ) 

516 

517 

518class InsertFakesTask(PipelineTask, CmdLineTask): 

519 """Insert fake objects into images. 

520 

521 Add fake stars and galaxies to the given image, read in through the dataRef. Galaxy parameters are read in 

522 from the specified file and then modelled using galsim. 

523 

524 `InsertFakesTask` has five functions that make images of the fake sources and then add them to the 

525 image. 

526 

527 `addPixCoords` 

528 Use the WCS information to add the pixel coordinates of each source. 

529 `mkFakeGalsimGalaxies` 

530 Use Galsim to make fake double sersic galaxies for each set of galaxy parameters in the input file. 

531 `mkFakeStars` 

532 Use the PSF information from the image to make a fake star using the magnitude information from the 

533 input file. 

534 `cleanCat` 

535 Remove rows of the input fake catalog which have half light radius, of either the bulge or the disk, 

536 that are 0. Also removes rows that have Sersic index outside of galsim's allowed paramters. If 

537 the config option sourceSelectionColName is set then this function limits the catalog of input fakes 

538 to only those which are True in this column. 

539 `addFakeSources` 

540 Add the fake sources to the image. 

541 

542 """ 

543 

544 _DefaultName = "insertFakes" 

545 ConfigClass = InsertFakesConfig 

546 

547 def runDataRef(self, dataRef): 

548 """Read in/write out the required data products and add fake sources to the deepCoadd. 

549 

550 Parameters 

551 ---------- 

552 dataRef : `lsst.daf.persistence.butlerSubset.ButlerDataRef` 

553 Data reference defining the image to have fakes added to it 

554 Used to access the following data products: 

555 deepCoadd 

556 """ 

557 

558 self.log.info("Adding fakes to: tract: %d, patch: %s, filter: %s", 

559 dataRef.dataId["tract"], dataRef.dataId["patch"], dataRef.dataId["filter"]) 

560 

561 # To do: should it warn when asked to insert variable sources into the coadd 

562 

563 if self.config.fakeType == "static": 

564 fakeCat = dataRef.get("deepCoadd_fakeSourceCat").toDataFrame() 

565 # To do: DM-16254, the read and write of the fake catalogs will be changed once the new pipeline 

566 # task structure for ref cats is in place. 

567 self.fakeSourceCatType = "deepCoadd_fakeSourceCat" 

568 else: 

569 fakeCat = Table.read(self.config.fakeType).to_pandas() 

570 

571 coadd = dataRef.get("deepCoadd") 

572 wcs = coadd.getWcs() 

573 photoCalib = coadd.getPhotoCalib() 

574 

575 imageWithFakes = self.run(fakeCat, coadd, wcs, photoCalib) 

576 

577 dataRef.put(imageWithFakes.imageWithFakes, "fakes_deepCoadd") 

578 

579 def runQuantum(self, butlerQC, inputRefs, outputRefs): 

580 inputs = butlerQC.get(inputRefs) 

581 inputs["wcs"] = inputs["image"].getWcs() 

582 inputs["photoCalib"] = inputs["image"].getPhotoCalib() 

583 

584 outputs = self.run(**inputs) 

585 butlerQC.put(outputs, outputRefs) 

586 

587 @classmethod 

588 def _makeArgumentParser(cls): 

589 parser = pipeBase.ArgumentParser(name=cls._DefaultName) 

590 parser.add_id_argument(name="--id", datasetType="deepCoadd", 

591 help="data IDs for the deepCoadd, e.g. --id tract=12345 patch=1,2 filter=r", 

592 ContainerClass=ExistingCoaddDataIdContainer) 

593 return parser 

594 

595 def run(self, fakeCat, image, wcs, photoCalib): 

596 """Add fake sources to an image. 

597 

598 Parameters 

599 ---------- 

600 fakeCat : `pandas.core.frame.DataFrame` 

601 The catalog of fake sources to be input 

602 image : `lsst.afw.image.exposure.exposure.ExposureF` 

603 The image into which the fake sources should be added 

604 wcs : `lsst.afw.geom.SkyWcs` 

605 WCS to use to add fake sources 

606 photoCalib : `lsst.afw.image.photoCalib.PhotoCalib` 

607 Photometric calibration to be used to calibrate the fake sources 

608 

609 Returns 

610 ------- 

611 resultStruct : `lsst.pipe.base.struct.Struct` 

612 contains : image : `lsst.afw.image.exposure.exposure.ExposureF` 

613 

614 Notes 

615 ----- 

616 Adds pixel coordinates for each source to the fakeCat and removes objects with bulge or disk half 

617 light radius = 0 (if ``config.doCleanCat = True``). 

618 

619 Adds the ``Fake`` mask plane to the image which is then set by `addFakeSources` to mark where fake 

620 sources have been added. Uses the information in the ``fakeCat`` to make fake galaxies (using galsim) 

621 and fake stars, using the PSF models from the PSF information for the image. These are then added to 

622 the image and the image with fakes included returned. 

623 

624 The galsim galaxies are made using a double sersic profile, one for the bulge and one for the disk, 

625 this is then convolved with the PSF at that point. 

626 """ 

627 # Attach overriding wcs and photoCalib to image, but retain originals 

628 # so we can reset at the end. 

629 origWcs = image.getWcs() 

630 origPhotoCalib = image.getPhotoCalib() 

631 image.setWcs(wcs) 

632 image.setPhotoCalib(photoCalib) 

633 

634 band = image.getFilterLabel().bandLabel 

635 fakeCat = self._standardizeColumns(fakeCat, band) 

636 

637 fakeCat = self.addPixCoords(fakeCat, image) 

638 fakeCat = self.trimFakeCat(fakeCat, image) 

639 

640 if len(fakeCat) > 0: 

641 if not self.config.insertImages: 

642 if isinstance(fakeCat[self.config.sourceType].iloc[0], str): 

643 galCheckVal = "galaxy" 

644 starCheckVal = "star" 

645 elif isinstance(fakeCat[self.config.sourceType].iloc[0], bytes): 

646 galCheckVal = b"galaxy" 

647 starCheckVal = b"star" 

648 elif isinstance(fakeCat[self.config.sourceType].iloc[0], (int, float)): 

649 galCheckVal = 1 

650 starCheckVal = 0 

651 else: 

652 raise TypeError( 

653 "sourceType column does not have required type, should be str, bytes or int" 

654 ) 

655 if self.config.doCleanCat: 

656 fakeCat = self.cleanCat(fakeCat, starCheckVal) 

657 

658 generator = self._generateGSObjectsFromCatalog(image, fakeCat, galCheckVal, starCheckVal) 

659 else: 

660 generator = self._generateGSObjectsFromImages(image, fakeCat) 

661 _add_fake_sources(image, generator, calibFluxRadius=self.config.calibFluxRadius, logger=self.log) 

662 elif len(fakeCat) == 0 and self.config.doProcessAllDataIds: 

663 self.log.warning("No fakes found for this dataRef; processing anyway.") 

664 image.mask.addMaskPlane("FAKE") 

665 else: 

666 raise RuntimeError("No fakes found for this dataRef.") 

667 

668 # restore original exposure WCS and photoCalib 

669 image.setWcs(origWcs) 

670 image.setPhotoCalib(origPhotoCalib) 

671 

672 resultStruct = pipeBase.Struct(imageWithFakes=image) 

673 

674 return resultStruct 

675 

676 def _standardizeColumns(self, fakeCat, band): 

677 """Use config variables to 'standardize' the expected columns and column 

678 names in the input catalog. 

679 

680 Parameters 

681 ---------- 

682 fakeCat : `pandas.core.frame.DataFrame` 

683 The catalog of fake sources to be input 

684 band : `str` 

685 Label for the current band being processed. 

686 

687 Returns 

688 ------- 

689 outCat : `pandas.core.frame.DataFrame` 

690 The standardized catalog of fake sources 

691 """ 

692 cfg = self.config 

693 replace_dict = {} 

694 

695 def add_to_replace_dict(new_name, depr_name, std_name): 

696 if new_name in fakeCat.columns: 

697 replace_dict[new_name] = std_name 

698 elif depr_name in fakeCat.columns: 

699 replace_dict[depr_name] = std_name 

700 else: 

701 raise ValueError(f"Could not determine column for {std_name}.") 

702 

703 # Prefer new config variables over deprecated config variables. 

704 # RA, dec, and mag are always required. Do these first 

705 for new_name, depr_name, std_name in [ 

706 (cfg.ra_col, cfg.raColName, 'ra'), 

707 (cfg.dec_col, cfg.decColName, 'dec'), 

708 (cfg.mag_col%band, cfg.magVar%band, 'mag') 

709 ]: 

710 add_to_replace_dict(new_name, depr_name, std_name) 

711 # Only handle bulge/disk params if not injecting images 

712 if not cfg.insertImages: 

713 for new_name, depr_name, std_name in [ 

714 (cfg.bulge_n_col, cfg.nBulge, 'bulge_n'), 

715 (cfg.bulge_pa_col, cfg.paBulge, 'bulge_pa'), 

716 (cfg.disk_n_col, cfg.nDisk, 'disk_n'), 

717 (cfg.disk_pa_col, cfg.paDisk, 'disk_pa'), 

718 ]: 

719 add_to_replace_dict(new_name, depr_name, std_name) 

720 

721 if cfg.doSubSelectSources: 

722 add_to_replace_dict( 

723 cfg.select_col, 

724 cfg.sourceSelectionColName, 

725 'select' 

726 ) 

727 fakeCat = fakeCat.rename(columns=replace_dict, copy=False) 

728 

729 # Handling the half-light radius and axis-ratio are trickier, since we 

730 # moved from expecting (HLR, a, b) to expecting (semimajor, axis_ratio). 

731 # Just handle these manually. 

732 if not cfg.insertImages: 

733 if ( 

734 cfg.bulge_semimajor_col in fakeCat.columns 

735 and cfg.bulge_axis_ratio_col in fakeCat.columns 

736 ): 

737 fakeCat = fakeCat.rename( 

738 columns={ 

739 cfg.bulge_semimajor_col: 'bulge_semimajor', 

740 cfg.bulge_axis_ratio_col: 'bulge_axis_ratio', 

741 cfg.disk_semimajor_col: 'disk_semimajor', 

742 cfg.disk_axis_ratio_col: 'disk_axis_ratio', 

743 }, 

744 copy=False 

745 ) 

746 elif ( 

747 cfg.bulgeHLR in fakeCat.columns 

748 and cfg.aBulge in fakeCat.columns 

749 and cfg.bBulge in fakeCat.columns 

750 ): 

751 fakeCat['bulge_axis_ratio'] = ( 

752 fakeCat[cfg.bBulge]/fakeCat[cfg.aBulge] 

753 ) 

754 fakeCat['bulge_semimajor'] = ( 

755 fakeCat[cfg.bulgeHLR]/np.sqrt(fakeCat['bulge_axis_ratio']) 

756 ) 

757 fakeCat['disk_axis_ratio'] = ( 

758 fakeCat[cfg.bDisk]/fakeCat[cfg.aDisk] 

759 ) 

760 fakeCat['disk_semimajor'] = ( 

761 fakeCat[cfg.diskHLR]/np.sqrt(fakeCat['disk_axis_ratio']) 

762 ) 

763 else: 

764 raise ValueError( 

765 "Could not determine columns for half-light radius and " 

766 "axis ratio." 

767 ) 

768 

769 # Process the bulge/disk flux ratio if possible. 

770 if cfg.bulge_disk_flux_ratio_col in fakeCat.columns: 

771 fakeCat = fakeCat.rename( 

772 columns={ 

773 cfg.bulge_disk_flux_ratio_col: 'bulge_disk_flux_ratio' 

774 }, 

775 copy=False 

776 ) 

777 else: 

778 fakeCat['bulge_disk_flux_ratio'] = 1.0 

779 

780 return fakeCat 

781 

782 def _generateGSObjectsFromCatalog(self, exposure, fakeCat, galCheckVal, starCheckVal): 

783 """Process catalog to generate `galsim.GSObject` s. 

784 

785 Parameters 

786 ---------- 

787 exposure : `lsst.afw.image.exposure.exposure.ExposureF` 

788 The exposure into which the fake sources should be added 

789 fakeCat : `pandas.core.frame.DataFrame` 

790 The catalog of fake sources to be input 

791 galCheckVal : `str`, `bytes` or `int` 

792 The value that is set in the sourceType column to specifiy an object is a galaxy. 

793 starCheckVal : `str`, `bytes` or `int` 

794 The value that is set in the sourceType column to specifiy an object is a star. 

795 

796 Yields 

797 ------ 

798 gsObjects : `generator` 

799 A generator of tuples of `lsst.geom.SpherePoint` and `galsim.GSObject`. 

800 """ 

801 wcs = exposure.getWcs() 

802 photoCalib = exposure.getPhotoCalib() 

803 

804 self.log.info("Making %d objects for insertion", len(fakeCat)) 

805 

806 for (index, row) in fakeCat.iterrows(): 

807 ra = row['ra'] 

808 dec = row['dec'] 

809 skyCoord = SpherePoint(ra, dec, radians) 

810 xy = wcs.skyToPixel(skyCoord) 

811 

812 try: 

813 flux = photoCalib.magnitudeToInstFlux(row['mag'], xy) 

814 except LogicError: 

815 continue 

816 

817 sourceType = row[self.config.sourceType] 

818 if sourceType == galCheckVal: 

819 # GalSim convention: HLR = sqrt(a * b) = a * sqrt(b / a) 

820 bulge_gs_HLR = row['bulge_semimajor']*np.sqrt(row['bulge_axis_ratio']) 

821 bulge = galsim.Sersic(n=row['bulge_n'], half_light_radius=bulge_gs_HLR) 

822 bulge = bulge.shear(q=row['bulge_axis_ratio'], beta=((90 - row['bulge_pa'])*galsim.degrees)) 

823 

824 disk_gs_HLR = row['disk_semimajor']*np.sqrt(row['disk_axis_ratio']) 

825 disk = galsim.Sersic(n=row['disk_n'], half_light_radius=disk_gs_HLR) 

826 disk = disk.shear(q=row['disk_axis_ratio'], beta=((90 - row['disk_pa'])*galsim.degrees)) 

827 

828 gal = bulge*row['bulge_disk_flux_ratio'] + disk 

829 gal = gal.withFlux(flux) 

830 

831 yield skyCoord, gal 

832 elif sourceType == starCheckVal: 

833 star = galsim.DeltaFunction() 

834 star = star.withFlux(flux) 

835 yield skyCoord, star 

836 else: 

837 raise TypeError(f"Unknown sourceType {sourceType}") 

838 

839 def _generateGSObjectsFromImages(self, exposure, fakeCat): 

840 """Process catalog to generate `galsim.GSObject` s. 

841 

842 Parameters 

843 ---------- 

844 exposure : `lsst.afw.image.exposure.exposure.ExposureF` 

845 The exposure into which the fake sources should be added 

846 fakeCat : `pandas.core.frame.DataFrame` 

847 The catalog of fake sources to be input 

848 

849 Yields 

850 ------ 

851 gsObjects : `generator` 

852 A generator of tuples of `lsst.geom.SpherePoint` and `galsim.GSObject`. 

853 """ 

854 band = exposure.getFilterLabel().bandLabel 

855 wcs = exposure.getWcs() 

856 photoCalib = exposure.getPhotoCalib() 

857 

858 self.log.info("Processing %d fake images", len(fakeCat)) 

859 

860 for (index, row) in fakeCat.iterrows(): 

861 ra = row['ra'] 

862 dec = row['dec'] 

863 skyCoord = SpherePoint(ra, dec, radians) 

864 xy = wcs.skyToPixel(skyCoord) 

865 

866 try: 

867 flux = photoCalib.magnitudeToInstFlux(row['mag'], xy) 

868 except LogicError: 

869 continue 

870 

871 imFile = row[band+"imFilename"] 

872 try: 

873 imFile = imFile.decode("utf-8") 

874 except AttributeError: 

875 pass 

876 imFile = imFile.strip() 

877 im = galsim.fits.read(imFile, read_header=True) 

878 

879 if self.config.fits_alignment == "wcs": 

880 # galsim.fits.read will always attach a WCS to its output. If it 

881 # can't find a WCS in the FITS header, then it defaults to 

882 # scale = 1.0 arcsec / pix. So if that's the scale, then we 

883 # need to check if it was explicitly set or if it's just the 

884 # default. If it's just the default then we should raise an 

885 # exception. 

886 if _isWCSGalsimDefault(im.wcs, im.header): 

887 raise RuntimeError( 

888 f"Cannot find WCS in input FITS file {imFile}" 

889 ) 

890 elif self.config.fits_alignment == "pixel": 

891 # Here we need to set im.wcs to the local WCS at the target 

892 # position. 

893 linWcs = wcs.linearizePixelToSky(skyCoord, geom.arcseconds) 

894 mat = linWcs.getMatrix() 

895 im.wcs = galsim.JacobianWCS( 

896 mat[0, 0], mat[0, 1], mat[1, 0], mat[1, 1] 

897 ) 

898 else: 

899 raise ValueError( 

900 f"Unknown fits_alignment type {self.config.fits_alignment}" 

901 ) 

902 

903 obj = galsim.InterpolatedImage(im, calculate_stepk=False) 

904 obj = obj.withFlux(flux) 

905 yield skyCoord, obj 

906 

907 def processImagesForInsertion(self, fakeCat, wcs, psf, photoCalib, band, pixelScale): 

908 """Process images from files into the format needed for insertion. 

909 

910 Parameters 

911 ---------- 

912 fakeCat : `pandas.core.frame.DataFrame` 

913 The catalog of fake sources to be input 

914 wcs : `lsst.afw.geom.skyWcs.skyWcs.SkyWc` 

915 WCS to use to add fake sources 

916 psf : `lsst.meas.algorithms.coaddPsf.coaddPsf.CoaddPsf` or 

917 `lsst.meas.extensions.psfex.psfexPsf.PsfexPsf` 

918 The PSF information to use to make the PSF images 

919 photoCalib : `lsst.afw.image.photoCalib.PhotoCalib` 

920 Photometric calibration to be used to calibrate the fake sources 

921 band : `str` 

922 The filter band that the observation was taken in. 

923 pixelScale : `float` 

924 The pixel scale of the image the sources are to be added to. 

925 

926 Returns 

927 ------- 

928 galImages : `list` 

929 A list of tuples of `lsst.afw.image.exposure.exposure.ExposureF` and 

930 `lsst.geom.Point2D` of their locations. 

931 For sources labelled as galaxy. 

932 starImages : `list` 

933 A list of tuples of `lsst.afw.image.exposure.exposure.ExposureF` and 

934 `lsst.geom.Point2D` of their locations. 

935 For sources labelled as star. 

936 

937 Notes 

938 ----- 

939 The input fakes catalog needs to contain the absolute path to the image in the 

940 band that is being used to add images to. It also needs to have the R.A. and 

941 declination of the fake source in radians and the sourceType of the object. 

942 """ 

943 galImages = [] 

944 starImages = [] 

945 

946 self.log.info("Processing %d fake images", len(fakeCat)) 

947 

948 for (imFile, sourceType, mag, x, y) in zip(fakeCat[band + "imFilename"].array, 

949 fakeCat["sourceType"].array, 

950 fakeCat['mag'].array, 

951 fakeCat["x"].array, fakeCat["y"].array): 

952 

953 im = afwImage.ImageF.readFits(imFile) 

954 

955 xy = geom.Point2D(x, y) 

956 

957 # We put these two PSF calculations within this same try block so that we catch cases 

958 # where the object's position is outside of the image. 

959 try: 

960 correctedFlux = psf.computeApertureFlux(self.config.calibFluxRadius, xy) 

961 psfKernel = psf.computeKernelImage(xy).getArray() 

962 psfKernel /= correctedFlux 

963 

964 except InvalidParameterError: 

965 self.log.info("%s at %0.4f, %0.4f outside of image", sourceType, x, y) 

966 continue 

967 

968 psfIm = galsim.InterpolatedImage(galsim.Image(psfKernel), scale=pixelScale) 

969 galsimIm = galsim.InterpolatedImage(galsim.Image(im.array), scale=pixelScale) 

970 convIm = galsim.Convolve([galsimIm, psfIm]) 

971 

972 try: 

973 outIm = convIm.drawImage(scale=pixelScale, method="real_space").array 

974 except (galsim.errors.GalSimFFTSizeError, MemoryError): 

975 continue 

976 

977 imSum = np.sum(outIm) 

978 divIm = outIm/imSum 

979 

980 try: 

981 flux = photoCalib.magnitudeToInstFlux(mag, xy) 

982 except LogicError: 

983 flux = 0 

984 

985 imWithFlux = flux*divIm 

986 

987 if sourceType == b"galaxy": 

988 galImages.append((afwImage.ImageF(imWithFlux), xy)) 

989 if sourceType == b"star": 

990 starImages.append((afwImage.ImageF(imWithFlux), xy)) 

991 

992 return galImages, starImages 

993 

994 def addPixCoords(self, fakeCat, image): 

995 

996 """Add pixel coordinates to the catalog of fakes. 

997 

998 Parameters 

999 ---------- 

1000 fakeCat : `pandas.core.frame.DataFrame` 

1001 The catalog of fake sources to be input 

1002 image : `lsst.afw.image.exposure.exposure.ExposureF` 

1003 The image into which the fake sources should be added 

1004 

1005 Returns 

1006 ------- 

1007 fakeCat : `pandas.core.frame.DataFrame` 

1008 """ 

1009 wcs = image.getWcs() 

1010 ras = fakeCat['ra'].values 

1011 decs = fakeCat['dec'].values 

1012 xs, ys = wcs.skyToPixelArray(ras, decs) 

1013 fakeCat["x"] = xs 

1014 fakeCat["y"] = ys 

1015 

1016 return fakeCat 

1017 

1018 def trimFakeCat(self, fakeCat, image): 

1019 """Trim the fake cat to about the size of the input image. 

1020 

1021 `fakeCat` must be processed with addPixCoords before using this method. 

1022 

1023 Parameters 

1024 ---------- 

1025 fakeCat : `pandas.core.frame.DataFrame` 

1026 The catalog of fake sources to be input 

1027 image : `lsst.afw.image.exposure.exposure.ExposureF` 

1028 The image into which the fake sources should be added 

1029 

1030 Returns 

1031 ------- 

1032 fakeCat : `pandas.core.frame.DataFrame` 

1033 The original fakeCat trimmed to the area of the image 

1034 """ 

1035 

1036 bbox = Box2D(image.getBBox()).dilatedBy(self.config.trimBuffer) 

1037 xs = fakeCat["x"].values 

1038 ys = fakeCat["y"].values 

1039 

1040 isContained = xs >= bbox.minX 

1041 isContained &= xs <= bbox.maxX 

1042 isContained &= ys >= bbox.minY 

1043 isContained &= ys <= bbox.maxY 

1044 

1045 return fakeCat[isContained] 

1046 

1047 def mkFakeGalsimGalaxies(self, fakeCat, band, photoCalib, pixelScale, psf, image): 

1048 """Make images of fake galaxies using GalSim. 

1049 

1050 Parameters 

1051 ---------- 

1052 band : `str` 

1053 pixelScale : `float` 

1054 psf : `lsst.meas.extensions.psfex.psfexPsf.PsfexPsf` 

1055 The PSF information to use to make the PSF images 

1056 fakeCat : `pandas.core.frame.DataFrame` 

1057 The catalog of fake sources to be input 

1058 photoCalib : `lsst.afw.image.photoCalib.PhotoCalib` 

1059 Photometric calibration to be used to calibrate the fake sources 

1060 

1061 Yields 

1062 ------- 

1063 galImages : `generator` 

1064 A generator of tuples of `lsst.afw.image.exposure.exposure.ExposureF` and 

1065 `lsst.geom.Point2D` of their locations. 

1066 

1067 Notes 

1068 ----- 

1069 

1070 Fake galaxies are made by combining two sersic profiles, one for the bulge and one for the disk. Each 

1071 component has an individual sersic index (n), a, b and position angle (PA). The combined profile is 

1072 then convolved with the PSF at the specified x, y position on the image. 

1073 

1074 The names of the columns in the ``fakeCat`` are configurable and are the column names from the 

1075 University of Washington simulations database as default. For more information see the doc strings 

1076 attached to the config options. 

1077 

1078 See mkFakeStars doc string for an explanation of calibration to instrumental flux. 

1079 """ 

1080 

1081 self.log.info("Making %d fake galaxy images", len(fakeCat)) 

1082 

1083 for (index, row) in fakeCat.iterrows(): 

1084 xy = geom.Point2D(row["x"], row["y"]) 

1085 

1086 # We put these two PSF calculations within this same try block so that we catch cases 

1087 # where the object's position is outside of the image. 

1088 try: 

1089 correctedFlux = psf.computeApertureFlux(self.config.calibFluxRadius, xy) 

1090 psfKernel = psf.computeKernelImage(xy).getArray() 

1091 psfKernel /= correctedFlux 

1092 

1093 except InvalidParameterError: 

1094 self.log.info("Galaxy at %0.4f, %0.4f outside of image", row["x"], row["y"]) 

1095 continue 

1096 

1097 try: 

1098 flux = photoCalib.magnitudeToInstFlux(row['mag'], xy) 

1099 except LogicError: 

1100 flux = 0 

1101 

1102 # GalSim convention: HLR = sqrt(a * b) = a * sqrt(b / a) 

1103 bulge_gs_HLR = row['bulge_semimajor']*np.sqrt(row['bulge_axis_ratio']) 

1104 bulge = galsim.Sersic(n=row['bulge_n'], half_light_radius=bulge_gs_HLR) 

1105 bulge = bulge.shear(q=row['bulge_axis_ratio'], beta=((90 - row['bulge_pa'])*galsim.degrees)) 

1106 

1107 disk_gs_HLR = row['disk_semimajor']*np.sqrt(row['disk_axis_ratio']) 

1108 disk = galsim.Sersic(n=row['disk_n'], half_light_radius=disk_gs_HLR) 

1109 disk = disk.shear(q=row['disk_axis_ratio'], beta=((90 - row['disk_pa'])*galsim.degrees)) 

1110 

1111 gal = bulge*row['bulge_disk_flux_ratio'] + disk 

1112 gal = gal.withFlux(flux) 

1113 

1114 psfIm = galsim.InterpolatedImage(galsim.Image(psfKernel), scale=pixelScale) 

1115 gal = galsim.Convolve([gal, psfIm]) 

1116 try: 

1117 galIm = gal.drawImage(scale=pixelScale, method="real_space").array 

1118 except (galsim.errors.GalSimFFTSizeError, MemoryError): 

1119 continue 

1120 

1121 yield (afwImage.ImageF(galIm), xy) 

1122 

1123 def mkFakeStars(self, fakeCat, band, photoCalib, psf, image): 

1124 

1125 """Make fake stars based off the properties in the fakeCat. 

1126 

1127 Parameters 

1128 ---------- 

1129 band : `str` 

1130 psf : `lsst.meas.extensions.psfex.psfexPsf.PsfexPsf` 

1131 The PSF information to use to make the PSF images 

1132 fakeCat : `pandas.core.frame.DataFrame` 

1133 The catalog of fake sources to be input 

1134 image : `lsst.afw.image.exposure.exposure.ExposureF` 

1135 The image into which the fake sources should be added 

1136 photoCalib : `lsst.afw.image.photoCalib.PhotoCalib` 

1137 Photometric calibration to be used to calibrate the fake sources 

1138 

1139 Yields 

1140 ------- 

1141 starImages : `generator` 

1142 A generator of tuples of `lsst.afw.image.ImageF` of fake stars and 

1143 `lsst.geom.Point2D` of their locations. 

1144 

1145 Notes 

1146 ----- 

1147 To take a given magnitude and translate to the number of counts in the image 

1148 we use photoCalib.magnitudeToInstFlux, which returns the instrumental flux for the 

1149 given calibration radius used in the photometric calibration step. 

1150 Thus `calibFluxRadius` should be set to this same radius so that we can normalize 

1151 the PSF model to the correct instrumental flux within calibFluxRadius. 

1152 """ 

1153 

1154 self.log.info("Making %d fake star images", len(fakeCat)) 

1155 

1156 for (index, row) in fakeCat.iterrows(): 

1157 xy = geom.Point2D(row["x"], row["y"]) 

1158 

1159 # We put these two PSF calculations within this same try block so that we catch cases 

1160 # where the object's position is outside of the image. 

1161 try: 

1162 correctedFlux = psf.computeApertureFlux(self.config.calibFluxRadius, xy) 

1163 starIm = psf.computeImage(xy) 

1164 starIm /= correctedFlux 

1165 

1166 except InvalidParameterError: 

1167 self.log.info("Star at %0.4f, %0.4f outside of image", row["x"], row["y"]) 

1168 continue 

1169 

1170 try: 

1171 flux = photoCalib.magnitudeToInstFlux(row['mag'], xy) 

1172 except LogicError: 

1173 flux = 0 

1174 

1175 starIm *= flux 

1176 yield ((starIm.convertF(), xy)) 

1177 

1178 def cleanCat(self, fakeCat, starCheckVal): 

1179 """Remove rows from the fakes catalog which have HLR = 0 for either the buldge or disk component, 

1180 also remove galaxies that have Sersic index outside the galsim min and max 

1181 allowed (0.3 <= n <= 6.2). 

1182 

1183 Parameters 

1184 ---------- 

1185 fakeCat : `pandas.core.frame.DataFrame` 

1186 The catalog of fake sources to be input 

1187 starCheckVal : `str`, `bytes` or `int` 

1188 The value that is set in the sourceType column to specifiy an object is a star. 

1189 

1190 Returns 

1191 ------- 

1192 fakeCat : `pandas.core.frame.DataFrame` 

1193 The input catalog of fake sources but with the bad objects removed 

1194 """ 

1195 

1196 rowsToKeep = (((fakeCat['bulge_semimajor'] != 0.0) & (fakeCat['disk_semimajor'] != 0.0)) 

1197 | (fakeCat[self.config.sourceType] == starCheckVal)) 

1198 numRowsNotUsed = len(fakeCat) - len(np.where(rowsToKeep)[0]) 

1199 self.log.info("Removing %d rows with HLR = 0 for either the bulge or disk", numRowsNotUsed) 

1200 fakeCat = fakeCat[rowsToKeep] 

1201 

1202 minN = galsim.Sersic._minimum_n 

1203 maxN = galsim.Sersic._maximum_n 

1204 rowsWithGoodSersic = (((fakeCat['bulge_n'] >= minN) & (fakeCat['bulge_n'] <= maxN) 

1205 & (fakeCat['disk_n'] >= minN) & (fakeCat['disk_n'] <= maxN)) 

1206 | (fakeCat[self.config.sourceType] == starCheckVal)) 

1207 numRowsNotUsed = len(fakeCat) - len(np.where(rowsWithGoodSersic)[0]) 

1208 self.log.info("Removing %d rows of galaxies with nBulge or nDisk outside of %0.2f <= n <= %0.2f", 

1209 numRowsNotUsed, minN, maxN) 

1210 fakeCat = fakeCat[rowsWithGoodSersic] 

1211 

1212 if self.config.doSubSelectSources: 

1213 numRowsNotUsed = len(fakeCat) - len(fakeCat['select']) 

1214 self.log.info("Removing %d rows which were not designated as template sources", numRowsNotUsed) 

1215 fakeCat = fakeCat[fakeCat['select']] 

1216 

1217 return fakeCat 

1218 

1219 def addFakeSources(self, image, fakeImages, sourceType): 

1220 """Add the fake sources to the given image 

1221 

1222 Parameters 

1223 ---------- 

1224 image : `lsst.afw.image.exposure.exposure.ExposureF` 

1225 The image into which the fake sources should be added 

1226 fakeImages : `typing.Iterator` [`tuple` ['lsst.afw.image.ImageF`, `lsst.geom.Point2d`]] 

1227 An iterator of tuples that contains (or generates) images of fake sources, 

1228 and the locations they are to be inserted at. 

1229 sourceType : `str` 

1230 The type (star/galaxy) of fake sources input 

1231 

1232 Returns 

1233 ------- 

1234 image : `lsst.afw.image.exposure.exposure.ExposureF` 

1235 

1236 Notes 

1237 ----- 

1238 Uses the x, y information in the ``fakeCat`` to position an image of the fake interpolated onto the 

1239 pixel grid of the image. Sets the ``FAKE`` mask plane for the pixels added with the fake source. 

1240 """ 

1241 

1242 imageBBox = image.getBBox() 

1243 imageMI = image.maskedImage 

1244 

1245 for (fakeImage, xy) in fakeImages: 

1246 X0 = xy.getX() - fakeImage.getWidth()/2 + 0.5 

1247 Y0 = xy.getY() - fakeImage.getHeight()/2 + 0.5 

1248 self.log.debug("Adding fake source at %d, %d", xy.getX(), xy.getY()) 

1249 if sourceType == "galaxy": 

1250 interpFakeImage = afwMath.offsetImage(fakeImage, X0, Y0, "lanczos3") 

1251 else: 

1252 interpFakeImage = fakeImage 

1253 

1254 interpFakeImBBox = interpFakeImage.getBBox() 

1255 interpFakeImBBox.clip(imageBBox) 

1256 

1257 if interpFakeImBBox.getArea() > 0: 

1258 imageMIView = imageMI[interpFakeImBBox] 

1259 clippedFakeImage = interpFakeImage[interpFakeImBBox] 

1260 clippedFakeImageMI = afwImage.MaskedImageF(clippedFakeImage) 

1261 clippedFakeImageMI.mask.set(self.bitmask) 

1262 imageMIView += clippedFakeImageMI 

1263 

1264 return image 

1265 

1266 def _getMetadataName(self): 

1267 """Disable metadata writing""" 

1268 return None