Coverage for tests/test_cell_coadd.py: 75%
341 statements
« prev ^ index » next coverage.py v7.16.1, created at 2026-09-23 10:30 +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 copy
15import dataclasses
16import json
17import os
18import pickle
19from pathlib import Path
20from typing import Any
22import numpy as np
23import pytest
25from lsst.images import YX, BoundsError, Box, Interval, MaskPlane, get_legacy_deep_coadd_mask_planes
26from lsst.images.cells import (
27 CellCoadd,
28 CellGrid,
29 CellGridBounds,
30 CellIJ,
31 CellPointSpreadFunctionSerializationModel,
32 CoaddProvenance,
33 PatchDefinition,
34)
35from lsst.images.describe import DescribeOptions, Report
36from lsst.images.fields import ChebyshevField
37from lsst.images.fits import FitsCompressionOptions
38from lsst.images.serialization import (
39 JsonRef,
40 class_for_schema,
41 parameterize_tree,
42 read_archive,
43)
44from lsst.images.tests import (
45 DP2_COADD_DATA_ID,
46 DP2_COADD_MISSING_CELL,
47 RoundtripFits,
48 RoundtripJson,
49 RoundtripNdf,
50 assert_cell_coadds_equal,
51 assert_images_equal,
52 assert_masked_images_equal,
53 assert_psfs_equal,
54 check_bounds_contains_broadcasting,
55 compare_cell_coadd_to_legacy,
56 compare_masked_image_to_legacy,
57 compare_psf_to_legacy,
58 compare_sky_projection_to_legacy_wcs,
59 current_fixture_path,
60 reset_afw_mask_planes, # noqa: F401
61)
63try:
64 import h5py # noqa: F401
66 HAVE_H5PY = True
67except ImportError:
68 HAVE_H5PY = False
70try:
71 import lsst.afw.image # noqa: F401
72 from lsst.cell_coadds import MultipleCellCoadd as LegacyMultipleCellCoadd
74 HAVE_LEGACY = True
75except ImportError:
76 HAVE_LEGACY = False
77 type LegacyMultipleCellCoadd = Any # type: ignore[no-redef]
79EXTERNAL_DATA_DIR = os.environ.get("TESTDATA_IMAGES_DIR", None)
80FIXTURE_DIR = Path(__file__).parent / "data" / "schemas"
82skip_no_h5py = pytest.mark.skipif(not HAVE_H5PY, reason="h5py is not installed")
83skip_no_legacy = pytest.mark.skipif(not HAVE_LEGACY, reason="lsst.afw (etc) could not be imported.")
86@dataclasses.dataclass
87class _LegacyTestData:
88 """A struct holding test data loaded from EXTERNAL_DATA_DIR."""
90 filename: str
91 tract_bbox: Box
92 legacy_cell_coadd: LegacyMultipleCellCoadd
93 cell_coadd: CellCoadd
94 plane_map: dict[str, MaskPlane] = dataclasses.field(default_factory=get_legacy_deep_coadd_mask_planes)
96 def make_psf_points(self, bbox: Box | None = None) -> YX[np.ndarray]:
97 """Create random PSF sample points within the given bbox."""
98 if bbox is None:
99 bbox = self.cell_coadd.bbox
100 rng = np.random.default_rng(44)
101 xc, yc = np.meshgrid(
102 np.arange(
103 bbox.x.start + self.cell_coadd.grid.cell_shape.x * 0.5,
104 bbox.x.stop,
105 self.cell_coadd.grid.cell_shape.x,
106 ),
107 np.arange(
108 bbox.y.start + self.cell_coadd.grid.cell_shape.y * 0.5,
109 bbox.y.stop,
110 self.cell_coadd.grid.cell_shape.y,
111 ),
112 )
113 return YX(
114 y=yc.ravel() + rng.uniform(-0.4, 0.4, size=yc.size),
115 x=xc.ravel() + rng.uniform(-0.4, 0.4, size=xc.size),
116 )
119@pytest.fixture
120def legacy_test_data(reset_afw_mask_planes: None) -> _LegacyTestData: # noqa: F811
121 """Return a struct of CellCoadd loaded from legacy test data.
123 Skips if ``TESTDATA_IMAGES_DIR`` is not set or if ``lsst.cell_coadds``
124 cannot be imported.
125 """
126 if EXTERNAL_DATA_DIR is None: 126 ↛ 128line 126 didn't jump to line 128 because the condition on line 126 was always true
127 pytest.skip("TESTDATA_IMAGES_DIR is not in the environment.")
128 try:
129 from lsst.cell_coadds import MultipleCellCoadd
130 except ImportError:
131 pytest.skip("lsst.cell_coadds could not be imported.")
132 filename = os.path.join(EXTERNAL_DATA_DIR, "dp2", "legacy", "deep_coadd_cell_predetection.fits")
133 plane_map = get_legacy_deep_coadd_mask_planes()
134 legacy_cell_coadd = MultipleCellCoadd.read_fits(filename)
135 with open(os.path.join(EXTERNAL_DATA_DIR, "dp2", "legacy", "skyMap.pickle"), "rb") as stream:
136 skymap = pickle.load(stream)
137 cell_coadd = CellCoadd.from_legacy_cell_coadd(
138 legacy_cell_coadd,
139 plane_map=plane_map,
140 tract_info=skymap[DP2_COADD_DATA_ID["tract"]],
141 )
142 return _LegacyTestData(
143 filename=filename,
144 tract_bbox=Box.from_legacy(skymap[DP2_COADD_DATA_ID["tract"]].getBBox()),
145 legacy_cell_coadd=legacy_cell_coadd,
146 cell_coadd=cell_coadd,
147 )
150@pytest.fixture
151def minified_cell_coadd() -> CellCoadd:
152 """Return a tiny CellCoadd from JSON data stored in this package."""
153 path = current_fixture_path(FIXTURE_DIR, "cell_coadd", variant="as_shipped")
154 return read_archive(path, CellCoadd)
157def make_subbox(full_bbox: Box) -> Box:
158 """Make a box that's useful for nontrivial subimage tests.
160 This box only overlaps (but does not fully cover) the middle 2 (of 4)
161 cells in y, while covering exactly the last column of cells in x. It does
162 not cover the missing cell.
163 """
164 return Box.factory[
165 full_bbox.y.start + 252 : full_bbox.y.stop - 175,
166 full_bbox.x.stop - 150 : full_bbox.x.stop,
167 ]
170def test_cell_coadd_repr_str_pinned(minified_cell_coadd: CellCoadd) -> None:
171 """Pin the exact str and repr output of a CellCoadd.
173 A coadd takes component types no string can express, so repr is the
174 descriptive form rather than a constructor call. Both come from the
175 report, as they do for the sibling `VisitImage`.
176 """
177 assert str(minified_cell_coadd) == "CellCoadd([y=48:60, x=36:48], tract=9813)"
178 assert repr(minified_cell_coadd) == "<CellCoadd([y=48:60, x=36:48], tract=9813)>"
181def report_fields(report: Report) -> dict[str, Any]:
182 """Return a mapping from report field label to field value."""
183 return {field.label: field.value for field in report.fields}
186def make_test_provenance() -> CoaddProvenance:
187 """Return a provenance with three input images taken over two nights,
188 whose contributions cover three cells unevenly.
189 """
190 inputs = CoaddProvenance.make_empty_input_table(3)
191 inputs["instrument"] = "LSSTCam"
192 inputs["physical_filter"] = "r_57"
193 inputs["visit"] = [101, 102, 103]
194 inputs["detector"] = [1, 2, 3]
195 inputs["day_obs"] = [20250520, 20250520, 20250521]
196 contributions = CoaddProvenance.make_empty_contribution_table(6)
197 contributions["cell_i"] = [0, 0, 0, 1, 1, 2]
198 contributions["cell_j"] = 0
199 return CoaddProvenance(inputs=inputs, contributions=contributions)
202def test_coadd_provenance_repr_str_pinned(minified_cell_coadd: CellCoadd) -> None:
203 """Pin the str and repr of a CoaddProvenance.
205 Neither table can be expressed in an eval-able string, so repr is the
206 descriptive form, and both report only what len() can answer.
207 """
208 provenance = minified_cell_coadd.provenance
209 assert str(provenance) == "CoaddProvenance(6 input images)"
210 assert repr(provenance) == "<CoaddProvenance(6 input images)>"
213def test_coadd_provenance_brief_report_skips_column_scans(minified_cell_coadd: CellCoadd) -> None:
214 """The brief report carries no fields, so repr and str never scan a
215 column of either table.
216 """
217 report = minified_cell_coadd.provenance._describe(DescribeOptions(brief=True))
218 assert not report.fields
219 assert report.to_repr() == "<CoaddProvenance(6 input images)>"
222def test_coadd_provenance_report_summarizes_its_tables(minified_cell_coadd: CellCoadd) -> None:
223 """The expanded report summarizes both tables instead of rendering them."""
224 report = minified_cell_coadd.provenance.describe()
225 assert report.type_name == "CoaddProvenance"
226 assert report_fields(report) == {
227 "instrument": "LSSTCam",
228 "physical_filter": "r_57",
229 "input images": "6 from 2 visits",
230 "day_obs": "20250520",
231 "cells": "3 with contributions",
232 "per cell": "2 input images",
233 }
234 assert list(report_fields(report)) == [
235 "instrument",
236 "physical_filter",
237 "input images",
238 "day_obs",
239 "cells",
240 "per cell",
241 ]
242 # Rich renders an astropy table as plain text and never consults its
243 # _repr_html_, so the report tabulates nothing.
244 assert not report.tables
245 assert "input images" in report._repr_html_()
248def test_coadd_provenance_report_counts_cells_against_bounds(minified_cell_coadd: CellCoadd) -> None:
249 """Given the image's cells, the report states coverage as a fraction."""
250 report = minified_cell_coadd.provenance._describe(bounds=minified_cell_coadd.bounds)
251 assert report_fields(report)["cells"] == "3 of 3 with contributions"
254def test_coadd_provenance_report_shows_partial_cell_coverage() -> None:
255 """Cells with image data but no contribution rows show up in the ratio."""
256 grid = CellGrid(bbox=Box.from_shape((30, 30)), cell_shape=YX(10, 10))
257 bounds = CellGridBounds(grid=grid, bbox=Box.factory[0:30, 0:30])
258 report = make_test_provenance()._describe(bounds=bounds)
259 assert report_fields(report)["cells"] == "3 of 9 with contributions"
262def test_coadd_provenance_report_ranges_over_nights_and_cells() -> None:
263 """Several visits, several nights and uneven cell coverage give the
264 ranged forms of each field.
265 """
266 fields = report_fields(make_test_provenance().describe())
267 assert fields["input images"] == "3 from 3 visits"
268 assert fields["day_obs"] == "20250520 - 20250521"
269 assert fields["cells"] == "3 with contributions"
270 assert fields["per cell"] == "1 - 3 input images (median 2)"
273def test_empty_coadd_provenance_cannot_be_constructed() -> None:
274 """Provenance is only meaningful when it records an input image."""
275 with pytest.raises(ValueError, match="at least one input"):
276 CoaddProvenance(
277 inputs=CoaddProvenance.make_empty_input_table(0),
278 contributions=CoaddProvenance.make_empty_contribution_table(0),
279 )
282def test_coadd_provenance_subset_without_contributions_is_none() -> None:
283 """Subsetting to cells that no image contributed to yields no provenance
284 rather than an empty one.
285 """
286 inputs = CoaddProvenance.make_empty_input_table(1)
287 inputs["instrument"] = "LSSTCam"
288 inputs["physical_filter"] = "r_57"
289 inputs["visit"] = [101]
290 inputs["detector"] = [1]
291 inputs["day_obs"] = [20250520]
292 contributions = CoaddProvenance.make_empty_contribution_table(1)
293 contributions["cell_i"] = [0]
294 contributions["cell_j"] = [0]
295 contributions["instrument"] = "LSSTCam"
296 contributions["visit"] = [101]
297 contributions["detector"] = [1]
298 provenance = CoaddProvenance(inputs=inputs, contributions=contributions)
300 populated = provenance[CellIJ(0, 0)]
301 assert populated is not None
302 assert len(populated.inputs) == 1
303 assert len(populated.contributions) == 1
305 assert provenance[CellIJ(1, 1)] is None
308def test_coadd_provenance_report_counts_one_input_image_in_the_singular() -> None:
309 """Pin the singular forms of the counted phrases.
311 One input image is what the standalone ``coadd_provenance`` file holds,
312 and what a single-cell subset of a coadd holds, so this is the case the
313 ``describe`` command line subcommand shows most often.
314 """
315 inputs = CoaddProvenance.make_empty_input_table(1)
316 inputs["instrument"] = "LSSTCam"
317 inputs["physical_filter"] = "r_57"
318 inputs["visit"] = [101]
319 inputs["detector"] = [1]
320 inputs["day_obs"] = [20250520]
321 provenance = CoaddProvenance(
322 inputs=inputs, contributions=CoaddProvenance.make_empty_contribution_table(1)
323 )
324 assert str(provenance) == "CoaddProvenance(1 input image)"
325 assert repr(provenance) == "<CoaddProvenance(1 input image)>"
326 fields = report_fields(provenance.describe())
327 assert fields["input images"] == "1 from 1 visit"
328 assert fields["per cell"] == "1 input image"
331def test_cell_coadd_report_accounts_for_every_component(minified_cell_coadd: CellCoadd) -> None:
332 """Nothing the coadd carries is absent from its report."""
333 report = minified_cell_coadd.describe()
334 fields = report_fields(report)
335 assert fields["mask_fractions"] == "rejected"
336 assert fields["noise_realizations"] == "1 image (dtype float32)"
337 assert fields["aperture_corrections"] == "3 fields"
338 # Provenance is a child, so it needs no field line of its own.
339 assert "provenance" not in fields
340 assert list(report.children) == [
341 "image",
342 "mask",
343 "variance",
344 "sky_projection",
345 "psf",
346 "provenance",
347 "backgrounds",
348 ]
349 provenance = report.children["provenance"]
350 assert provenance.type_name == "CoaddProvenance"
351 # The coadd passed its cells down, so coverage is stated as a fraction.
352 assert report_fields(provenance)["cells"] == "3 of 3 with contributions"
355def test_cell_coadd_report_counts_aperture_corrections_only(minified_cell_coadd: CellCoadd) -> None:
356 """A coadd can carry dozens of aperture corrections with long names, so
357 the report counts them and never names them, detail or not.
358 """
359 plain = report_fields(minified_cell_coadd.describe())["aperture_corrections"]
360 assert plain == "3 fields"
361 detailed = report_fields(minified_cell_coadd.describe(detail=True))["aperture_corrections"]
362 assert detailed == "3 fields" # No additional detail
365def test_cell_coadd_report_states_absent_provenance(minified_cell_coadd: CellCoadd) -> None:
366 """A coadd with no provenance says so, while components that are merely
367 empty stay out of the report.
368 """
369 bare = CellCoadd(
370 minified_cell_coadd.image,
371 mask=minified_cell_coadd.mask,
372 variance=minified_cell_coadd.variance,
373 sky_projection=minified_cell_coadd.sky_projection,
374 band=minified_cell_coadd.band,
375 psf=minified_cell_coadd.psf,
376 patch=minified_cell_coadd.patch,
377 )
378 report = bare.describe()
379 fields = report_fields(report)
380 assert fields["provenance"] == "none"
381 assert "provenance" not in report.children
382 assert "mask_fractions" not in fields
383 assert "noise_realizations" not in fields
384 assert "aperture_corrections" not in fields
387def test_cell_coadd_copy_preserves_absent_optional_state(minified_cell_coadd: CellCoadd) -> None:
388 """Copying must support coadds without patch or provenance records."""
389 minified_cell_coadd._patch = None
390 minified_cell_coadd._provenance = None
392 copied = minified_cell_coadd.copy()
394 assert copied._patch is None
395 assert copied._provenance is None
396 assert_masked_images_equal(copied, minified_cell_coadd)
397 with pytest.raises(AttributeError, match="no patch"):
398 copied.patch
399 with pytest.raises(AttributeError, match="no provenance"):
400 copied.provenance
403def test_cell_grid_patch_str_uses_clean_geometry() -> None:
404 """CellGrid, PatchDefinition and CellGridBounds str drop the
405 Interval/YX/Box/CellIJ wrappers.
407 The report renders field values with str, so these must use the compact
408 geometry forms rather than pydantic's default field-by-field repr.
409 """
410 grid = CellGrid(bbox=Box.from_shape((100, 200)), cell_shape=YX(10, 20))
411 assert str(grid) == "[y=0:100, x=0:200], cell_shape=(y=10, x=20)"
413 patch = PatchDefinition(id=73, index=YX(7, 3), inner_bbox=Box.factory[1:3, 2:4], cells=grid)
414 assert str(patch) == (
415 "id=73, index=(y=7, x=3), inner_bbox=[y=1:3, x=2:4], "
416 "cells=([y=0:100, x=0:200], cell_shape=(y=10, x=20))"
417 )
418 # repr stays as the pydantic default so it remains eval-ish and distinct.
419 assert "Interval(" in repr(patch)
420 assert "YX(" in repr(patch)
422 bounds = CellGridBounds(grid=grid, bbox=Box.factory[0:40, 0:60])
423 assert str(bounds) == "[y=0:40, x=0:60] in grid ([y=0:100, x=0:200], cell_shape=(y=10, x=20))"
424 bounds_missing = CellGridBounds(
425 grid=grid, bbox=Box.factory[0:40, 0:60], missing=frozenset({CellIJ(i=1, j=1), CellIJ(i=0, j=2)})
426 )
427 assert str(bounds_missing) == (
428 "[y=0:40, x=0:60] in grid ([y=0:100, x=0:200], cell_shape=(y=10, x=20)), "
429 "missing={(i=0, j=2), (i=1, j=1)}"
430 )
431 assert "Interval(" in repr(bounds_missing)
434def test_cell_shape_accepts_both_spellings() -> None:
435 """Verify both on-disk spellings of an XY pair still validate.
437 Shipped cell_coadd 1.0.0 files carry the array spelling, so this stays
438 readable until a major bump retires the as_shipped fixture. A focused
439 check here gives a two-line failure instead of a 52 KB fixture failing
440 to read.
441 """
442 tree_cls = class_for_schema("cell_psf")
443 assert tree_cls is not None
444 model = parameterize_tree(tree_cls, JsonRef)
445 fixture = json.loads(current_fixture_path(FIXTURE_DIR, "cell_psf").read_text())
446 assert fixture["bounds"]["grid"]["cell_shape"] == {"y": 4, "x": 4}
448 legacy = copy.deepcopy(fixture)
449 legacy["bounds"]["grid"]["cell_shape"] = [4, 4]
450 tree = model.model_validate(legacy)
451 assert isinstance(tree, CellPointSpreadFunctionSerializationModel)
452 assert tree.bounds.grid.cell_shape.y == 4
453 assert tree.bounds.grid.cell_shape.x == 4
456def test_cell_psf_rejects_index_below_bounds(minified_cell_coadd: CellCoadd) -> None:
457 """An index below the PSF bounds must not wrap around its array."""
458 start = minified_cell_coadd.psf.bounds.subgrid_start
459 with pytest.raises(BoundsError, match="out of bounds"):
460 minified_cell_coadd.psf[CellIJ(i=start.i - 1, j=start.j)]
463def test_from_legacy(legacy_test_data: _LegacyTestData) -> None:
464 """Test constructing a CellCoadd by converting a legacy
465 ``MultipleCellCoadd``.
466 """
467 assert legacy_test_data.cell_coadd.bounds.missing == {CellIJ(**DP2_COADD_MISSING_CELL)}
468 assert legacy_test_data.cell_coadd.bbox == Box.factory[12900:13500, 9600:10050]
469 compare_cell_coadd_to_legacy(
470 legacy_test_data.cell_coadd,
471 legacy_test_data.legacy_cell_coadd,
472 tract_bbox=legacy_test_data.tract_bbox,
473 plane_map=legacy_test_data.plane_map,
474 psf_points=legacy_test_data.make_psf_points(),
475 )
478def test_roundtrip(legacy_test_data: _LegacyTestData) -> None:
479 """Test that a CellCoadd roundtrips through FITS."""
480 with RoundtripFits(legacy_test_data.cell_coadd, "CellCoadd") as roundtrip:
481 # Check a subimage read (no component arg — does not trigger a skip).
482 subbox = Box.factory[
483 legacy_test_data.cell_coadd.bbox.y.start + 252 : legacy_test_data.cell_coadd.bbox.y.stop - 175,
484 legacy_test_data.cell_coadd.bbox.x.stop - 150 : legacy_test_data.cell_coadd.bbox.x.stop,
485 ]
486 subimage = roundtrip.get(bbox=subbox)
487 assert_masked_images_equal(subimage, legacy_test_data.cell_coadd[subbox], expect_view=False)
488 with roundtrip.inspect() as fits:
489 for extname in ["IMAGE", "MASK", "VARIANCE", "MASK_FRACTIONS/REJECTED"] + [
490 f"NOISE_REALIZATIONS/{n}" for n in range(len(legacy_test_data.cell_coadd.noise_realizations))
491 ]:
492 assert fits[extname].header["ZTILE1"] == legacy_test_data.cell_coadd.grid.cell_shape.x
493 assert fits[extname].header["ZTILE2"] == legacy_test_data.cell_coadd.grid.cell_shape.y
494 # Fixture self-consistency: bbox and missing-cell set are as expected.
495 assert legacy_test_data.cell_coadd.bounds.missing == {CellIJ(**DP2_COADD_MISSING_CELL)}
496 assert legacy_test_data.cell_coadd.bbox == Box.factory[12900:13500, 9600:10050]
497 # Full round-trip fidelity.
498 assert_cell_coadds_equal(roundtrip.result, legacy_test_data.cell_coadd, expect_view=False)
499 compare_cell_coadd_to_legacy(
500 roundtrip.result,
501 legacy_test_data.legacy_cell_coadd,
502 tract_bbox=legacy_test_data.tract_bbox,
503 plane_map=legacy_test_data.plane_map,
504 psf_points=legacy_test_data.make_psf_points(),
505 )
508def test_roundtrip_components(legacy_test_data: _LegacyTestData) -> None:
509 """Test component and subimage reads.
511 This test will be skipped if `lsst.daf.butler` is not available instead of
512 falling back to non-butler I/O, which is why we don't want to merge it
513 with `test_roundtrip`.
514 """
515 with RoundtripFits(legacy_test_data.cell_coadd, "CellCoadd") as roundtrip:
516 subbox = make_subbox(legacy_test_data.cell_coadd.bbox)
517 subpsf = roundtrip.get("psf", bbox=subbox)
518 assert subpsf.bounds.bbox == Box(
519 y=Interval.factory[
520 legacy_test_data.cell_coadd.bbox.y.start + 150 : legacy_test_data.cell_coadd.bbox.y.stop - 150
521 ],
522 x=subbox.x,
523 )
524 assert_psfs_equal(
525 subpsf,
526 legacy_test_data.cell_coadd.psf,
527 points=legacy_test_data.make_psf_points(subbox),
528 )
529 assert roundtrip.get("bbox") == legacy_test_data.cell_coadd.bbox
530 alternates = {
531 k: roundtrip.get(k)
532 for k in [
533 "sky_projection",
534 "image",
535 "mask",
536 "variance",
537 "masked_image",
538 "psf",
539 "aperture_corrections",
540 "provenance",
541 "backgrounds",
542 "bbox",
543 ]
544 }
545 # Read all the components at once.
546 all_components = roundtrip.get("components")
547 assert set(all_components) == set(alternates) - {"masked_image"}
548 assert all_components["bbox"] == alternates["bbox"]
549 assert_psfs_equal(all_components["psf"], alternates["psf"])
550 assert_images_equal(all_components["image"], alternates["image"])
552 backgrounds = roundtrip.get("backgrounds")
553 assert backgrounds.keys() == set()
554 assert backgrounds.subtracted is None
556 compare_cell_coadd_to_legacy(
557 roundtrip.result,
558 legacy_test_data.legacy_cell_coadd,
559 tract_bbox=legacy_test_data.tract_bbox,
560 plane_map=legacy_test_data.plane_map,
561 alternates=alternates,
562 psf_points=legacy_test_data.make_psf_points(),
563 )
566def test_fits_compression(legacy_test_data: _LegacyTestData) -> None:
567 """Test lossy FITS compression produces the expected headers."""
568 with RoundtripFits(
569 legacy_test_data.cell_coadd,
570 storage_class="CellCoadd",
571 recipe="lossy16",
572 compression_options={
573 "image": FitsCompressionOptions.LOSSY,
574 "variance": FitsCompressionOptions.LOSSY,
575 },
576 ) as roundtrip:
577 with roundtrip.inspect() as fits:
578 for extname in ["IMAGE", "MASK", "VARIANCE", "MASK_FRACTIONS/REJECTED"] + [
579 f"NOISE_REALIZATIONS/{n}" for n in range(len(legacy_test_data.cell_coadd.noise_realizations))
580 ]:
581 assert fits[extname].header["ZTILE1"] == legacy_test_data.cell_coadd.grid.cell_shape.x
582 assert fits[extname].header["ZTILE2"] == legacy_test_data.cell_coadd.grid.cell_shape.y
583 if extname == "MASK" or extname.startswith("MASK_FRACTIONS"):
584 assert fits[extname].header["ZCMPTYPE"] == "GZIP_2"
585 else:
586 assert fits[extname].header["ZCMPTYPE"] == "RICE_1"
587 assert fits[extname].header["ZQUANTIZ"] == "SUBTRACTIVE_DITHER_2"
590def test_json_roundtrip(legacy_test_data: _LegacyTestData) -> None:
591 """Verify a CellCoadd round-trips correctly through the JSON archive."""
592 with RoundtripJson(legacy_test_data.cell_coadd) as roundtrip:
593 pass
594 assert_cell_coadds_equal(roundtrip.result, legacy_test_data.cell_coadd, expect_view=False)
597def test_to_legacy_cell_coadd(legacy_test_data: _LegacyTestData) -> None:
598 """Verify converting a CellCoadd back into a legacy MultipleCellCoadd."""
599 legacy_cell_coadd = legacy_test_data.cell_coadd.to_legacy_cell_coadd()
600 compare_cell_coadd_to_legacy(
601 legacy_test_data.cell_coadd,
602 legacy_cell_coadd,
603 tract_bbox=legacy_test_data.tract_bbox,
604 plane_map=legacy_test_data.plane_map,
605 psf_points=legacy_test_data.make_psf_points(),
606 )
607 with pytest.raises(
608 ValueError, match="MultipleCellCoadd requires its bounding box to lie on the cell grid."
609 ):
610 legacy_test_data.cell_coadd[make_subbox(legacy_test_data.cell_coadd.bbox)].to_legacy_cell_coadd()
613@skip_no_legacy
614def test_to_legacy(legacy_test_data: _LegacyTestData) -> None:
615 """Test converting a CellCoadd back into a legacy Exposure."""
616 legacy_exposure = legacy_test_data.cell_coadd.to_legacy()
617 assert legacy_exposure.getFilter().bandLabel == legacy_test_data.cell_coadd.band
618 assert Box.from_legacy(legacy_exposure.getBBox()) == legacy_test_data.cell_coadd.bbox
619 compare_masked_image_to_legacy(
620 legacy_test_data.cell_coadd,
621 legacy_exposure.maskedImage,
622 plane_map=legacy_test_data.plane_map,
623 expect_view=True,
624 )
625 compare_psf_to_legacy(
626 legacy_test_data.cell_coadd.psf,
627 legacy_exposure.getPsf(),
628 points=legacy_test_data.make_psf_points(),
629 expect_legacy_raise_on_out_of_bounds=True,
630 )
631 compare_sky_projection_to_legacy_wcs(
632 legacy_test_data.cell_coadd.sky_projection,
633 legacy_exposure.getWcs(),
634 legacy_test_data.cell_coadd.sky_projection.pixel_frame,
635 subimage_bbox=legacy_test_data.cell_coadd.bbox,
636 is_fits=True,
637 )
638 subbox = make_subbox(legacy_test_data.cell_coadd.bbox)
639 compare_masked_image_to_legacy(
640 legacy_test_data.cell_coadd[subbox],
641 legacy_test_data.cell_coadd[subbox].to_legacy().maskedImage,
642 plane_map=legacy_test_data.plane_map,
643 expect_view=True,
644 )
647@skip_no_h5py
648def test_ndf_roundtrip(legacy_test_data: _LegacyTestData) -> None:
649 """Test that CellCoadd round-trips through NDF."""
650 with RoundtripNdf(legacy_test_data.cell_coadd, "CellCoadd") as roundtrip:
651 assert_cell_coadds_equal(roundtrip.result, legacy_test_data.cell_coadd, expect_view=False)
654# The float32 image pixels are of order 10 nJy, as is the test background
655# below, so background arithmetic is only reproducible to a few ULPs at that
656# scale.
657BACKGROUND_ATOL = 1e-5
660def _add_gradient_background(coadd: CellCoadd, name: str = "pretty") -> ChebyshevField:
661 """Add a non-constant background over the coadd's full bbox and return the
662 field that was added.
663 """
664 field = ChebyshevField(
665 coadd.bbox,
666 np.array([[10.0, 2.0, 0.5], [1.0, 0.75, 0.0], [0.25, 0.0, 0.0]]),
667 unit=coadd.unit,
668 )
669 coadd.backgrounds.add(name, field, description="Gradient background for tests.")
670 return field
673def _make_background_subbox(coadd: CellCoadd) -> Box:
674 """Make a box that trims a different number of pixels from each side of the
675 coadd, so a background rendered over the wrong box cannot match by chance.
676 """
677 return Box.factory[
678 coadd.bbox.y.start + 2 : coadd.bbox.y.stop - 3,
679 coadd.bbox.x.start + 1 : coadd.bbox.x.stop - 2,
680 ]
683def test_apply_background_subimage(minified_cell_coadd: CellCoadd) -> None:
684 """Test that applying a background to a subimage subtracts only the
685 portion of the background model that overlaps the subimage.
686 """
687 field = _add_gradient_background(minified_cell_coadd)
688 subbox = _make_background_subbox(minified_cell_coadd)
689 subimage = minified_cell_coadd[subbox]
690 original = subimage.image.array.copy()
692 subimage.apply_background("pretty")
694 assert subimage.backgrounds.subtracted is not None
695 assert subimage.backgrounds.subtracted.name == "pretty"
696 assert subimage.image.bbox == subbox
697 expected = field.render(subbox, dtype=subimage.image.array.dtype).array
698 np.testing.assert_allclose(original - subimage.image.array, expected, atol=BACKGROUND_ATOL)
701def test_restore_background_subimage(minified_cell_coadd: CellCoadd) -> None:
702 """Test that restoring the original background on a subimage adds back
703 only the portion of the background model that overlaps the subimage.
704 """
705 _add_gradient_background(minified_cell_coadd)
706 subbox = _make_background_subbox(minified_cell_coadd)
707 subimage = minified_cell_coadd[subbox]
708 original = subimage.image.array.copy()
709 subimage.apply_background("pretty")
711 subimage.apply_background(None)
713 assert subimage.backgrounds.subtracted is None
714 np.testing.assert_allclose(subimage.image.array, original, atol=BACKGROUND_ATOL)
717def test_failed_background_switch_preserves_pixels_and_state(
718 minified_cell_coadd: CellCoadd,
719) -> None:
720 """An invalid replacement must not restore the current background."""
721 _add_gradient_background(minified_cell_coadd)
722 minified_cell_coadd.apply_background("pretty")
723 before = minified_cell_coadd.image.array.copy()
725 with pytest.raises(KeyError):
726 minified_cell_coadd.apply_background("missing")
728 np.testing.assert_array_equal(minified_cell_coadd.image.array, before)
729 assert minified_cell_coadd.backgrounds.subtracted is not None
730 assert minified_cell_coadd.backgrounds.subtracted.name == "pretty"
733def test_switch_background_applies_replacement(minified_cell_coadd: CellCoadd) -> None:
734 """Switching models restores the old background and subtracts the new."""
735 field = _add_gradient_background(minified_cell_coadd)
736 minified_cell_coadd.backgrounds.add("other", field * 2.0)
737 original = minified_cell_coadd.image.array.copy()
738 minified_cell_coadd.apply_background("pretty")
740 minified_cell_coadd.apply_background("other")
742 expected = (
743 original
744 - (field * 2.0).render(minified_cell_coadd.bbox, dtype=minified_cell_coadd.image.array.dtype).array
745 )
746 np.testing.assert_allclose(minified_cell_coadd.image.array, expected, atol=BACKGROUND_ATOL)
747 assert minified_cell_coadd.backgrounds.subtracted is not None
748 assert minified_cell_coadd.backgrounds.subtracted.name == "other"
751def test_apply_background_after_bounded_read(minified_cell_coadd: CellCoadd) -> None:
752 """Test that a background can be applied to a coadd read with a bbox
753 parameter.
754 """
755 field = _add_gradient_background(minified_cell_coadd)
756 subbox = _make_background_subbox(minified_cell_coadd)
757 with RoundtripJson(minified_cell_coadd, "CellCoadd") as roundtrip:
758 subimage = roundtrip.get(bbox=subbox)
759 original = subimage.image.array.copy()
761 subimage.apply_background("pretty")
763 assert subimage.image.bbox == subbox
764 expected = field.render(subbox, dtype=subimage.image.array.dtype).array
765 np.testing.assert_allclose(original - subimage.image.array, expected, atol=BACKGROUND_ATOL)
768def test_cell_grid_bounds_contains_broadcasting(minified_cell_coadd: CellCoadd) -> None:
769 """Test that CellGridBounds.contains broadcasts like a numpy ufunc."""
770 assert minified_cell_coadd.bounds.missing, "fixture should retain a missing cell"
771 check_bounds_contains_broadcasting(minified_cell_coadd.bounds)
774def test_intersection_bounds_contains_broadcasting(minified_cell_coadd: CellCoadd) -> None:
775 """Test that IntersectionBounds.contains broadcasts like a numpy ufunc."""
776 # Clip the CellGridBounds with a Box offset by 1 pixel on each side so it
777 # does not snap to any cell boundary, forcing a lazy IntersectionBounds.
778 bounds = minified_cell_coadd.bounds
779 clip = Box.factory[
780 bounds.bbox.y.start + 1 : bounds.bbox.y.stop - 1,
781 bounds.bbox.x.start + 1 : bounds.bbox.x.stop - 1,
782 ]
783 check_bounds_contains_broadcasting(bounds.intersection(clip))