Coverage for tests/test_cell_coadd.py: 69%
276 statements
« prev ^ index » next coverage.py v7.15.2, created at 2026-07-30 08:38 +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 dataclasses
15import os
16import pickle
17from typing import Any
19import numpy as np
20import pytest
22from lsst.images import YX, Box, Interval, MaskPlane, get_legacy_deep_coadd_mask_planes
23from lsst.images.cells import CellCoadd, CellGrid, CellGridBounds, CellIJ, CoaddProvenance, PatchDefinition
24from lsst.images.describe import DescribeOptions, Report
25from lsst.images.fields import ChebyshevField
26from lsst.images.fits import FitsCompressionOptions
27from lsst.images.serialization import read_archive
28from lsst.images.tests import (
29 DP2_COADD_DATA_ID,
30 DP2_COADD_MISSING_CELL,
31 RoundtripFits,
32 RoundtripJson,
33 RoundtripNdf,
34 assert_cell_coadds_equal,
35 assert_images_equal,
36 assert_masked_images_equal,
37 assert_psfs_equal,
38 check_bounds_contains_broadcasting,
39 compare_cell_coadd_to_legacy,
40 compare_masked_image_to_legacy,
41 compare_psf_to_legacy,
42 compare_sky_projection_to_legacy_wcs,
43)
45try:
46 import h5py # noqa: F401
48 HAVE_H5PY = True
49except ImportError:
50 HAVE_H5PY = False
52try:
53 import lsst.afw.image # noqa: F401
54 from lsst.cell_coadds import MultipleCellCoadd as LegacyMultipleCellCoadd
56 HAVE_LEGACY = True
57except ImportError:
58 HAVE_LEGACY = False
59 type LegacyMultipleCellCoadd = Any # type: ignore[no-redef]
61EXTERNAL_DATA_DIR = os.environ.get("TESTDATA_IMAGES_DIR", None)
62LOCAL_DATA_DIR = os.path.join(os.path.dirname(__file__), "data")
64skip_no_h5py = pytest.mark.skipif(not HAVE_H5PY, reason="h5py is not installed")
65skip_no_legacy = pytest.mark.skipif(not HAVE_LEGACY, reason="lsst.afw (etc) could not be imported.")
68@dataclasses.dataclass
69class _LegacyTestData:
70 """A struct holding test data loaded from EXTERNAL_DATA_DIR."""
72 filename: str
73 tract_bbox: Box
74 legacy_cell_coadd: LegacyMultipleCellCoadd
75 cell_coadd: CellCoadd
76 plane_map: dict[str, MaskPlane] = dataclasses.field(default_factory=get_legacy_deep_coadd_mask_planes)
78 def make_psf_points(self, bbox: Box | None = None) -> YX[np.ndarray]:
79 """Create random PSF sample points within the given bbox."""
80 if bbox is None:
81 bbox = self.cell_coadd.bbox
82 rng = np.random.default_rng(44)
83 xc, yc = np.meshgrid(
84 np.arange(
85 bbox.x.start + self.cell_coadd.grid.cell_shape.x * 0.5,
86 bbox.x.stop,
87 self.cell_coadd.grid.cell_shape.x,
88 ),
89 np.arange(
90 bbox.y.start + self.cell_coadd.grid.cell_shape.y * 0.5,
91 bbox.y.stop,
92 self.cell_coadd.grid.cell_shape.y,
93 ),
94 )
95 return YX(
96 y=yc.ravel() + rng.uniform(-0.4, 0.4, size=yc.size),
97 x=xc.ravel() + rng.uniform(-0.4, 0.4, size=xc.size),
98 )
101@pytest.fixture(scope="session")
102def legacy_test_data() -> _LegacyTestData:
103 """Return a struct of CellCoadd loaded from legacy test data.
105 Skips if ``TESTDATA_IMAGES_DIR`` is not set or if ``lsst.cell_coadds``
106 cannot be imported.
107 """
108 if EXTERNAL_DATA_DIR is None: 108 ↛ 110line 108 didn't jump to line 110 because the condition on line 108 was always true
109 pytest.skip("TESTDATA_IMAGES_DIR is not in the environment.")
110 try:
111 from lsst.cell_coadds import MultipleCellCoadd
112 except ImportError:
113 pytest.skip("lsst.cell_coadds could not be imported.")
114 filename = os.path.join(EXTERNAL_DATA_DIR, "dp2", "legacy", "deep_coadd_cell_predetection.fits")
115 plane_map = get_legacy_deep_coadd_mask_planes()
116 legacy_cell_coadd = MultipleCellCoadd.read_fits(filename)
117 with open(os.path.join(EXTERNAL_DATA_DIR, "dp2", "legacy", "skyMap.pickle"), "rb") as stream:
118 skymap = pickle.load(stream)
119 cell_coadd = CellCoadd.from_legacy_cell_coadd(
120 legacy_cell_coadd,
121 plane_map=plane_map,
122 tract_info=skymap[DP2_COADD_DATA_ID["tract"]],
123 )
124 return _LegacyTestData(
125 filename=filename,
126 tract_bbox=Box.from_legacy(skymap[DP2_COADD_DATA_ID["tract"]].getBBox()),
127 legacy_cell_coadd=legacy_cell_coadd,
128 cell_coadd=cell_coadd,
129 )
132@pytest.fixture
133def minified_cell_coadd() -> CellCoadd:
134 """Return a tiny CellCoadd from JSON data stored in this package."""
135 path = os.path.join(LOCAL_DATA_DIR, "schema_v1", "legacy", "cell_coadd.json")
136 return read_archive(path, CellCoadd)
139def make_subbox(full_bbox: Box) -> Box:
140 """Make a box that's useful for nontrivial subimage tests.
142 This box only overlaps (but does not fully cover) the middle 2 (of 4)
143 cells in y, while covering exactly the last column of cells in x. It does
144 not cover the missing cell.
145 """
146 return Box.factory[
147 full_bbox.y.start + 252 : full_bbox.y.stop - 175,
148 full_bbox.x.stop - 150 : full_bbox.x.stop,
149 ]
152def test_cell_coadd_repr_str_pinned(minified_cell_coadd: CellCoadd) -> None:
153 """Pin the exact str and repr output of a CellCoadd.
155 A coadd takes component types no string can express, so repr is the
156 descriptive form rather than a constructor call. Both come from the
157 report, as they do for the sibling `VisitImage`.
158 """
159 assert str(minified_cell_coadd) == "CellCoadd([y=48:60, x=36:48], tract=9813)"
160 assert repr(minified_cell_coadd) == "<CellCoadd([y=48:60, x=36:48], tract=9813)>"
163def report_fields(report: Report) -> dict[str, Any]:
164 """Return a mapping from report field label to field value."""
165 return {field.label: field.value for field in report.fields}
168def make_test_provenance() -> CoaddProvenance:
169 """Return a provenance with three input images taken over two nights,
170 whose contributions cover three cells unevenly.
171 """
172 inputs = CoaddProvenance.make_empty_input_table(3)
173 inputs["instrument"] = "LSSTCam"
174 inputs["physical_filter"] = "r_57"
175 inputs["visit"] = [101, 102, 103]
176 inputs["detector"] = [1, 2, 3]
177 inputs["day_obs"] = [20250520, 20250520, 20250521]
178 contributions = CoaddProvenance.make_empty_contribution_table(6)
179 contributions["cell_i"] = [0, 0, 0, 1, 1, 2]
180 contributions["cell_j"] = 0
181 return CoaddProvenance(inputs=inputs, contributions=contributions)
184def test_coadd_provenance_repr_str_pinned(minified_cell_coadd: CellCoadd) -> None:
185 """Pin the str and repr of a CoaddProvenance.
187 Neither table can be expressed in an eval-able string, so repr is the
188 descriptive form, and both report only what len() can answer.
189 """
190 provenance = minified_cell_coadd.provenance
191 assert str(provenance) == "CoaddProvenance(6 input images)"
192 assert repr(provenance) == "<CoaddProvenance(6 input images)>"
195def test_coadd_provenance_brief_report_skips_column_scans(minified_cell_coadd: CellCoadd) -> None:
196 """The brief report carries no fields, so repr and str never scan a
197 column of either table.
198 """
199 report = minified_cell_coadd.provenance._describe(DescribeOptions(brief=True))
200 assert not report.fields
201 assert report.to_repr() == "<CoaddProvenance(6 input images)>"
204def test_coadd_provenance_report_summarizes_its_tables(minified_cell_coadd: CellCoadd) -> None:
205 """The expanded report summarizes both tables instead of rendering them."""
206 report = minified_cell_coadd.provenance.describe()
207 assert report.type_name == "CoaddProvenance"
208 assert report_fields(report) == {
209 "instrument": "LSSTCam",
210 "physical_filter": "r_57",
211 "input images": "6 from 2 visits",
212 "day_obs": "20250520",
213 "cells": "3 with contributions",
214 "per cell": "2 input images",
215 }
216 assert list(report_fields(report)) == [
217 "instrument",
218 "physical_filter",
219 "input images",
220 "day_obs",
221 "cells",
222 "per cell",
223 ]
224 # Rich renders an astropy table as plain text and never consults its
225 # _repr_html_, so the report tabulates nothing.
226 assert not report.tables
227 assert "input images" in report._repr_html_()
230def test_coadd_provenance_report_counts_cells_against_bounds(minified_cell_coadd: CellCoadd) -> None:
231 """Given the image's cells, the report states coverage as a fraction."""
232 report = minified_cell_coadd.provenance._describe(bounds=minified_cell_coadd.bounds)
233 assert report_fields(report)["cells"] == "3 of 3 with contributions"
236def test_coadd_provenance_report_shows_partial_cell_coverage() -> None:
237 """Cells with image data but no contribution rows show up in the ratio."""
238 grid = CellGrid(bbox=Box.from_shape((30, 30)), cell_shape=YX(10, 10))
239 bounds = CellGridBounds(grid=grid, bbox=Box.factory[0:30, 0:30])
240 report = make_test_provenance()._describe(bounds=bounds)
241 assert report_fields(report)["cells"] == "3 of 9 with contributions"
244def test_coadd_provenance_report_ranges_over_nights_and_cells() -> None:
245 """Several visits, several nights and uneven cell coverage give the
246 ranged forms of each field.
247 """
248 fields = report_fields(make_test_provenance().describe())
249 assert fields["input images"] == "3 from 3 visits"
250 assert fields["day_obs"] == "20250520 - 20250521"
251 assert fields["cells"] == "3 with contributions"
252 assert fields["per cell"] == "1 - 3 input images (median 2)"
255def test_coadd_provenance_report_handles_empty_tables() -> None:
256 """A provenance with no rows describes without raising."""
257 provenance = CoaddProvenance(
258 inputs=CoaddProvenance.make_empty_input_table(0),
259 contributions=CoaddProvenance.make_empty_contribution_table(0),
260 )
261 assert str(provenance) == "CoaddProvenance(no input images)"
262 report = provenance.describe()
263 assert report_fields(report) == {"input images": "none", "cells": "none"}
264 report.__rich__()
267def test_coadd_provenance_report_counts_one_input_image_in_the_singular() -> None:
268 """Pin the singular forms of the counted phrases.
270 One input image is what the standalone ``coadd_provenance`` file holds,
271 and what a single-cell subset of a coadd holds, so this is the case the
272 ``describe`` command line subcommand shows most often.
273 """
274 inputs = CoaddProvenance.make_empty_input_table(1)
275 inputs["instrument"] = "LSSTCam"
276 inputs["physical_filter"] = "r_57"
277 inputs["visit"] = [101]
278 inputs["detector"] = [1]
279 inputs["day_obs"] = [20250520]
280 provenance = CoaddProvenance(
281 inputs=inputs, contributions=CoaddProvenance.make_empty_contribution_table(1)
282 )
283 assert str(provenance) == "CoaddProvenance(1 input image)"
284 assert repr(provenance) == "<CoaddProvenance(1 input image)>"
285 fields = report_fields(provenance.describe())
286 assert fields["input images"] == "1 from 1 visit"
287 assert fields["per cell"] == "1 input image"
290def test_cell_coadd_report_accounts_for_every_component(minified_cell_coadd: CellCoadd) -> None:
291 """Nothing the coadd carries is absent from its report."""
292 report = minified_cell_coadd.describe()
293 fields = report_fields(report)
294 assert fields["mask_fractions"] == "rejected"
295 assert fields["noise_realizations"] == "1 image (dtype float32)"
296 assert fields["aperture_corrections"] == "3 fields"
297 # Provenance is a child, so it needs no field line of its own.
298 assert "provenance" not in fields
299 assert list(report.children) == [
300 "image",
301 "mask",
302 "variance",
303 "sky_projection",
304 "psf",
305 "provenance",
306 "backgrounds",
307 ]
308 provenance = report.children["provenance"]
309 assert provenance.type_name == "CoaddProvenance"
310 # The coadd passed its cells down, so coverage is stated as a fraction.
311 assert report_fields(provenance)["cells"] == "3 of 3 with contributions"
314def test_cell_coadd_report_counts_aperture_corrections_only(minified_cell_coadd: CellCoadd) -> None:
315 """A coadd can carry dozens of aperture corrections with long names, so
316 the report counts them and never names them, detail or not.
317 """
318 plain = report_fields(minified_cell_coadd.describe())["aperture_corrections"]
319 assert plain == "3 fields"
320 detailed = report_fields(minified_cell_coadd.describe(detail=True))["aperture_corrections"]
321 assert detailed == "3 fields" # No additional detail
324def test_cell_coadd_report_states_absent_provenance(minified_cell_coadd: CellCoadd) -> None:
325 """A coadd with no provenance says so, while components that are merely
326 empty stay out of the report.
327 """
328 bare = CellCoadd(
329 minified_cell_coadd.image,
330 mask=minified_cell_coadd.mask,
331 variance=minified_cell_coadd.variance,
332 sky_projection=minified_cell_coadd.sky_projection,
333 band=minified_cell_coadd.band,
334 psf=minified_cell_coadd.psf,
335 patch=minified_cell_coadd.patch,
336 )
337 report = bare.describe()
338 fields = report_fields(report)
339 assert fields["provenance"] == "none"
340 assert "provenance" not in report.children
341 assert "mask_fractions" not in fields
342 assert "noise_realizations" not in fields
343 assert "aperture_corrections" not in fields
346def test_cell_grid_patch_str_uses_clean_geometry() -> None:
347 """CellGrid, PatchDefinition and CellGridBounds str drop the
348 Interval/YX/Box/CellIJ wrappers.
350 The report renders field values with str, so these must use the compact
351 geometry forms rather than pydantic's default field-by-field repr.
352 """
353 grid = CellGrid(bbox=Box.from_shape((100, 200)), cell_shape=YX(10, 20))
354 assert str(grid) == "[y=0:100, x=0:200], cell_shape=(y=10, x=20)"
356 patch = PatchDefinition(id=73, index=YX(7, 3), inner_bbox=Box.factory[1:3, 2:4], cells=grid)
357 assert str(patch) == (
358 "id=73, index=(y=7, x=3), inner_bbox=[y=1:3, x=2:4], "
359 "cells=([y=0:100, x=0:200], cell_shape=(y=10, x=20))"
360 )
361 # repr stays as the pydantic default so it remains eval-ish and distinct.
362 assert "Interval(" in repr(patch)
363 assert "YX(" in repr(patch)
365 bounds = CellGridBounds(grid=grid, bbox=Box.factory[0:40, 0:60])
366 assert str(bounds) == "[y=0:40, x=0:60] in grid ([y=0:100, x=0:200], cell_shape=(y=10, x=20))"
367 bounds_missing = CellGridBounds(
368 grid=grid, bbox=Box.factory[0:40, 0:60], missing=frozenset({CellIJ(i=1, j=1), CellIJ(i=0, j=2)})
369 )
370 assert str(bounds_missing) == (
371 "[y=0:40, x=0:60] in grid ([y=0:100, x=0:200], cell_shape=(y=10, x=20)), "
372 "missing={(i=0, j=2), (i=1, j=1)}"
373 )
374 assert "Interval(" in repr(bounds_missing)
377def test_from_legacy(legacy_test_data: _LegacyTestData) -> None:
378 """Test constructing a CellCoadd by converting a legacy
379 ``MultipleCellCoadd``.
380 """
381 assert legacy_test_data.cell_coadd.bounds.missing == {CellIJ(**DP2_COADD_MISSING_CELL)}
382 assert legacy_test_data.cell_coadd.bbox == Box.factory[12900:13500, 9600:10050]
383 compare_cell_coadd_to_legacy(
384 legacy_test_data.cell_coadd,
385 legacy_test_data.legacy_cell_coadd,
386 tract_bbox=legacy_test_data.tract_bbox,
387 plane_map=legacy_test_data.plane_map,
388 psf_points=legacy_test_data.make_psf_points(),
389 )
392def test_roundtrip(legacy_test_data: _LegacyTestData) -> None:
393 """Test that a CellCoadd roundtrips through FITS."""
394 with RoundtripFits(legacy_test_data.cell_coadd, "CellCoadd") as roundtrip:
395 # Check a subimage read (no component arg — does not trigger a skip).
396 subbox = Box.factory[
397 legacy_test_data.cell_coadd.bbox.y.start + 252 : legacy_test_data.cell_coadd.bbox.y.stop - 175,
398 legacy_test_data.cell_coadd.bbox.x.stop - 150 : legacy_test_data.cell_coadd.bbox.x.stop,
399 ]
400 subimage = roundtrip.get(bbox=subbox)
401 assert_masked_images_equal(subimage, legacy_test_data.cell_coadd[subbox], expect_view=False)
402 with roundtrip.inspect() as fits:
403 for extname in ["IMAGE", "MASK", "VARIANCE", "MASK_FRACTIONS/REJECTED"] + [
404 f"NOISE_REALIZATIONS/{n}" for n in range(len(legacy_test_data.cell_coadd.noise_realizations))
405 ]:
406 assert fits[extname].header["ZTILE1"] == legacy_test_data.cell_coadd.grid.cell_shape.x
407 assert fits[extname].header["ZTILE2"] == legacy_test_data.cell_coadd.grid.cell_shape.y
408 # Fixture self-consistency: bbox and missing-cell set are as expected.
409 assert legacy_test_data.cell_coadd.bounds.missing == {CellIJ(**DP2_COADD_MISSING_CELL)}
410 assert legacy_test_data.cell_coadd.bbox == Box.factory[12900:13500, 9600:10050]
411 # Full round-trip fidelity.
412 assert_cell_coadds_equal(roundtrip.result, legacy_test_data.cell_coadd, expect_view=False)
413 compare_cell_coadd_to_legacy(
414 roundtrip.result,
415 legacy_test_data.legacy_cell_coadd,
416 tract_bbox=legacy_test_data.tract_bbox,
417 plane_map=legacy_test_data.plane_map,
418 psf_points=legacy_test_data.make_psf_points(),
419 )
422def test_roundtrip_components(legacy_test_data: _LegacyTestData) -> None:
423 """Test component and subimage reads.
425 This test will be skipped if `lsst.daf.butler` is not available instead of
426 falling back to non-butler I/O, which is why we don't want to merge it
427 with `test_roundtrip`.
428 """
429 with RoundtripFits(legacy_test_data.cell_coadd, "CellCoadd") as roundtrip:
430 subbox = make_subbox(legacy_test_data.cell_coadd.bbox)
431 subpsf = roundtrip.get("psf", bbox=subbox)
432 assert subpsf.bounds.bbox == Box(
433 y=Interval.factory[
434 legacy_test_data.cell_coadd.bbox.y.start + 150 : legacy_test_data.cell_coadd.bbox.y.stop - 150
435 ],
436 x=subbox.x,
437 )
438 assert_psfs_equal(
439 subpsf,
440 legacy_test_data.cell_coadd.psf,
441 points=legacy_test_data.make_psf_points(subbox),
442 )
443 assert roundtrip.get("bbox") == legacy_test_data.cell_coadd.bbox
444 alternates = {
445 k: roundtrip.get(k)
446 for k in [
447 "sky_projection",
448 "image",
449 "mask",
450 "variance",
451 "masked_image",
452 "psf",
453 "aperture_corrections",
454 "provenance",
455 "backgrounds",
456 "bbox",
457 ]
458 }
459 # Read all the components at once.
460 all_components = roundtrip.get("components")
461 assert set(all_components) == set(alternates) - {"masked_image"}
462 assert all_components["bbox"] == alternates["bbox"]
463 assert_psfs_equal(all_components["psf"], alternates["psf"])
464 assert_images_equal(all_components["image"], alternates["image"])
466 backgrounds = roundtrip.get("backgrounds")
467 assert backgrounds.keys() == set()
468 assert backgrounds.subtracted is None
470 compare_cell_coadd_to_legacy(
471 roundtrip.result,
472 legacy_test_data.legacy_cell_coadd,
473 tract_bbox=legacy_test_data.tract_bbox,
474 plane_map=legacy_test_data.plane_map,
475 alternates=alternates,
476 psf_points=legacy_test_data.make_psf_points(),
477 )
480def test_fits_compression(legacy_test_data: _LegacyTestData) -> None:
481 """Test lossy FITS compression produces the expected headers."""
482 with RoundtripFits(
483 legacy_test_data.cell_coadd,
484 storage_class="CellCoadd",
485 recipe="lossy16",
486 compression_options={
487 "image": FitsCompressionOptions.LOSSY,
488 "variance": FitsCompressionOptions.LOSSY,
489 },
490 ) as roundtrip:
491 with roundtrip.inspect() as fits:
492 for extname in ["IMAGE", "MASK", "VARIANCE", "MASK_FRACTIONS/REJECTED"] + [
493 f"NOISE_REALIZATIONS/{n}" for n in range(len(legacy_test_data.cell_coadd.noise_realizations))
494 ]:
495 assert fits[extname].header["ZTILE1"] == legacy_test_data.cell_coadd.grid.cell_shape.x
496 assert fits[extname].header["ZTILE2"] == legacy_test_data.cell_coadd.grid.cell_shape.y
497 if extname == "MASK" or extname.startswith("MASK_FRACTIONS"):
498 assert fits[extname].header["ZCMPTYPE"] == "GZIP_2"
499 else:
500 assert fits[extname].header["ZCMPTYPE"] == "RICE_1"
501 assert fits[extname].header["ZQUANTIZ"] == "SUBTRACTIVE_DITHER_2"
504def test_json_roundtrip(legacy_test_data: _LegacyTestData) -> None:
505 """Verify a CellCoadd round-trips correctly through the JSON archive."""
506 with RoundtripJson(legacy_test_data.cell_coadd) as roundtrip:
507 pass
508 assert_cell_coadds_equal(roundtrip.result, legacy_test_data.cell_coadd, expect_view=False)
511def test_to_legacy_cell_coadd(legacy_test_data: _LegacyTestData) -> None:
512 """Verify converting a CellCoadd back into a legacy MultipleCellCoadd."""
513 legacy_cell_coadd = legacy_test_data.cell_coadd.to_legacy_cell_coadd()
514 compare_cell_coadd_to_legacy(
515 legacy_test_data.cell_coadd,
516 legacy_cell_coadd,
517 tract_bbox=legacy_test_data.tract_bbox,
518 plane_map=legacy_test_data.plane_map,
519 psf_points=legacy_test_data.make_psf_points(),
520 )
521 with pytest.raises(
522 ValueError, match="MultipleCellCoadd requires its bounding box to lie on the cell grid."
523 ):
524 legacy_test_data.cell_coadd[make_subbox(legacy_test_data.cell_coadd.bbox)].to_legacy_cell_coadd()
527@skip_no_legacy
528def test_to_legacy(legacy_test_data: _LegacyTestData) -> None:
529 """Test converting a CellCoadd back into a legacy Exposure."""
530 legacy_exposure = legacy_test_data.cell_coadd.to_legacy()
531 assert legacy_exposure.getFilter().bandLabel == legacy_test_data.cell_coadd.band
532 assert Box.from_legacy(legacy_exposure.getBBox()) == legacy_test_data.cell_coadd.bbox
533 compare_masked_image_to_legacy(
534 legacy_test_data.cell_coadd,
535 legacy_exposure.maskedImage,
536 plane_map=legacy_test_data.plane_map,
537 expect_view=True,
538 )
539 compare_psf_to_legacy(
540 legacy_test_data.cell_coadd.psf,
541 legacy_exposure.getPsf(),
542 points=legacy_test_data.make_psf_points(),
543 expect_legacy_raise_on_out_of_bounds=True,
544 )
545 compare_sky_projection_to_legacy_wcs(
546 legacy_test_data.cell_coadd.sky_projection,
547 legacy_exposure.getWcs(),
548 legacy_test_data.cell_coadd.sky_projection.pixel_frame,
549 subimage_bbox=legacy_test_data.cell_coadd.bbox,
550 is_fits=True,
551 )
552 subbox = make_subbox(legacy_test_data.cell_coadd.bbox)
553 compare_masked_image_to_legacy(
554 legacy_test_data.cell_coadd[subbox],
555 legacy_test_data.cell_coadd[subbox].to_legacy().maskedImage,
556 plane_map=legacy_test_data.plane_map,
557 expect_view=True,
558 )
561@skip_no_h5py
562def test_ndf_roundtrip(legacy_test_data: _LegacyTestData) -> None:
563 """Test that CellCoadd round-trips through NDF."""
564 with RoundtripNdf(legacy_test_data.cell_coadd, "CellCoadd") as roundtrip:
565 assert_cell_coadds_equal(roundtrip.result, legacy_test_data.cell_coadd, expect_view=False)
568# The float32 image pixels are of order 10 nJy, as is the test background
569# below, so background arithmetic is only reproducible to a few ULPs at that
570# scale.
571BACKGROUND_ATOL = 1e-5
574def _add_gradient_background(coadd: CellCoadd, name: str = "pretty") -> ChebyshevField:
575 """Add a non-constant background over the coadd's full bbox and return the
576 field that was added.
577 """
578 field = ChebyshevField(
579 coadd.bbox,
580 np.array([[10.0, 2.0, 0.5], [1.0, 0.75, 0.0], [0.25, 0.0, 0.0]]),
581 unit=coadd.unit,
582 )
583 coadd.backgrounds.add(name, field, description="Gradient background for tests.")
584 return field
587def _make_background_subbox(coadd: CellCoadd) -> Box:
588 """Make a box that trims a different number of pixels from each side of the
589 coadd, so a background rendered over the wrong box cannot match by chance.
590 """
591 return Box.factory[
592 coadd.bbox.y.start + 2 : coadd.bbox.y.stop - 3,
593 coadd.bbox.x.start + 1 : coadd.bbox.x.stop - 2,
594 ]
597def test_apply_background_subimage(minified_cell_coadd: CellCoadd) -> None:
598 """Test that applying a background to a subimage subtracts only the
599 portion of the background model that overlaps the subimage.
600 """
601 field = _add_gradient_background(minified_cell_coadd)
602 subbox = _make_background_subbox(minified_cell_coadd)
603 subimage = minified_cell_coadd[subbox]
604 original = subimage.image.array.copy()
606 subimage.apply_background("pretty")
608 assert subimage.backgrounds.subtracted is not None
609 assert subimage.backgrounds.subtracted.name == "pretty"
610 assert subimage.image.bbox == subbox
611 expected = field.render(subbox, dtype=subimage.image.array.dtype).array
612 np.testing.assert_allclose(original - subimage.image.array, expected, atol=BACKGROUND_ATOL)
615def test_restore_background_subimage(minified_cell_coadd: CellCoadd) -> None:
616 """Test that restoring the original background on a subimage adds back
617 only the portion of the background model that overlaps the subimage.
618 """
619 _add_gradient_background(minified_cell_coadd)
620 subbox = _make_background_subbox(minified_cell_coadd)
621 subimage = minified_cell_coadd[subbox]
622 original = subimage.image.array.copy()
623 subimage.apply_background("pretty")
625 subimage.apply_background(None)
627 assert subimage.backgrounds.subtracted is None
628 np.testing.assert_allclose(subimage.image.array, original, atol=BACKGROUND_ATOL)
631def test_apply_background_after_bounded_read(minified_cell_coadd: CellCoadd) -> None:
632 """Test that a background can be applied to a coadd read with a bbox
633 parameter.
634 """
635 field = _add_gradient_background(minified_cell_coadd)
636 subbox = _make_background_subbox(minified_cell_coadd)
637 with RoundtripJson(minified_cell_coadd, "CellCoadd") as roundtrip:
638 subimage = roundtrip.get(bbox=subbox)
639 original = subimage.image.array.copy()
641 subimage.apply_background("pretty")
643 assert subimage.image.bbox == subbox
644 expected = field.render(subbox, dtype=subimage.image.array.dtype).array
645 np.testing.assert_allclose(original - subimage.image.array, expected, atol=BACKGROUND_ATOL)
648def test_cell_grid_bounds_contains_broadcasting(minified_cell_coadd: CellCoadd) -> None:
649 """Test that CellGridBounds.contains broadcasts like a numpy ufunc."""
650 assert minified_cell_coadd.bounds.missing, "fixture should retain a missing cell"
651 check_bounds_contains_broadcasting(minified_cell_coadd.bounds)
654def test_intersection_bounds_contains_broadcasting(minified_cell_coadd: CellCoadd) -> None:
655 """Test that IntersectionBounds.contains broadcasts like a numpy ufunc."""
656 # Clip the CellGridBounds with a Box offset by 1 pixel on each side so it
657 # does not snap to any cell boundary, forcing a lazy IntersectionBounds.
658 bounds = minified_cell_coadd.bounds
659 clip = Box.factory[
660 bounds.bbox.y.start + 1 : bounds.bbox.y.stop - 1,
661 bounds.bbox.x.start + 1 : bounds.bbox.x.stop - 1,
662 ]
663 check_bounds_contains_broadcasting(bounds.intersection(clip))