Coverage for tests/test_coadds.py: 19%
270 statements
« prev ^ index » next coverage.py v7.11.0, created at 2025-11-07 02:32 -0800
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/>.
22import unittest
23from collections.abc import Iterable, Mapping
24from itertools import product
26import numpy as np
27from frozendict import frozendict
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 Polygon, Quadrupole
34from lsst.afw.image import ExposureF, ImageF
35from lsst.cell_coadds import (
36 CellCoaddFitsReader,
37 CellIdentifiers,
38 CoaddInputs,
39 CoaddUnits,
40 CommonComponents,
41 ExplodedCoadd,
42 MultipleCellCoadd,
43 ObservationIdentifiers,
44 OwnedImagePlanes,
45 PatchIdentifiers,
46 SingleCellCoadd,
47 StitchedCoadd,
48 UniformGrid,
49)
50from lsst.meas.algorithms import SingleGaussianPsf
51from lsst.skymap import Index2D
54class BaseMultipleCellCoaddTestCase(lsst.utils.tests.TestCase):
55 """A base class that provides a common set of methods."""
57 psf_size: int
58 psf_sigmas: Mapping[Index2D, float]
59 border_size: int
60 inner_size: int
61 outer_size: int
62 test_positions: Iterable[tuple[geom.Point2D, Index2D]]
63 exposures: Mapping[Index2D, ExposureF]
64 multiple_cell_coadd: MultipleCellCoadd
66 @classmethod
67 def setUpClass(cls) -> None:
68 """Set up a multiple cell coadd with 2x2 cells."""
69 np.random.seed(42)
70 data_id = test_utils.generate_data_id()
71 common = CommonComponents(
72 units=CoaddUnits.nJy, # units here are arbitrary.
73 wcs=test_utils.generate_wcs(),
74 band=data_id["band"],
75 identifiers=PatchIdentifiers.from_data_id(data_id),
76 visit_polygons={ # This is arbitrary, just to check it if it goes through FITS persistence.
77 ObservationIdentifiers(
78 instrument="dummy",
79 physical_filter="dummy-I",
80 visit=12345,
81 detector=67,
82 day_obs=20000101,
83 ): Polygon([geom.Point2D(0, 0), geom.Point2D(1, 0), geom.Point2D(1, 1), geom.Point2D(0, 1)])
84 },
85 )
87 cls.nx, cls.ny = 3, 2
88 cls.psf_sigmas = {
89 Index2D(x=0, y=0): 1.2,
90 Index2D(x=0, y=1): 0.7,
91 Index2D(x=1, y=0): 0.9,
92 Index2D(x=1, y=1): 1.1,
93 Index2D(x=2, y=0): 1.3,
94 Index2D(x=2, y=1): 0.8,
95 }
97 cls.border_size = 5
98 # In practice, we expect this to be squares.
99 # To check for any possible x, y swap, we use different values.
100 cls.psf_size_x, cls.psf_size_y = 21, 23
101 cls.inner_size_x, cls.inner_size_y = 17, 15
102 cls.outer_size_x = cls.inner_size_x + 2 * cls.border_size
103 cls.outer_size_y = cls.inner_size_y + 2 * cls.border_size
104 # The origin should not be at (0, 0) for robust testing.
105 cls.x0, cls.y0 = 5, 2
106 cls.n_noise_realizations = 2
108 patch_outer_bbox = geom.Box2I(
109 geom.Point2I(cls.x0, cls.y0), geom.Extent2I(cls.nx * cls.inner_size_x, cls.ny * cls.inner_size_y)
110 )
111 patch_outer_bbox.grow(cls.border_size)
113 # Add one star and one galaxy per quadrant.
114 # The mapping of positions to cell indices assume inner_size = (17, 15)
115 # and border_size = 5. If that is changed, these values need an update.
116 sources = (
117 # flux, centroid, shape
118 (1000.0, geom.Point2D(cls.x0 + 6.3, cls.y0 + 7.2), None),
119 (2500.0, geom.Point2D(cls.x0 + 16.8, cls.y0 + 18.3), None),
120 (1500.0, geom.Point2D(cls.x0 + 21.2, cls.y0 + 5.1), None),
121 (3200.0, geom.Point2D(cls.x0 + 16.1, cls.y0 + 23.9), None),
122 (1800.0, geom.Point2D(cls.x0 + 44.7, cls.y0 + 8.9), None),
123 (2100.0, geom.Point2D(cls.x0 + 34.1, cls.y0 + 19.2), None),
124 (900.0, geom.Point2D(cls.x0 + 9.1, cls.y0 + 13.9), Quadrupole(2.5, 1.5, 0.8)),
125 (1250.0, geom.Point2D(cls.x0 + 19.3, cls.y0 + 11.2), Quadrupole(1.5, 2.5, 0.75)),
126 (2100.0, geom.Point2D(cls.x0 + 5.1, cls.y0 + 21.2), Quadrupole(1.7, 1.9, 0.05)),
127 (2800.0, geom.Point2D(cls.x0 + 24.1, cls.y0 + 19.2), Quadrupole(1.9, 1.7, 0.1)),
128 (2350.0, geom.Point2D(cls.x0 + 40.3, cls.y0 + 13.9), Quadrupole(1.8, 1.8, -0.4)),
129 (4999.0, geom.Point2D(cls.x0 + 45.8, cls.y0 + 22.0), Quadrupole(1.6, 1.2, 0.2)),
130 )
132 # The test points are chosen to cover various corner cases assuming
133 # inner_size = (17, 15) and border_size = 5. If that is changed, the
134 # test points should be updated to not fall outside the coadd and still
135 # cover the description in the inline comments. This is checked below.
136 test_points = (
137 geom.Point2D(cls.x0 + 5, cls.y0 + 4), # inner point in lower left
138 geom.Point2D(cls.x0 + 6, cls.y0 + 24), # inner point in upper left
139 geom.Point2D(cls.x0 + 25.2, cls.y0 + 7.8), # inner point in lower middle
140 geom.Point2D(cls.x0 + 23, cls.y0 + 22), # inner point in upper middle
141 geom.Point2D(cls.x0 + 39, cls.y0 + 9.4), # inner point in lower right
142 geom.Point2D(cls.x0 + 44, cls.y0 + 24), # inner point in upper right
143 # Some points that lie on the border
144 geom.Point2D(cls.x0 + 33, cls.y0 + 24), # inner point in upper right
145 geom.Point2D(cls.x0 + 46, cls.y0 + 0), # inner point in lower right
146 geom.Point2D(cls.x0 + 19, cls.y0 + 16), # inner point in upper middle
147 geom.Point2D(cls.x0 + 17, cls.y0 + 8), # inner point in lower middle
148 geom.Point2D(cls.x0 + 0, cls.y0 + 29), # inner point in upper left
149 geom.Point2D(cls.x0 + 0, cls.y0 + 0), # inner point in lower left
150 )
151 # A list of (point, cell_index) pairs.
152 cls.test_positions = [
153 (
154 point,
155 Index2D(
156 x=int((point.getX() - cls.x0) // cls.inner_size_x),
157 y=int((point.getY() - cls.y0) // cls.inner_size_y),
158 ),
159 )
160 for point in test_points
161 ]
163 schema = lsst.meas.base.tests.TestDataset.makeMinimalSchema()
165 single_cell_coadds = []
166 cls.exposures = dict.fromkeys(cls.psf_sigmas.keys())
168 for x in range(cls.nx):
169 for y in range(cls.ny):
170 identifiers = CellIdentifiers(
171 cell=Index2D(x=x, y=y),
172 skymap=common.identifiers.skymap,
173 tract=common.identifiers.tract,
174 patch=common.identifiers.patch,
175 band=common.identifiers.band,
176 )
178 outer_bbox = geom.Box2I(
179 geom.Point2I(cls.x0 + x * cls.inner_size_x, cls.y0 + y * cls.inner_size_y),
180 geom.Extent2I(cls.inner_size_x, cls.inner_size_y),
181 )
182 outer_bbox.grow(cls.border_size)
184 dataset = lsst.meas.base.tests.TestDataset(
185 patch_outer_bbox, psfSigma=cls.psf_sigmas[identifiers.cell]
186 )
188 for inst_flux, position, shape in sources:
189 dataset.addSource(inst_flux, position, shape)
191 # Create a spatially varying variance plane.
192 variance = ImageF(
193 # np.random.uniform returns an array with x-y flipped.
194 np.random.uniform(
195 0.8,
196 1.2,
197 (
198 cls.ny * cls.inner_size_y + 2 * cls.border_size,
199 cls.nx * cls.inner_size_x + 2 * cls.border_size,
200 ),
201 ).astype(np.float32),
202 xy0=outer_bbox.getMin(),
203 )
204 noise_realizations = [
205 ImageF(
206 np.random.normal(
207 0.0,
208 1.0,
209 (
210 cls.ny * cls.inner_size_y + 2 * cls.border_size,
211 cls.nx * cls.inner_size_x + 2 * cls.border_size,
212 ),
213 ).astype(np.float32)
214 * variance.array,
215 xy0=outer_bbox.getMin(),
216 )
217 for _ in range(cls.n_noise_realizations)
218 ]
219 mask_fraction = ImageF(
220 np.random.uniform(
221 0.0,
222 1.0,
223 (
224 cls.ny * cls.inner_size_y + 2 * cls.border_size,
225 cls.nx * cls.inner_size_x + 2 * cls.border_size,
226 ),
227 ).astype(np.float32),
228 xy0=outer_bbox.getMin(),
229 )
230 exposure, _ = dataset.realize(variance.getArray() ** 0.5, schema, randomSeed=123456789)
231 cls.exposures[identifiers.cell] = exposure
232 exposure = exposure[outer_bbox]
233 image_plane = OwnedImagePlanes(
234 image=exposure.image,
235 variance=exposure.variance,
236 mask=exposure.mask,
237 noise_realizations=noise_realizations,
238 mask_fractions=mask_fraction,
239 )
240 aperture_correction_map = frozendict(
241 base_GaussianFlux_instFlux=0.9 * (x + 1),
242 base_GaussianFlux_instFluxErr=0.01 * (y + 1),
243 base_PsfFlux_instFlux=0.8 * (y + 2),
244 base_PsfFlux_instFluxErr=0.02 * (x + 2),
245 )
247 single_cell_coadds.append(
248 SingleCellCoadd(
249 outer=image_plane,
250 psf=SingleGaussianPsf(
251 cls.psf_size_x, cls.psf_size_y, cls.psf_sigmas[Index2D(x=x, y=y)]
252 ).computeKernelImage(outer_bbox.getCenter()),
253 inner_bbox=geom.Box2I(
254 geom.Point2I(cls.x0 + x * cls.inner_size_x, cls.y0 + y * cls.inner_size_y),
255 geom.Extent2I(cls.inner_size_x, cls.inner_size_y),
256 ),
257 inputs={
258 ObservationIdentifiers(
259 instrument="dummy",
260 physical_filter="dummy-I",
261 visit=12345,
262 detector=67,
263 day_obs=20000101,
264 ): CoaddInputs(
265 True,
266 1.0,
267 1.0,
268 Quadrupole(
269 cls.psf_sigmas[Index2D(x=x, y=y)] ** 2,
270 cls.psf_sigmas[Index2D(x=x, y=y)] ** 2,
271 ),
272 False,
273 )
274 },
275 common=common,
276 identifiers=identifiers,
277 aperture_correction_map=aperture_correction_map,
278 )
279 )
281 grid_bbox = geom.Box2I(
282 geom.Point2I(cls.x0, cls.y0), geom.Extent2I(cls.nx * cls.inner_size_x, cls.ny * cls.inner_size_y)
283 )
284 grid = UniformGrid.from_bbox_shape(grid_bbox, Index2D(x=cls.nx, y=cls.ny), padding=cls.border_size)
286 for test_point, test_index in cls.test_positions:
287 assert test_index == grid.index(
288 geom.Point2I(test_point)
289 ), f"Test point {test_point} is not in cell {test_index}."
291 cls.multiple_cell_coadd = MultipleCellCoadd(
292 single_cell_coadds,
293 grid=grid,
294 outer_cell_size=geom.Extent2I(cls.outer_size_x, cls.outer_size_y),
295 inner_bbox=None,
296 common=common,
297 psf_image_size=geom.Extent2I(cls.psf_size_x, cls.psf_size_y),
298 )
300 @classmethod
301 def tearDownClass(cls) -> None: # noqa: D102
302 # Docstring inherited
303 del cls.multiple_cell_coadd
304 del cls.exposures
305 super().tearDownClass()
307 def assertMultipleCellCoaddsEqual(self, mcc1: MultipleCellCoadd, mcc2: MultipleCellCoadd) -> None:
308 """Check the equality of two instances of `MultipleCellCoadd`.
310 Parameters
311 ----------
312 mcc1 : `MultipleCellCoadd`
313 The MultipleCellCoadd created by reading a FITS file.
314 mcc2 : `MultipleCellCoadd`
315 The reference MultipleCellCoadd for comparison.
316 """
317 self.assertEqual(mcc1.band, mcc2.band)
318 self.assertEqual(mcc1.identifiers, mcc2.identifiers)
319 self.assertEqual(mcc1.inner_bbox, mcc2.inner_bbox)
320 self.assertEqual(mcc1.outer_bbox, mcc2.outer_bbox)
321 self.assertEqual(mcc1.outer_cell_size, mcc2.outer_cell_size)
322 self.assertEqual(mcc1.grid, mcc2.grid)
323 self.assertEqual(mcc1.mask_fraction_names, mcc2.mask_fraction_names)
324 self.assertEqual(mcc1.n_noise_realizations, mcc2.n_noise_realizations)
325 self.assertEqual(mcc1.psf_image_size, mcc2.psf_image_size)
326 self.assertEqual(mcc1.units, mcc2.units)
327 self.assertEqual(mcc1.wcs.getFitsMetadata().toString(), mcc2.wcs.getFitsMetadata().toString())
329 # Check that the individual cells are identical.
330 self.assertEqual(mcc1.cells.keys(), mcc2.cells.keys())
331 for idx in mcc1.cells.keys(): # noqa: SIM118
332 self.assertImagesEqual(mcc1.cells[idx].outer.image, mcc2.cells[idx].outer.image)
333 self.assertMasksEqual(mcc1.cells[idx].outer.mask, mcc2.cells[idx].outer.mask)
334 self.assertImagesEqual(mcc1.cells[idx].outer.variance, mcc2.cells[idx].outer.variance)
335 self.assertImagesEqual(mcc1.cells[idx].outer.mask_fractions, mcc2.cells[idx].outer.mask_fractions)
336 self.assertEqual(
337 len(mcc1.cells[idx].outer.noise_realizations), len(mcc2.cells[idx].outer.noise_realizations)
338 )
339 for noise_id in range(len(mcc1.cells[idx].outer.noise_realizations)):
340 with self.subTest(noise_id):
341 self.assertImagesEqual(
342 mcc1.cells[idx].outer.noise_realizations[noise_id],
343 mcc2.cells[idx].outer.noise_realizations[noise_id],
344 )
345 self.assertImagesEqual(mcc1.cells[idx].psf_image, mcc2.cells[idx].psf_image)
346 self.assertEqual(mcc1.cells[idx].inputs, mcc2.cells[idx].inputs)
348 self.assertEqual(mcc1.cells[idx].band, mcc1.band)
349 self.assertEqual(mcc1.cells[idx].common, mcc1.common)
350 self.assertEqual(mcc1.cells[idx].units, mcc2.units)
351 self.assertEqual(mcc1.cells[idx].wcs, mcc1.wcs)
352 self.assertEqual(mcc1.cells[idx].aperture_correction_map, mcc2.cells[idx].aperture_correction_map)
354 # Identifiers differ because of the ``cell`` component.
355 # Check the other attributes within the identifiers.
356 for attr in ("skymap", "tract", "patch", "band"):
357 self.assertEqual(getattr(mcc1.cells[idx].identifiers, attr), getattr(mcc1.identifiers, attr))
360class MultipleCellCoaddTestCase(BaseMultipleCellCoaddTestCase):
361 """Test the construction and interfaces of MultipleCellCoadd."""
363 def test_fits(self):
364 """Test that we can write a coadd to a FITS file and read it."""
365 with lsst.utils.tests.getTempFilePath(".fits") as filename:
366 self.multiple_cell_coadd.write_fits(filename)
367 mcc1 = MultipleCellCoadd.read_fits(filename) # Test the readFits method.
369 # Test the reader class.
370 reader = CellCoaddFitsReader(filename)
371 mcc2 = reader.readAsMultipleCellCoadd()
373 wcs = reader.readWcs()
375 self.assertMultipleCellCoaddsEqual(mcc1, self.multiple_cell_coadd)
376 self.assertMultipleCellCoaddsEqual(mcc2, self.multiple_cell_coadd)
377 # By transititve property of equality, mcc1 == mcc2.
379 self.assertEqual(self.multiple_cell_coadd.band, self.multiple_cell_coadd.common.band)
380 for poly in self.multiple_cell_coadd.common.visit_polygons.values():
381 self.assertEqual(len(poly), 4)
382 self.assertEqual(
383 wcs.getFitsMetadata().toString(), self.multiple_cell_coadd.wcs.getFitsMetadata().toString()
384 )
386 def test_visit_count(self):
387 """Test the visit_count method."""
388 # Since we don't simulate coaddition from multiple warps, the cells are
389 # all going to have just a single visit.
390 for cellId, singleCellCoadd in self.multiple_cell_coadd.cells.items():
391 with self.subTest(x=cellId.x, y=cellId.y):
392 self.assertEqual(singleCellCoadd.visit_count, 1)
394 def test_aperture_correction(self):
395 """Test the aperture correction values are what we expect."""
396 for cellId, single_cell_coadd in self.multiple_cell_coadd.cells.items():
397 with self.subTest(x=cellId.x, y=cellId.y):
398 ap_corr_map = single_cell_coadd.aperture_correction_map
399 self.assertIsNotNone(ap_corr_map)
400 self.assertEqual(ap_corr_map["base_GaussianFlux_instFlux"], 0.9 * (cellId.x + 1))
401 self.assertEqual(ap_corr_map["base_GaussianFlux_instFluxErr"], 0.01 * (cellId.y + 1))
402 self.assertEqual(ap_corr_map["base_PsfFlux_instFlux"], 0.8 * (cellId.y + 2))
403 self.assertEqual(ap_corr_map["base_PsfFlux_instFluxErr"], 0.02 * (cellId.x + 2))
406class ExplodedCoaddTestCase(BaseMultipleCellCoaddTestCase):
407 """Test the construction and methods of an ExplodedCoadd instance."""
409 exploded_coadd: ExplodedCoadd
411 @classmethod
412 def setUpClass(cls) -> None: # noqa: D102
413 # Docstring inherited
414 super().setUpClass()
415 cls.exploded_coadd = cls.multiple_cell_coadd.explode()
417 @classmethod
418 def tearDownClass(cls) -> None: # noqa: D102
419 # Docstring inherited
420 del cls.exploded_coadd
421 super().tearDownClass()
423 def test_exploded_psf_image(self):
424 """Show that psf_image sizes are absurd."""
425 self.assertEqual(
426 self.exploded_coadd.psf_image.getBBox().getDimensions(),
427 geom.Extent2I(self.nx * self.psf_size_x, self.ny * self.psf_size_y),
428 )
429 for pad_psfs_with in (-999, -4, 0, 4, 8, 21, 40, 100):
430 exploded_coadd = self.multiple_cell_coadd.explode(pad_psfs_with=pad_psfs_with)
431 self.assertEqual(
432 exploded_coadd.psf_image.getBBox().getDimensions(),
433 geom.Extent2I(self.nx * self.outer_size_x, self.ny * self.outer_size_y),
434 )
436 def test_asMaskedImage(self):
437 """Test the asMaskedImage method for an ExplodedCoadd object."""
438 masked_image = self.exploded_coadd.asMaskedImage()
439 masked_image.setXY0(self.multiple_cell_coadd.outer_bbox.getMin())
440 base_bbox = self.multiple_cell_coadd.cells.first.outer.bbox
441 for cell_x, cell_y in product(range(self.nx), range(self.ny)):
442 bbox = base_bbox.shiftedBy(geom.Extent2I(cell_x * self.outer_size_x, cell_y * self.outer_size_y))
443 with self.subTest(cell_x=cell_x, cell_y=cell_y):
444 self.assertMaskedImagesEqual(
445 masked_image[bbox],
446 self.multiple_cell_coadd.cells[Index2D(cell_x, cell_y)].outer.asMaskedImage(),
447 )
450class StitchedCoaddTestCase(BaseMultipleCellCoaddTestCase):
451 """Test the construction and methods of a StitchedCoadd instance."""
453 stitched_coadd: StitchedCoadd
455 @classmethod
456 def setUpClass(cls) -> None: # noqa: D102
457 # Docstring inherited
458 super().setUpClass()
459 cls.stitched_coadd = cls.multiple_cell_coadd.stitch()
461 @classmethod
462 def tearDownClass(cls) -> None: # noqa: D102
463 # Docstring inherited
464 del cls.stitched_coadd
465 super().tearDownClass()
467 def test_computeBBox(self):
468 """Test the computeBBox method for a StitchedPsf object."""
469 stitched_psf = self.stitched_coadd.psf
471 psf_bbox = geom.Box2I(
472 geom.Point2I(-(self.psf_size_x // 2), -(self.psf_size_y // 2)),
473 geom.Extent2I(self.psf_size_x, self.psf_size_y),
474 )
476 for position, _ in self.test_positions:
477 bbox = stitched_psf.computeBBox(position)
478 self.assertEqual(bbox, psf_bbox)
480 def test_computeShape(self):
481 """Test the computeShape method for a StitchedPsf object."""
482 stitched_psf = self.stitched_coadd.psf
483 for position, cell_index in self.test_positions:
484 psf_shape = stitched_psf.computeShape(position) # check we can compute shape
485 self.assertIsNot(psf_shape.getIxx(), np.nan)
486 self.assertIsNot(psf_shape.getIyy(), np.nan)
487 self.assertIsNot(psf_shape.getIxy(), np.nan)
489 # Moments measured from pixellated images are significantly
490 # underestimated for small PSFs.
491 if self.psf_sigmas[cell_index] >= 1.0:
492 self.assertAlmostEqual(psf_shape.getIxx(), self.psf_sigmas[cell_index] ** 2, delta=1e-3)
493 self.assertAlmostEqual(psf_shape.getIyy(), self.psf_sigmas[cell_index] ** 2, delta=1e-3)
494 self.assertAlmostEqual(psf_shape.getIxy(), 0.0)
496 def test_computeKernelImage(self):
497 """Test the computeKernelImage method for a StitchedPsf object."""
498 stitched_psf = self.stitched_coadd.psf
499 psf_bbox = geom.Box2I(
500 geom.Point2I(-(self.psf_size_x // 2), -(self.psf_size_y // 2)),
501 geom.Extent2I(self.psf_size_x, self.psf_size_y),
502 )
504 for position, cell_index in self.test_positions:
505 image1 = stitched_psf.computeKernelImage(position)
506 image2 = SingleGaussianPsf(
507 self.psf_size_x, self.psf_size_y, self.psf_sigmas[cell_index]
508 ).computeKernelImage(position)
509 self.assertImagesEqual(image1, image2)
510 self.assertEqual(image1.getBBox(), psf_bbox)
512 def test_computeImage(self):
513 """Test the computeImage method for a StitchedPsf object."""
514 stitched_psf = self.stitched_coadd.psf
515 psf_extent = geom.Extent2I(self.psf_size_x, self.psf_size_y)
517 for position, cell_index in self.test_positions:
518 image1 = stitched_psf.computeImage(position)
519 image2 = SingleGaussianPsf(
520 self.psf_size_x, self.psf_size_y, self.psf_sigmas[cell_index]
521 ).computeImage(position)
522 self.assertImagesEqual(image1, image2)
523 self.assertEqual(image1.getBBox().getDimensions(), psf_extent)
525 def test_computeImage_computeKernelImage(self):
526 """Test that computeImage called at integer points gives the same
527 result as calling computeKernelImage.
528 """
529 stitched_psf = self.stitched_coadd.psf
530 for position, _cell_index in self.test_positions:
531 pos = geom.Point2D(geom.Point2I(position)) # round to integer
532 image1 = stitched_psf.computeKernelImage(pos)
533 image2 = stitched_psf.computeImage(pos)
534 self.assertImagesEqual(image1, image2)
536 def test_computeApetureFlux(self):
537 """Test the computeApertureFlux method for a StitchedPsf object."""
538 stitched_psf = self.stitched_coadd.psf
539 for position, cell_index in self.test_positions:
540 flux1sigma = stitched_psf.computeApertureFlux(self.psf_sigmas[cell_index], position=position)
541 self.assertAlmostEqual(flux1sigma, 0.39, delta=5e-2)
543 flux3sigma = stitched_psf.computeApertureFlux(
544 3.0 * self.psf_sigmas[cell_index], position=position
545 )
546 self.assertAlmostEqual(flux3sigma, 0.97, delta=2e-2)
548 def test_asExposure(self):
549 """Test the asExposure method for a StitchedCoadd object."""
550 exposure = self.stitched_coadd.asExposure()
552 # Check that the bounding box is correct.
553 bbox = exposure.getBBox()
554 self.assertEqual(bbox.getWidth(), self.inner_size_x * self.nx + 2 * self.border_size)
555 self.assertEqual(bbox.getHeight(), self.inner_size_y * self.ny + 2 * self.border_size)
557 for y in range(self.ny):
558 for x in range(self.nx):
559 bbox = geom.Box2I(
560 geom.Point2I(self.x0 + x * self.inner_size_x, self.y0 + y * self.inner_size_y),
561 geom.Extent2I(self.inner_size_x, self.inner_size_y),
562 )
563 index = Index2D(x=x, y=y)
564 self.assertMaskedImagesEqual(exposure[bbox], self.exposures[index][bbox])
566 self.assertMaskedImagesEqual(self.stitched_coadd.asExposure(noise_index=None), exposure)
567 for noise_index in range(self.n_noise_realizations):
568 noise_exposure = self.stitched_coadd.asExposure(noise_index=noise_index)
570 self.assertImagesEqual(noise_exposure.variance, exposure.variance)
571 self.assertImagesEqual(noise_exposure.mask, exposure.mask)
573 for y in range(self.ny):
574 for x in range(self.nx):
575 bbox = geom.Box2I(
576 geom.Point2I(self.x0 + x * self.inner_size_x, self.y0 + y * self.inner_size_y),
577 geom.Extent2I(self.inner_size_x, self.inner_size_y),
578 )
579 index = Index2D(x=x, y=y)
580 self.assertImagesEqual(
581 noise_exposure.image[bbox],
582 self.multiple_cell_coadd.cells[index].outer.noise_realizations[noise_index][bbox],
583 )
585 with self.assertRaises(ValueError):
586 self.stitched_coadd.asExposure(noise_index=self.n_noise_realizations)
588 def test_aperture_correction(self):
589 """Test the aperture correction values are what we expect."""
590 ap_corr_map = self.stitched_coadd.ap_corr_map
591 algorithm_names = ap_corr_map.keys()
593 for (position, cellId), algorithm_name in product(self.test_positions, algorithm_names):
594 with self.subTest(x=cellId.x, y=cellId.y, algorithm_name=algorithm_name):
595 field_name = algorithm_name
596 self.assertEqual(
597 ap_corr_map[field_name].evaluate(position),
598 self.multiple_cell_coadd.cells[cellId].aperture_correction_map[field_name],
599 )
601 def test_inputs(self):
602 """Test that the inputs are populated correctly on stitching."""
603 inputs = self.stitched_coadd.ccds
604 for position, cellId in self.test_positions:
605 with self.subTest(x=cellId.x, y=cellId.y):
606 self.assertEqual(
607 inputs.subsetContaining(position),
608 tuple(self.multiple_cell_coadd.cells[cellId].inputs.keys()),
609 )
610 self.assertEqual(len(self.stitched_coadd.visits), 1)
612 def test_borders(self):
613 """Test that the borders are populated correctly on stitching."""
614 mi = self.stitched_coadd.asMaskedImage()
616 # Check that the borders are not all zero.
617 with self.assertRaises(AssertionError):
618 np.testing.assert_array_equal(mi.image.array[: self.border_size, :], 0)
619 with self.assertRaises(AssertionError):
620 np.testing.assert_array_equal(mi.image.array[-self.border_size :, :], 0)
621 with self.assertRaises(AssertionError):
622 np.testing.assert_array_equal(mi.image.array[:, : self.border_size], 0)
623 with self.assertRaises(AssertionError):
624 np.testing.assert_array_equal(mi.image.array[:, -self.border_size :], 0)
626 def test_fits(self):
627 """Test that we can write an Exposure with StitchedPsf to a FITS file
628 and read it.
629 """
630 write_exposure = self.stitched_coadd.asExposure()
631 with lsst.utils.tests.getTempFilePath(".fits") as filename:
632 write_exposure.writeFits(filename)
633 read_exposure = ExposureF.readFits(filename) # Test the readFits method.
635 # Test that the image planes are identical.
636 self.assertMaskedImagesEqual(read_exposure, write_exposure)
638 # Test the PSF images in the StitchedPsf.
639 for index in write_exposure.psf.images.indices():
640 self.assertImagesEqual(read_exposure.psf.images[index], write_exposure.psf.images[index])
642 # Test that the WCSs are equal.
643 self.assertEqual(
644 read_exposure.wcs.getFitsMetadata().toString(),
645 write_exposure.wcs.getFitsMetadata().toString(),
646 )
649class TestMemory(lsst.utils.tests.MemoryTestCase):
650 """Check for resource/memory leaks."""
653def setup_module(module): # noqa: D103
654 lsst.utils.tests.init()
657if __name__ == "__main__": 657 ↛ 658line 657 didn't jump to line 658 because the condition on line 657 was never true
658 lsst.utils.tests.init()
659 unittest.main()