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