Coverage for tests/test_cell_coadd.py: 28%
102 statements
« prev ^ index » next coverage.py v7.14.1, created at 2026-06-12 07:50 +0000
1# This file is part of lsst-images.
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# Use of this source code is governed by a 3-clause BSD-style
10# license that can be found in the LICENSE file.
12from __future__ import annotations
14import os
15import pickle
16import unittest
17from typing import Any
19import numpy as np
21from lsst.images import YX, Box, Interval, get_legacy_deep_coadd_mask_planes
22from lsst.images.cells import CellCoadd, CellIJ
23from lsst.images.fits import FitsCompressionOptions
24from lsst.images.tests import (
25 DP2_COADD_DATA_ID,
26 DP2_COADD_MISSING_CELL,
27 RoundtripFits,
28 RoundtripJson,
29 RoundtripNdf,
30 assert_cell_coadds_equal,
31 assert_masked_images_equal,
32 assert_psfs_equal,
33 compare_cell_coadd_to_legacy,
34 compare_masked_image_to_legacy,
35 compare_psf_to_legacy,
36 compare_sky_projection_to_legacy_wcs,
37)
39try:
40 import h5py # noqa: F401
42 HAVE_H5PY = True
43except ImportError:
44 HAVE_H5PY = False
46DATA_DIR = os.environ.get("TESTDATA_IMAGES_DIR", None)
49@unittest.skipUnless(DATA_DIR is not None, "TESTDATA_IMAGES_DIR is not in the environment.")
50class CellCoaddTestCase(unittest.TestCase):
51 """Tests for the CellCoadd class and its many component classes."""
53 @classmethod
54 def setUpClass(cls) -> None:
55 assert DATA_DIR is not None, "Guaranteed by decorator."
56 cls.filename = os.path.join(DATA_DIR, "dp2", "legacy", "deep_coadd_cell_predetection.fits")
57 cls.plane_map = get_legacy_deep_coadd_mask_planes()
58 cls.missing_cell = CellIJ(**DP2_COADD_MISSING_CELL)
59 try:
60 from lsst.cell_coadds import MultipleCellCoadd
62 cls.legacy_cell_coadd = MultipleCellCoadd.read_fits(cls.filename)
63 except ImportError:
64 raise unittest.SkipTest("lsst.cell_coadds could not be imported.") from None
65 with open(os.path.join(DATA_DIR, "dp2", "legacy", "skyMap.pickle"), "rb") as stream:
66 cls.skymap = pickle.load(stream)
67 cls.cell_coadd = CellCoadd.from_legacy(
68 cls.legacy_cell_coadd,
69 plane_map=cls.plane_map,
70 tract_info=cls.skymap[DP2_COADD_DATA_ID["tract"]],
71 )
73 def make_psf_points(self, bbox: Box) -> YX[np.ndarray]:
74 """Make arrays of points to test PSFs at, given a bbox that is assumed
75 to be snapped to the cell_coadd grid.
76 """
77 xc, yc = np.meshgrid(
78 np.arange(
79 bbox.x.start + self.cell_coadd.grid.cell_shape.x * 0.5,
80 bbox.x.stop,
81 self.cell_coadd.grid.cell_shape.x,
82 ),
83 np.arange(
84 bbox.y.start + self.cell_coadd.grid.cell_shape.y * 0.5,
85 bbox.y.stop,
86 self.cell_coadd.grid.cell_shape.y,
87 ),
88 )
89 return YX(
90 y=yc.ravel() + self.rng.uniform(-0.4, 0.4, size=yc.size),
91 x=xc.ravel() + self.rng.uniform(-0.4, 0.4, size=xc.size),
92 )
94 def setUp(self) -> None:
95 self.rng = np.random.default_rng(44)
96 self.psf_points = self.make_psf_points(self.cell_coadd.bbox)
98 def test_from_legacy(self) -> None:
99 """Test constructing a CellCoadd by converting a legacy
100 lsst.cell_coadds.MultipleCellCoadd.
101 """
102 self.assertEqual(self.cell_coadd.bounds.missing, {self.missing_cell})
103 self.assertEqual(self.cell_coadd.bbox, Box.factory[12900:13500, 9600:10050])
104 compare_cell_coadd_to_legacy(
105 self,
106 self.cell_coadd,
107 self.legacy_cell_coadd,
108 tract_bbox=Box.from_legacy(self.skymap[DP2_COADD_DATA_ID["tract"]].getBBox()),
109 plane_map=self.plane_map,
110 psf_points=self.psf_points,
111 )
113 def test_roundtrip(self) -> None:
114 """Test serializing a CellCoadd and reading it back in, including
115 subimage and component reads.
116 """
117 with RoundtripFits(self, self.cell_coadd, "CellCoadd") as roundtrip:
118 # Check a subimage read. The subbox only overlaps (but does not
119 # fully cover) the middle 2 (of 4) cells in y, while covering
120 # exactly the last column of cells in x. It does not cover the
121 # missing cell.
122 subbox = Box.factory[
123 self.cell_coadd.bbox.y.start + 252 : self.cell_coadd.bbox.y.stop - 175,
124 self.cell_coadd.bbox.x.stop - 150 : self.cell_coadd.bbox.x.stop,
125 ]
126 subimage = roundtrip.get(bbox=subbox)
127 assert_masked_images_equal(self, subimage, self.cell_coadd[subbox], expect_view=False)
128 alternates: dict[str, Any] = {}
129 with self.subTest():
130 subpsf = roundtrip.get("psf", bbox=subbox)
131 self.assertEqual(
132 subpsf.bounds.bbox,
133 Box(
134 y=Interval.factory[
135 self.cell_coadd.bbox.y.start + 150 : self.cell_coadd.bbox.y.stop - 150
136 ],
137 x=subbox.x,
138 ),
139 )
140 assert_psfs_equal(self, subpsf, self.cell_coadd.psf, points=self.make_psf_points(subbox))
141 self.assertEqual(roundtrip.get("bbox"), self.cell_coadd.bbox)
142 alternates = {
143 k: roundtrip.get(k)
144 for k in [
145 "sky_projection",
146 "image",
147 "mask",
148 "variance",
149 "psf",
150 "aperture_corrections",
151 "provenance",
152 ]
153 }
154 with self.subTest():
155 backgrounds = roundtrip.get("backgrounds")
156 self.assertEqual(backgrounds.keys(), set())
157 self.assertIsNone(backgrounds.subtracted)
158 with roundtrip.inspect() as fits:
159 for extname in ["IMAGE", "MASK", "VARIANCE", "MASK_FRACTIONS/REJECTED"] + [
160 f"NOISE_REALIZATIONS/{n}" for n in range(len(self.cell_coadd.noise_realizations))
161 ]:
162 self.assertEqual(fits[extname].header["ZTILE1"], self.cell_coadd.grid.cell_shape.x)
163 self.assertEqual(fits[extname].header["ZTILE2"], self.cell_coadd.grid.cell_shape.y)
164 # Fixture self-consistency: bbox and missing-cell set are what setUp
165 # claims they are.
166 self.assertEqual(self.cell_coadd.bounds.missing, {self.missing_cell})
167 self.assertEqual(self.cell_coadd.bbox, Box.factory[12900:13500, 9600:10050])
168 # Full round-trip fidelity, including background contents.
169 assert_cell_coadds_equal(self, roundtrip.result, self.cell_coadd, expect_view=False)
170 compare_cell_coadd_to_legacy(
171 self,
172 roundtrip.result,
173 self.legacy_cell_coadd,
174 tract_bbox=Box.from_legacy(self.skymap[DP2_COADD_DATA_ID["tract"]].getBBox()),
175 plane_map=self.plane_map,
176 alternates=alternates,
177 psf_points=self.psf_points,
178 )
180 def test_fits_compression(self) -> None:
181 """Test writing with quantized FITS compression."""
182 with RoundtripFits(
183 self,
184 self.cell_coadd,
185 storage_class="CellCoadd",
186 recipe="lossy16",
187 compression_options={
188 "image": FitsCompressionOptions.LOSSY,
189 "variance": FitsCompressionOptions.LOSSY,
190 },
191 ) as roundtrip:
192 with roundtrip.inspect() as fits:
193 for extname in ["IMAGE", "MASK", "VARIANCE", "MASK_FRACTIONS/REJECTED"] + [
194 f"NOISE_REALIZATIONS/{n}" for n in range(len(self.cell_coadd.noise_realizations))
195 ]:
196 with self.subTest(extname=extname):
197 self.assertEqual(fits[extname].header["ZTILE1"], self.cell_coadd.grid.cell_shape.x)
198 self.assertEqual(fits[extname].header["ZTILE2"], self.cell_coadd.grid.cell_shape.y)
199 if extname == "MASK" or extname.startswith("MASK_FRACTIONS"):
200 self.assertEqual(fits[extname].header["ZCMPTYPE"], "GZIP_2")
201 else:
202 self.assertEqual(fits[extname].header["ZCMPTYPE"], "RICE_1")
203 self.assertEqual(fits[extname].header["ZQUANTIZ"], "SUBTRACTIVE_DITHER_2")
205 def test_fits_json_consistency(self) -> None:
206 """FITS and JSON backends produce equal CellCoadds on round-trip."""
207 with (
208 RoundtripFits(self, self.cell_coadd) as fits_rt,
209 RoundtripJson(self, self.cell_coadd) as json_rt,
210 ):
211 assert_cell_coadds_equal(self, self.cell_coadd, fits_rt.result, expect_view=False)
212 assert_cell_coadds_equal(self, self.cell_coadd, json_rt.result, expect_view=False)
213 assert_cell_coadds_equal(self, fits_rt.result, json_rt.result, expect_view=False)
215 def test_to_legacy(self) -> None:
216 """Test converting a CellCoadd back into a legacy MultipleCellCoadd."""
217 legacy_cell_coadd = self.cell_coadd.to_legacy()
218 compare_cell_coadd_to_legacy(
219 self,
220 self.cell_coadd,
221 legacy_cell_coadd,
222 tract_bbox=Box.from_legacy(self.skymap[DP2_COADD_DATA_ID["tract"]].getBBox()),
223 plane_map=self.plane_map,
224 psf_points=self.psf_points,
225 )
227 def test_to_legacy_exposure(self) -> None:
228 """Test converting a CellCoadd back into a legacy Exposure."""
229 legacy_exposure = self.cell_coadd.to_legacy_exposure()
231 self.assertEqual(legacy_exposure.getFilter().bandLabel, self.cell_coadd.band)
232 self.assertEqual(Box.from_legacy(legacy_exposure.getBBox()), self.cell_coadd.bbox)
233 compare_masked_image_to_legacy(
234 self, self.cell_coadd, legacy_exposure.maskedImage, plane_map=self.plane_map, expect_view=True
235 )
236 compare_psf_to_legacy(
237 self,
238 self.cell_coadd.psf,
239 legacy_exposure.getPsf(),
240 points=self.psf_points,
241 expect_legacy_raise_on_out_of_bounds=True,
242 )
243 compare_sky_projection_to_legacy_wcs(
244 self,
245 self.cell_coadd.sky_projection,
246 legacy_exposure.getWcs(),
247 self.cell_coadd.sky_projection.pixel_frame,
248 subimage_bbox=self.cell_coadd.bbox,
249 is_fits=True,
250 )
252 @unittest.skipUnless(HAVE_H5PY, "h5py is not installed")
253 def test_round_trip_ndf(self) -> None:
254 """NDF round-trip for CellCoadd, exercising hoisted long-named arrays.
256 This test covers the HDS name-shrinker fix for noise_realizations.
257 """
258 with RoundtripNdf(self, self.cell_coadd, "CellCoadd") as roundtrip:
259 assert_cell_coadds_equal(self, roundtrip.result, self.cell_coadd, expect_view=False)
261 @unittest.skipUnless(HAVE_H5PY, "h5py is not installed")
262 def test_fits_ndf_consistency(self) -> None:
263 """FITS and NDF backends produce equal CellCoadds on round-trip."""
264 with (
265 RoundtripFits(self, self.cell_coadd) as fits_rt,
266 RoundtripNdf(self, self.cell_coadd) as ndf_rt,
267 ):
268 assert_cell_coadds_equal(self, self.cell_coadd, fits_rt.result, expect_view=False)
269 assert_cell_coadds_equal(self, self.cell_coadd, ndf_rt.result, expect_view=False)
270 assert_cell_coadds_equal(self, fits_rt.result, ndf_rt.result, expect_view=False)
273if __name__ == "__main__":
274 unittest.main()