Coverage for tests/test_coadds.py: 20%

260 statements  

« prev     ^ index     » next       coverage.py v7.10.6, created at 2025-10-11 02:08 -0700

1# This file is part of cell_coadds. 

2# 

3# Developed for the LSST Data Management System. 

4# This product includes software developed by the LSST Project 

5# (https://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 <https://www.gnu.org/licenses/>. 

21 

22import unittest 

23from collections.abc import Iterable, Mapping 

24from itertools import product 

25 

26import numpy as np 

27from frozendict import frozendict 

28 

29import lsst.cell_coadds.test_utils as test_utils 

30import lsst.geom as geom 

31import lsst.meas.base.tests 

32import lsst.utils.tests 

33from lsst.afw.geom import Quadrupole 

34from lsst.afw.image import ExposureF, ImageF 

35from lsst.cell_coadds import ( 

36 CellCoaddFitsReader, 

37 CellIdentifiers, 

38 CoaddUnits, 

39 CommonComponents, 

40 ExplodedCoadd, 

41 MultipleCellCoadd, 

42 ObservationIdentifiers, 

43 OwnedImagePlanes, 

44 PatchIdentifiers, 

45 SingleCellCoadd, 

46 StitchedCoadd, 

47 UniformGrid, 

48) 

49from lsst.meas.algorithms import SingleGaussianPsf 

50from lsst.skymap import Index2D 

51 

52 

53class BaseMultipleCellCoaddTestCase(lsst.utils.tests.TestCase): 

54 """A base class that provides a common set of methods.""" 

55 

56 psf_size: int 

57 psf_sigmas: Mapping[Index2D, float] 

58 border_size: int 

59 inner_size: int 

60 outer_size: int 

61 test_positions: Iterable[tuple[geom.Point2D, Index2D]] 

62 exposures: Mapping[Index2D, ExposureF] 

63 multiple_cell_coadd: MultipleCellCoadd 

64 

65 @classmethod 

66 def setUpClass(cls) -> None: 

67 """Set up a multiple cell coadd with 2x2 cells.""" 

68 np.random.seed(42) 

69 data_id = test_utils.generate_data_id() 

70 common = CommonComponents( 

71 units=CoaddUnits.nJy, # units here are arbitrary. 

72 wcs=test_utils.generate_wcs(), 

73 band=data_id["band"], 

74 identifiers=PatchIdentifiers.from_data_id(data_id), 

75 ) 

76 

77 cls.nx, cls.ny = 3, 2 

78 cls.psf_sigmas = { 

79 Index2D(x=0, y=0): 1.2, 

80 Index2D(x=0, y=1): 0.7, 

81 Index2D(x=1, y=0): 0.9, 

82 Index2D(x=1, y=1): 1.1, 

83 Index2D(x=2, y=0): 1.3, 

84 Index2D(x=2, y=1): 0.8, 

85 } 

86 

87 cls.border_size = 5 

88 # In practice, we expect this to be squares. 

89 # To check for any possible x, y swap, we use different values. 

90 cls.psf_size_x, cls.psf_size_y = 21, 23 

91 cls.inner_size_x, cls.inner_size_y = 17, 15 

92 cls.outer_size_x = cls.inner_size_x + 2 * cls.border_size 

93 cls.outer_size_y = cls.inner_size_y + 2 * cls.border_size 

94 # The origin should not be at (0, 0) for robust testing. 

95 cls.x0, cls.y0 = 5, 2 

96 cls.n_noise_realizations = 2 

97 

98 patch_outer_bbox = geom.Box2I( 

99 geom.Point2I(cls.x0, cls.y0), geom.Extent2I(cls.nx * cls.inner_size_x, cls.ny * cls.inner_size_y) 

100 ) 

101 patch_outer_bbox.grow(cls.border_size) 

102 

103 # Add one star and one galaxy per quadrant. 

104 # The mapping of positions to cell indices assume inner_size = (17, 15) 

105 # and border_size = 5. If that is changed, these values need an update. 

106 sources = ( 

107 # flux, centroid, shape 

108 (1000.0, geom.Point2D(cls.x0 + 6.3, cls.y0 + 7.2), None), 

109 (2500.0, geom.Point2D(cls.x0 + 16.8, cls.y0 + 18.3), None), 

110 (1500.0, geom.Point2D(cls.x0 + 21.2, cls.y0 + 5.1), None), 

111 (3200.0, geom.Point2D(cls.x0 + 16.1, cls.y0 + 23.9), None), 

112 (1800.0, geom.Point2D(cls.x0 + 44.7, cls.y0 + 8.9), None), 

113 (2100.0, geom.Point2D(cls.x0 + 34.1, cls.y0 + 19.2), None), 

114 (900.0, geom.Point2D(cls.x0 + 9.1, cls.y0 + 13.9), Quadrupole(2.5, 1.5, 0.8)), 

115 (1250.0, geom.Point2D(cls.x0 + 19.3, cls.y0 + 11.2), Quadrupole(1.5, 2.5, 0.75)), 

116 (2100.0, geom.Point2D(cls.x0 + 5.1, cls.y0 + 21.2), Quadrupole(1.7, 1.9, 0.05)), 

117 (2800.0, geom.Point2D(cls.x0 + 24.1, cls.y0 + 19.2), Quadrupole(1.9, 1.7, 0.1)), 

118 (2350.0, geom.Point2D(cls.x0 + 40.3, cls.y0 + 13.9), Quadrupole(1.8, 1.8, -0.4)), 

119 (4999.0, geom.Point2D(cls.x0 + 45.8, cls.y0 + 22.0), Quadrupole(1.6, 1.2, 0.2)), 

120 ) 

121 

122 # The test points are chosen to cover various corner cases assuming 

123 # inner_size = (17, 15) and border_size = 5. If that is changed, the 

124 # test points should be updated to not fall outside the coadd and still 

125 # cover the description in the inline comments. This is checked below. 

126 test_points = ( 

127 geom.Point2D(cls.x0 + 5, cls.y0 + 4), # inner point in lower left 

128 geom.Point2D(cls.x0 + 6, cls.y0 + 24), # inner point in upper left 

129 geom.Point2D(cls.x0 + 25.2, cls.y0 + 7.8), # inner point in lower middle 

130 geom.Point2D(cls.x0 + 23, cls.y0 + 22), # inner point in upper middle 

131 geom.Point2D(cls.x0 + 39, cls.y0 + 9.4), # inner point in lower right 

132 geom.Point2D(cls.x0 + 44, cls.y0 + 24), # inner point in upper right 

133 # Some points that lie on the border 

134 geom.Point2D(cls.x0 + 33, cls.y0 + 24), # inner point in upper right 

135 geom.Point2D(cls.x0 + 46, cls.y0 + 0), # inner point in lower right 

136 geom.Point2D(cls.x0 + 19, cls.y0 + 16), # inner point in upper middle 

137 geom.Point2D(cls.x0 + 17, cls.y0 + 8), # inner point in lower middle 

138 geom.Point2D(cls.x0 + 0, cls.y0 + 29), # inner point in upper left 

139 geom.Point2D(cls.x0 + 0, cls.y0 + 0), # inner point in lower left 

140 ) 

141 # A list of (point, cell_index) pairs. 

142 cls.test_positions = [ 

143 ( 

144 point, 

145 Index2D( 

146 x=int((point.getX() - cls.x0) // cls.inner_size_x), 

147 y=int((point.getY() - cls.y0) // cls.inner_size_y), 

148 ), 

149 ) 

150 for point in test_points 

151 ] 

152 

153 schema = lsst.meas.base.tests.TestDataset.makeMinimalSchema() 

154 

155 single_cell_coadds = [] 

156 cls.exposures = dict.fromkeys(cls.psf_sigmas.keys()) 

157 

158 for x in range(cls.nx): 

159 for y in range(cls.ny): 

160 identifiers = CellIdentifiers( 

161 cell=Index2D(x=x, y=y), 

162 skymap=common.identifiers.skymap, 

163 tract=common.identifiers.tract, 

164 patch=common.identifiers.patch, 

165 band=common.identifiers.band, 

166 ) 

167 

168 outer_bbox = geom.Box2I( 

169 geom.Point2I(cls.x0 + x * cls.inner_size_x, cls.y0 + y * cls.inner_size_y), 

170 geom.Extent2I(cls.inner_size_x, cls.inner_size_y), 

171 ) 

172 outer_bbox.grow(cls.border_size) 

173 

174 dataset = lsst.meas.base.tests.TestDataset( 

175 patch_outer_bbox, psfSigma=cls.psf_sigmas[identifiers.cell] 

176 ) 

177 

178 for inst_flux, position, shape in sources: 

179 dataset.addSource(inst_flux, position, shape) 

180 

181 # Create a spatially varying variance plane. 

182 variance = ImageF( 

183 # np.random.uniform returns an array with x-y flipped. 

184 np.random.uniform( 

185 0.8, 

186 1.2, 

187 ( 

188 cls.ny * cls.inner_size_y + 2 * cls.border_size, 

189 cls.nx * cls.inner_size_x + 2 * cls.border_size, 

190 ), 

191 ).astype(np.float32), 

192 xy0=outer_bbox.getMin(), 

193 ) 

194 noise_realizations = [ 

195 ImageF( 

196 np.random.normal( 

197 0.0, 

198 1.0, 

199 ( 

200 cls.ny * cls.inner_size_y + 2 * cls.border_size, 

201 cls.nx * cls.inner_size_x + 2 * cls.border_size, 

202 ), 

203 ).astype(np.float32) 

204 * variance.array, 

205 xy0=outer_bbox.getMin(), 

206 ) 

207 for _ in range(cls.n_noise_realizations) 

208 ] 

209 mask_fraction = ImageF( 

210 np.random.uniform( 

211 0.0, 

212 1.0, 

213 ( 

214 cls.ny * cls.inner_size_y + 2 * cls.border_size, 

215 cls.nx * cls.inner_size_x + 2 * cls.border_size, 

216 ), 

217 ).astype(np.float32), 

218 xy0=outer_bbox.getMin(), 

219 ) 

220 exposure, _ = dataset.realize(variance.getArray() ** 0.5, schema, randomSeed=123456789) 

221 cls.exposures[identifiers.cell] = exposure 

222 exposure = exposure[outer_bbox] 

223 image_plane = OwnedImagePlanes( 

224 image=exposure.image, 

225 variance=exposure.variance, 

226 mask=exposure.mask, 

227 noise_realizations=noise_realizations, 

228 mask_fractions=mask_fraction, 

229 ) 

230 aperture_correction_map = frozendict( 

231 base_GaussianFlux_instFlux=0.9 * (x + 1), 

232 base_GaussianFlux_instFluxErr=0.01 * (y + 1), 

233 base_PsfFlux_instFlux=0.8 * (y + 2), 

234 base_PsfFlux_instFluxErr=0.02 * (x + 2), 

235 ) 

236 

237 single_cell_coadds.append( 

238 SingleCellCoadd( 

239 outer=image_plane, 

240 psf=SingleGaussianPsf( 

241 cls.psf_size_x, cls.psf_size_y, cls.psf_sigmas[Index2D(x=x, y=y)] 

242 ).computeKernelImage(outer_bbox.getCenter()), 

243 inner_bbox=geom.Box2I( 

244 geom.Point2I(cls.x0 + x * cls.inner_size_x, cls.y0 + y * cls.inner_size_y), 

245 geom.Extent2I(cls.inner_size_x, cls.inner_size_y), 

246 ), 

247 inputs=( 

248 ObservationIdentifiers( 

249 instrument="dummy", 

250 physical_filter="dummy-I", 

251 visit=12345, 

252 detector=67, 

253 day_obs=20000101, 

254 ), 

255 ), 

256 common=common, 

257 identifiers=identifiers, 

258 aperture_correction_map=aperture_correction_map, 

259 ) 

260 ) 

261 

262 grid_bbox = geom.Box2I( 

263 geom.Point2I(cls.x0, cls.y0), geom.Extent2I(cls.nx * cls.inner_size_x, cls.ny * cls.inner_size_y) 

264 ) 

265 grid = UniformGrid.from_bbox_shape(grid_bbox, Index2D(x=cls.nx, y=cls.ny), padding=cls.border_size) 

266 

267 for test_point, test_index in cls.test_positions: 

268 assert test_index == grid.index( 

269 geom.Point2I(test_point) 

270 ), f"Test point {test_point} is not in cell {test_index}." 

271 

272 cls.multiple_cell_coadd = MultipleCellCoadd( 

273 single_cell_coadds, 

274 grid=grid, 

275 outer_cell_size=geom.Extent2I(cls.outer_size_x, cls.outer_size_y), 

276 inner_bbox=None, 

277 common=common, 

278 psf_image_size=geom.Extent2I(cls.psf_size_x, cls.psf_size_y), 

279 ) 

280 

281 @classmethod 

282 def tearDownClass(cls) -> None: # noqa: D102 

283 # Docstring inherited 

284 del cls.multiple_cell_coadd 

285 del cls.exposures 

286 super().tearDownClass() 

287 

288 def assertMultipleCellCoaddsEqual(self, mcc1: MultipleCellCoadd, mcc2: MultipleCellCoadd) -> None: 

289 """Check the equality of two instances of `MultipleCellCoadd`. 

290 

291 Parameters 

292 ---------- 

293 mcc1 : `MultipleCellCoadd` 

294 The MultipleCellCoadd created by reading a FITS file. 

295 mcc2 : `MultipleCellCoadd` 

296 The reference MultipleCellCoadd for comparison. 

297 """ 

298 self.assertEqual(mcc1.band, mcc2.band) 

299 self.assertEqual(mcc1.identifiers, mcc2.identifiers) 

300 self.assertEqual(mcc1.inner_bbox, mcc2.inner_bbox) 

301 self.assertEqual(mcc1.outer_bbox, mcc2.outer_bbox) 

302 self.assertEqual(mcc1.outer_cell_size, mcc2.outer_cell_size) 

303 self.assertEqual(mcc1.grid, mcc2.grid) 

304 self.assertEqual(mcc1.mask_fraction_names, mcc2.mask_fraction_names) 

305 self.assertEqual(mcc1.n_noise_realizations, mcc2.n_noise_realizations) 

306 self.assertEqual(mcc1.psf_image_size, mcc2.psf_image_size) 

307 self.assertEqual(mcc1.units, mcc2.units) 

308 self.assertEqual(mcc1.wcs.getFitsMetadata().toString(), mcc2.wcs.getFitsMetadata().toString()) 

309 

310 # Check that the individual cells are identical. 

311 self.assertEqual(mcc1.cells.keys(), mcc2.cells.keys()) 

312 for idx in mcc1.cells.keys(): # noqa: SIM118 

313 self.assertImagesEqual(mcc1.cells[idx].outer.image, mcc2.cells[idx].outer.image) 

314 self.assertMasksEqual(mcc1.cells[idx].outer.mask, mcc2.cells[idx].outer.mask) 

315 self.assertImagesEqual(mcc1.cells[idx].outer.variance, mcc2.cells[idx].outer.variance) 

316 self.assertImagesEqual(mcc1.cells[idx].outer.mask_fractions, mcc2.cells[idx].outer.mask_fractions) 

317 self.assertEqual( 

318 len(mcc1.cells[idx].outer.noise_realizations), len(mcc2.cells[idx].outer.noise_realizations) 

319 ) 

320 for noise_id in range(len(mcc1.cells[idx].outer.noise_realizations)): 

321 with self.subTest(noise_id): 

322 self.assertImagesEqual( 

323 mcc1.cells[idx].outer.noise_realizations[noise_id], 

324 mcc2.cells[idx].outer.noise_realizations[noise_id], 

325 ) 

326 self.assertImagesEqual(mcc1.cells[idx].psf_image, mcc2.cells[idx].psf_image) 

327 self.assertEqual(mcc1.cells[idx].inputs, mcc2.cells[idx].inputs) 

328 

329 self.assertEqual(mcc1.cells[idx].band, mcc1.band) 

330 self.assertEqual(mcc1.cells[idx].common, mcc1.common) 

331 self.assertEqual(mcc1.cells[idx].units, mcc2.units) 

332 self.assertEqual(mcc1.cells[idx].wcs, mcc1.wcs) 

333 self.assertEqual(mcc1.cells[idx].aperture_correction_map, mcc2.cells[idx].aperture_correction_map) 

334 

335 # Identifiers differ because of the ``cell`` component. 

336 # Check the other attributes within the identifiers. 

337 for attr in ("skymap", "tract", "patch", "band"): 

338 self.assertEqual(getattr(mcc1.cells[idx].identifiers, attr), getattr(mcc1.identifiers, attr)) 

339 

340 

341class MultipleCellCoaddTestCase(BaseMultipleCellCoaddTestCase): 

342 """Test the construction and interfaces of MultipleCellCoadd.""" 

343 

344 def test_fits(self): 

345 """Test that we can write a coadd to a FITS file and read it.""" 

346 with lsst.utils.tests.getTempFilePath(".fits") as filename: 

347 self.multiple_cell_coadd.write_fits(filename) 

348 mcc1 = MultipleCellCoadd.read_fits(filename) # Test the readFits method. 

349 

350 # Test the reader class. 

351 reader = CellCoaddFitsReader(filename) 

352 mcc2 = reader.readAsMultipleCellCoadd() 

353 

354 wcs = reader.readWcs() 

355 

356 self.assertMultipleCellCoaddsEqual(mcc1, self.multiple_cell_coadd) 

357 self.assertMultipleCellCoaddsEqual(mcc2, self.multiple_cell_coadd) 

358 # By transititve property of equality, mcc1 == mcc2. 

359 

360 self.assertEqual(self.multiple_cell_coadd.band, self.multiple_cell_coadd.common.band) 

361 self.assertEqual( 

362 wcs.getFitsMetadata().toString(), self.multiple_cell_coadd.wcs.getFitsMetadata().toString() 

363 ) 

364 

365 def test_visit_count(self): 

366 """Test the visit_count method.""" 

367 # Since we don't simulate coaddition from multiple warps, the cells are 

368 # all going to have just a single visit. 

369 for cellId, singleCellCoadd in self.multiple_cell_coadd.cells.items(): 

370 with self.subTest(x=cellId.x, y=cellId.y): 

371 self.assertEqual(singleCellCoadd.visit_count, 1) 

372 

373 def test_aperture_correction(self): 

374 """Test the aperture correction values are what we expect.""" 

375 for cellId, single_cell_coadd in self.multiple_cell_coadd.cells.items(): 

376 with self.subTest(x=cellId.x, y=cellId.y): 

377 ap_corr_map = single_cell_coadd.aperture_correction_map 

378 self.assertIsNotNone(ap_corr_map) 

379 self.assertEqual(ap_corr_map["base_GaussianFlux_instFlux"], 0.9 * (cellId.x + 1)) 

380 self.assertEqual(ap_corr_map["base_GaussianFlux_instFluxErr"], 0.01 * (cellId.y + 1)) 

381 self.assertEqual(ap_corr_map["base_PsfFlux_instFlux"], 0.8 * (cellId.y + 2)) 

382 self.assertEqual(ap_corr_map["base_PsfFlux_instFluxErr"], 0.02 * (cellId.x + 2)) 

383 

384 

385class ExplodedCoaddTestCase(BaseMultipleCellCoaddTestCase): 

386 """Test the construction and methods of an ExplodedCoadd instance.""" 

387 

388 exploded_coadd: ExplodedCoadd 

389 

390 @classmethod 

391 def setUpClass(cls) -> None: # noqa: D102 

392 # Docstring inherited 

393 super().setUpClass() 

394 cls.exploded_coadd = cls.multiple_cell_coadd.explode() 

395 

396 @classmethod 

397 def tearDownClass(cls) -> None: # noqa: D102 

398 # Docstring inherited 

399 del cls.exploded_coadd 

400 super().tearDownClass() 

401 

402 def test_exploded_psf_image(self): 

403 """Show that psf_image sizes are absurd.""" 

404 self.assertEqual( 

405 self.exploded_coadd.psf_image.getBBox().getDimensions(), 

406 geom.Extent2I(self.nx * self.psf_size_x, self.ny * self.psf_size_y), 

407 ) 

408 for pad_psfs_with in (-999, -4, 0, 4, 8, 21, 40, 100): 

409 exploded_coadd = self.multiple_cell_coadd.explode(pad_psfs_with=pad_psfs_with) 

410 self.assertEqual( 

411 exploded_coadd.psf_image.getBBox().getDimensions(), 

412 geom.Extent2I(self.nx * self.outer_size_x, self.ny * self.outer_size_y), 

413 ) 

414 

415 def test_asMaskedImage(self): 

416 """Test the asMaskedImage method for an ExplodedCoadd object.""" 

417 masked_image = self.exploded_coadd.asMaskedImage() 

418 masked_image.setXY0(self.multiple_cell_coadd.outer_bbox.getMin()) 

419 base_bbox = self.multiple_cell_coadd.cells.first.outer.bbox 

420 for cell_x, cell_y in product(range(self.nx), range(self.ny)): 

421 bbox = base_bbox.shiftedBy(geom.Extent2I(cell_x * self.outer_size_x, cell_y * self.outer_size_y)) 

422 with self.subTest(cell_x=cell_x, cell_y=cell_y): 

423 self.assertMaskedImagesEqual( 

424 masked_image[bbox], 

425 self.multiple_cell_coadd.cells[Index2D(cell_x, cell_y)].outer.asMaskedImage(), 

426 ) 

427 

428 

429class StitchedCoaddTestCase(BaseMultipleCellCoaddTestCase): 

430 """Test the construction and methods of a StitchedCoadd instance.""" 

431 

432 stitched_coadd: StitchedCoadd 

433 

434 @classmethod 

435 def setUpClass(cls) -> None: # noqa: D102 

436 # Docstring inherited 

437 super().setUpClass() 

438 cls.stitched_coadd = cls.multiple_cell_coadd.stitch() 

439 

440 @classmethod 

441 def tearDownClass(cls) -> None: # noqa: D102 

442 # Docstring inherited 

443 del cls.stitched_coadd 

444 super().tearDownClass() 

445 

446 def test_computeBBox(self): 

447 """Test the computeBBox method for a StitchedPsf object.""" 

448 stitched_psf = self.stitched_coadd.psf 

449 

450 psf_bbox = geom.Box2I( 

451 geom.Point2I(-(self.psf_size_x // 2), -(self.psf_size_y // 2)), 

452 geom.Extent2I(self.psf_size_x, self.psf_size_y), 

453 ) 

454 

455 for position, _ in self.test_positions: 

456 bbox = stitched_psf.computeBBox(position) 

457 self.assertEqual(bbox, psf_bbox) 

458 

459 def test_computeShape(self): 

460 """Test the computeShape method for a StitchedPsf object.""" 

461 stitched_psf = self.stitched_coadd.psf 

462 for position, cell_index in self.test_positions: 

463 psf_shape = stitched_psf.computeShape(position) # check we can compute shape 

464 self.assertIsNot(psf_shape.getIxx(), np.nan) 

465 self.assertIsNot(psf_shape.getIyy(), np.nan) 

466 self.assertIsNot(psf_shape.getIxy(), np.nan) 

467 

468 # Moments measured from pixellated images are significantly 

469 # underestimated for small PSFs. 

470 if self.psf_sigmas[cell_index] >= 1.0: 

471 self.assertAlmostEqual(psf_shape.getIxx(), self.psf_sigmas[cell_index] ** 2, delta=1e-3) 

472 self.assertAlmostEqual(psf_shape.getIyy(), self.psf_sigmas[cell_index] ** 2, delta=1e-3) 

473 self.assertAlmostEqual(psf_shape.getIxy(), 0.0) 

474 

475 def test_computeKernelImage(self): 

476 """Test the computeKernelImage method for a StitchedPsf object.""" 

477 stitched_psf = self.stitched_coadd.psf 

478 psf_bbox = geom.Box2I( 

479 geom.Point2I(-(self.psf_size_x // 2), -(self.psf_size_y // 2)), 

480 geom.Extent2I(self.psf_size_x, self.psf_size_y), 

481 ) 

482 

483 for position, cell_index in self.test_positions: 

484 image1 = stitched_psf.computeKernelImage(position) 

485 image2 = SingleGaussianPsf( 

486 self.psf_size_x, self.psf_size_y, self.psf_sigmas[cell_index] 

487 ).computeKernelImage(position) 

488 self.assertImagesEqual(image1, image2) 

489 self.assertEqual(image1.getBBox(), psf_bbox) 

490 

491 def test_computeImage(self): 

492 """Test the computeImage method for a StitchedPsf object.""" 

493 stitched_psf = self.stitched_coadd.psf 

494 psf_extent = geom.Extent2I(self.psf_size_x, self.psf_size_y) 

495 

496 for position, cell_index in self.test_positions: 

497 image1 = stitched_psf.computeImage(position) 

498 image2 = SingleGaussianPsf( 

499 self.psf_size_x, self.psf_size_y, self.psf_sigmas[cell_index] 

500 ).computeImage(position) 

501 self.assertImagesEqual(image1, image2) 

502 self.assertEqual(image1.getBBox().getDimensions(), psf_extent) 

503 

504 def test_computeImage_computeKernelImage(self): 

505 """Test that computeImage called at integer points gives the same 

506 result as calling computeKernelImage. 

507 """ 

508 stitched_psf = self.stitched_coadd.psf 

509 for position, _cell_index in self.test_positions: 

510 pos = geom.Point2D(geom.Point2I(position)) # round to integer 

511 image1 = stitched_psf.computeKernelImage(pos) 

512 image2 = stitched_psf.computeImage(pos) 

513 self.assertImagesEqual(image1, image2) 

514 

515 def test_computeApetureFlux(self): 

516 """Test the computeApertureFlux method for a StitchedPsf object.""" 

517 stitched_psf = self.stitched_coadd.psf 

518 for position, cell_index in self.test_positions: 

519 flux1sigma = stitched_psf.computeApertureFlux(self.psf_sigmas[cell_index], position=position) 

520 self.assertAlmostEqual(flux1sigma, 0.39, delta=5e-2) 

521 

522 flux3sigma = stitched_psf.computeApertureFlux( 

523 3.0 * self.psf_sigmas[cell_index], position=position 

524 ) 

525 self.assertAlmostEqual(flux3sigma, 0.97, delta=2e-2) 

526 

527 def test_asExposure(self): 

528 """Test the asExposure method for a StitchedCoadd object.""" 

529 exposure = self.stitched_coadd.asExposure() 

530 

531 # Check that the bounding box is correct. 

532 bbox = exposure.getBBox() 

533 self.assertEqual(bbox.getWidth(), self.inner_size_x * self.nx + 2 * self.border_size) 

534 self.assertEqual(bbox.getHeight(), self.inner_size_y * self.ny + 2 * self.border_size) 

535 

536 for y in range(self.ny): 

537 for x in range(self.nx): 

538 bbox = geom.Box2I( 

539 geom.Point2I(self.x0 + x * self.inner_size_x, self.y0 + y * self.inner_size_y), 

540 geom.Extent2I(self.inner_size_x, self.inner_size_y), 

541 ) 

542 index = Index2D(x=x, y=y) 

543 self.assertImagesEqual(exposure.image[bbox], self.exposures[index].image[bbox]) 

544 self.assertImagesEqual(exposure.variance[bbox], self.exposures[index].variance[bbox]) 

545 self.assertImagesEqual(exposure.mask[bbox], self.exposures[index].mask[bbox]) 

546 

547 def test_aperture_correction(self): 

548 """Test the aperture correction values are what we expect.""" 

549 ap_corr_map = self.stitched_coadd.ap_corr_map 

550 algorithm_names = ap_corr_map.keys() 

551 

552 for (position, cellId), algorithm_name in product(self.test_positions, algorithm_names): 

553 with self.subTest(x=cellId.x, y=cellId.y, algorithm_name=algorithm_name): 

554 field_name = algorithm_name 

555 self.assertEqual( 

556 ap_corr_map[field_name].evaluate(position), 

557 self.multiple_cell_coadd.cells[cellId].aperture_correction_map[field_name], 

558 ) 

559 

560 def test_inputs(self): 

561 """Test that the inputs are populated correctly on stitching.""" 

562 inputs = self.stitched_coadd.ccds 

563 for position, cellId in self.test_positions: 

564 with self.subTest(x=cellId.x, y=cellId.y): 

565 self.assertEqual( 

566 inputs.subsetContaining(position), 

567 self.multiple_cell_coadd.cells[cellId].inputs, 

568 ) 

569 self.assertEqual(len(self.stitched_coadd.visits), 1) 

570 

571 def test_borders(self): 

572 """Test that the borders are populated correctly on stitching.""" 

573 mi = self.stitched_coadd.asMaskedImage() 

574 

575 # Check that the borders are not all zero. 

576 with self.assertRaises(AssertionError): 

577 np.testing.assert_array_equal(mi.image.array[: self.border_size, :], 0) 

578 with self.assertRaises(AssertionError): 

579 np.testing.assert_array_equal(mi.image.array[-self.border_size :, :], 0) 

580 with self.assertRaises(AssertionError): 

581 np.testing.assert_array_equal(mi.image.array[:, : self.border_size], 0) 

582 with self.assertRaises(AssertionError): 

583 np.testing.assert_array_equal(mi.image.array[:, -self.border_size :], 0) 

584 

585 def test_fits(self): 

586 """Test that we can write an Exposure with StitchedPsf to a FITS file 

587 and read it. 

588 """ 

589 write_exposure = self.stitched_coadd.asExposure() 

590 with lsst.utils.tests.getTempFilePath(".fits") as filename: 

591 write_exposure.writeFits(filename) 

592 read_exposure = ExposureF.readFits(filename) # Test the readFits method. 

593 

594 # Test that the image planes are identical. 

595 self.assertImagesEqual(read_exposure.image, write_exposure.image) 

596 self.assertImagesEqual(read_exposure.variance, write_exposure.variance) 

597 self.assertImagesEqual(read_exposure.mask, write_exposure.mask) 

598 

599 # Test the PSF images in the StitchedPsf. 

600 for index in write_exposure.psf.images.indices(): 

601 self.assertImagesEqual(read_exposure.psf.images[index], write_exposure.psf.images[index]) 

602 

603 # Test that the WCSs are equal. 

604 self.assertEqual( 

605 read_exposure.wcs.getFitsMetadata().toString(), 

606 write_exposure.wcs.getFitsMetadata().toString(), 

607 ) 

608 

609 

610class TestMemory(lsst.utils.tests.MemoryTestCase): 

611 """Check for resource/memory leaks.""" 

612 

613 

614def setup_module(module): # noqa: D103 

615 lsst.utils.tests.init() 

616 

617 

618if __name__ == "__main__": 618 ↛ 619line 618 didn't jump to line 619 because the condition on line 618 was never true

619 lsst.utils.tests.init() 

620 unittest.main()