Coverage for tests/test_transforms.py: 49%
180 statements
« prev ^ index » next coverage.py v7.14.3, created at 2026-06-29 15:33 -0700
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 functools
16import os
17import unittest
18from typing import Any, ClassVar
20import astropy.units as u
21import numpy as np
22import pydantic
24from lsst.images import (
25 ICRS,
26 Box,
27 CameraFrameSet,
28 CameraFrameSetSerializationModel,
29 DetectorFrame,
30 FocalPlaneFrame,
31 SkyProjection,
32 Transform,
33 TransformSerializationModel,
34)
35from lsst.images._transforms import _ast as astshim
36from lsst.images.fits import PointerModel
37from lsst.images.serialization import ArchiveTree, InputArchive, JsonRef, OutputArchive
38from lsst.images.tests import (
39 DP2_VISIT_DETECTOR_DATA_ID,
40 RoundtripFits,
41 RoundtripJson,
42 check_transform,
43 compare_sky_projection_to_legacy_wcs,
44 legacy_points_to_xy_array,
45)
47DATA_DIR = os.environ.get("TESTDATA_IMAGES_DIR", None)
50class TransformTestCase(unittest.TestCase):
51 """Tests for the Transform, SkyProjection, and FrameSet classes."""
53 def test_identity(self) -> None:
54 """Test an identity transform."""
55 frame = DetectorFrame(**DP2_VISIT_DETECTOR_DATA_ID, bbox=Box.factory[:5, :4])
56 xy = frame.bbox.meshgrid().map(np.ravel)
57 identity = Transform.identity(frame)
58 check_transform(self, identity, xy, xy, frame, frame)
59 self.assertEqual(identity.decompose(), [])
60 with RoundtripJson(self, identity) as roundtrip:
61 pass
62 check_transform(self, roundtrip.result, xy, xy, frame, frame)
64 def test_transform_equality(self) -> None:
65 """Test ``Transform.__eq__`` across all of its comparison branches."""
66 pixel_frame = DetectorFrame(**DP2_VISIT_DETECTOR_DATA_ID, bbox=Box.factory[:5, :4])
67 focal_plane = FocalPlaneFrame(instrument="LSSTCam", visit=1, unit=u.mm)
68 # A distinct frame for the in-frame and out-frame branches.
69 alt_frame = DetectorFrame(instrument="LSSTCam", visit=1, detector=12, bbox=Box.factory[:5, :4])
70 in_bounds = Box.factory[:5, :4]
71 out_bounds = Box.factory[:10, :8]
73 def make(
74 *,
75 in_frame: Any = pixel_frame,
76 out_frame: Any = focal_plane,
77 ast_mapping: astshim.Mapping | None = None,
78 in_bounds_: Box | None = in_bounds,
79 out_bounds_: Box | None = out_bounds,
80 components: Any = (),
81 ) -> Transform[Any, Any]:
82 return Transform(
83 in_frame,
84 out_frame,
85 ast_mapping if ast_mapping is not None else astshim.UnitMap(2),
86 in_bounds=in_bounds_,
87 out_bounds=out_bounds_,
88 components=components,
89 )
91 base = make()
93 # Identity short-circuit: an object is always equal to itself.
94 self.assertEqual(base, base)
96 # Two independently constructed but equivalent transforms are equal,
97 # and equality is symmetric.
98 self.assertEqual(base, make())
99 self.assertEqual(make(), base)
101 # Comparison against a non-Transform yields NotImplemented, so Python
102 # falls back to identity: the objects are unequal and ``!=`` is True.
103 self.assertFalse(base == "not a transform")
104 self.assertTrue(base != "not a transform")
105 self.assertNotEqual(base, None)
106 self.assertNotEqual(base, 42)
108 # Each remaining branch differs from ``base`` in exactly one attribute.
109 self.assertNotEqual(base, make(ast_mapping=astshim.ShiftMap([1.0, 2.0])))
110 self.assertNotEqual(base, make(in_bounds_=out_bounds))
111 self.assertNotEqual(base, make(out_bounds_=in_bounds))
112 self.assertNotEqual(base, make(in_frame=alt_frame))
113 self.assertNotEqual(base, make(out_frame=alt_frame))
114 self.assertNotEqual(base, make(components=[Transform.identity(alt_frame)]))
116 def test_sky_projection_equality(self) -> None:
117 """Test ``SkyProjection.__eq__`` across all of its comparison
118 branches.
119 """
120 pixel_frame = DetectorFrame(**DP2_VISIT_DETECTOR_DATA_ID, bbox=Box.factory[:5, :4])
122 # Check the two failure modes.
123 with self.assertRaises(ValueError):
124 SkyProjection(Transform(ICRS, ICRS, astshim.UnitMap(2)))
126 with self.assertRaises(ValueError):
127 SkyProjection(Transform(pixel_frame, pixel_frame, astshim.UnitMap(2)))
129 def make_pixel_to_sky(ast_mapping: astshim.Mapping | None = None) -> Transform[Any, Any]:
130 mapping = ast_mapping if ast_mapping is not None else astshim.UnitMap(2)
131 return Transform(pixel_frame, ICRS, mapping)
133 base = SkyProjection(make_pixel_to_sky())
135 # Identity short-circuit: an object is always equal to itself.
136 self.assertEqual(base, base)
138 # Two independently constructed but equivalent projections are equal.
139 self.assertEqual(base, SkyProjection(make_pixel_to_sky()))
141 # Comparison against a non-SkyProjection yields NotImplemented, so
142 # Python falls back to identity.
143 self.assertFalse(base == "not a projection")
144 self.assertTrue(base != "not a projection")
145 self.assertNotEqual(base, None)
147 # Differ only in the pixel-to-sky transform.
148 self.assertNotEqual(base, SkyProjection(make_pixel_to_sky(astshim.ShiftMap([1.0, 2.0]))))
150 # The fits_approximation branch: absent on ``base`` but present here.
151 with_approx = SkyProjection(
152 make_pixel_to_sky(), fits_approximation=make_pixel_to_sky(astshim.ShiftMap([0.1, 0.2]))
153 )
154 self.assertNotEqual(base, with_approx)
156 # Equal pixel-to-sky and equal fits_approximations are equal.
157 with_approx_again = SkyProjection(
158 make_pixel_to_sky(), fits_approximation=make_pixel_to_sky(astshim.ShiftMap([0.1, 0.2]))
159 )
160 self.assertEqual(with_approx, with_approx_again)
162 # Same pixel-to-sky transform but a different fits_approximation.
163 other_approx = SkyProjection(
164 make_pixel_to_sky(), fits_approximation=make_pixel_to_sky(astshim.ShiftMap([0.3, 0.4]))
165 )
166 self.assertNotEqual(with_approx, other_approx)
168 @unittest.skipUnless(DATA_DIR is not None, "TESTDATA_IMAGES_DIR is not in the environment.")
169 def test_camera(self) -> None:
170 """Test that we can:
172 - make a CameraFrameSet from the AST representation returned by afw;
173 - transform points and get the same result as afw;
174 - round-trip the CameraFrameSet through FITS serialization and still
175 do all of that;
176 - also roundtrip a Transform that can be obtained from the
177 CameraFrameSet, by referencing the mappings in the frame set.
179 This test is skipped if legacy modules cannot be imported.
181 This test provides coverage for the archive system's pointer and
182 frame-set reference machinery.
183 """
184 try:
185 from lsst.afw.cameraGeom import Camera
186 except ImportError:
187 raise unittest.SkipTest("'lsst.afw.cameraGeom' could not be imported.") from None
188 assert DATA_DIR is not None, "Guaranteed by decorator."
189 filename = os.path.join(DATA_DIR, "dp2", "legacy", "camera.fits")
190 legacy_camera = Camera.readFits(filename)
191 frame_set = CameraFrameSet.from_legacy(legacy_camera)
192 detector_id: int = DP2_VISIT_DETECTOR_DATA_ID["detector"]
193 self.compare_to_legacy_camera(legacy_camera, frame_set)
194 test_holder = FrameSetTestHolder(
195 frames=frame_set,
196 pixels_to_fp=frame_set[frame_set.detector(detector_id), frame_set.focal_plane()],
197 )
198 with RoundtripFits(self, test_holder) as roundtrip1:
199 self.assertEqual(len(roundtrip1.serialized.pixels_to_fp.frames), 2)
200 self.assertEqual(len(roundtrip1.serialized.pixels_to_fp.bounds), 2)
201 self.assertEqual(len(roundtrip1.serialized.pixels_to_fp.mappings), 1)
202 # Instead of storing the AST mapping directly, we should have
203 # stored a reference to the frame set:
204 self.assertIsInstance(roundtrip1.serialized.pixels_to_fp.mappings[0], PointerModel)
205 self.compare_to_legacy_camera(legacy_camera, roundtrip1.result.frames)
206 self.assertEqual(roundtrip1.result.pixels_to_fp.in_frame, frame_set.detector(detector_id))
207 self.assertEqual(roundtrip1.result.pixels_to_fp.out_frame, frame_set.focal_plane())
208 self.assertEqual(
209 roundtrip1.result.pixels_to_fp._ast_mapping.simplified().show(),
210 test_holder.pixels_to_fp._ast_mapping.simplified().show(),
211 )
212 with RoundtripJson(self, test_holder) as roundtrip2:
213 self.assertEqual(len(roundtrip2.serialized.pixels_to_fp.frames), 2)
214 self.assertEqual(len(roundtrip2.serialized.pixels_to_fp.bounds), 2)
215 self.assertEqual(len(roundtrip2.serialized.pixels_to_fp.mappings), 1)
216 # Instead of storing the AST mapping directly, we should have
217 # stored a reference to the frame set:
218 self.assertIsInstance(roundtrip2.serialized.pixels_to_fp.mappings[0], JsonRef)
219 raw_data = roundtrip2.inspect()
220 self.assertEqual(len(raw_data["indirect"]), 1)
221 self.assertEqual(raw_data["frames"], {"$ref": "#/indirect/0"})
222 self.compare_to_legacy_camera(legacy_camera, roundtrip2.result.frames)
223 self.assertEqual(roundtrip2.result.pixels_to_fp.in_frame, frame_set.detector(detector_id))
224 self.assertEqual(roundtrip2.result.pixels_to_fp.out_frame, frame_set.focal_plane())
225 self.assertEqual(
226 roundtrip2.result.pixels_to_fp._ast_mapping.simplified().show(),
227 test_holder.pixels_to_fp._ast_mapping.simplified().show(),
228 )
230 def compare_to_legacy_camera(self, legacy_camera: Any, frame_set: CameraFrameSet) -> None:
231 """Test the transforms extracted from a CameraFrameSet against the
232 legacy lsst.afw.cameraGeom implementations.
233 """
234 from lsst.afw.cameraGeom import FIELD_ANGLE, FOCAL_PLANE, PIXELS
235 from lsst.geom import Point2D
237 legacy_detector = legacy_camera[16]
238 pixel_legacy_points = [Point2D(50.0, 60.0), Point2D(801.2, 322.8), Point2D(33.5, 22.1)]
239 fp_legacy_points = [legacy_detector.transform(p, PIXELS, FOCAL_PLANE) for p in pixel_legacy_points]
240 fa_legacy_points = [legacy_detector.transform(p, PIXELS, FIELD_ANGLE) for p in pixel_legacy_points]
241 pixel_xy_array = legacy_points_to_xy_array(pixel_legacy_points)
242 fp_xy_array = legacy_points_to_xy_array(fp_legacy_points)
243 fa_xy_array = legacy_points_to_xy_array(fa_legacy_points)
244 # Test transforms extracted directly from the frame set.
245 pixel_to_fp = frame_set[frame_set.detector(16), frame_set.focal_plane()]
246 check_transform(
247 self, pixel_to_fp, pixel_xy_array, fp_xy_array, frame_set.detector(16), frame_set.focal_plane()
248 )
249 pixel_to_fa = frame_set[frame_set.detector(16), frame_set.field_angle()]
250 check_transform(
251 self, pixel_to_fa, pixel_xy_array, fa_xy_array, frame_set.detector(16), frame_set.field_angle()
252 )
253 fp_to_fa = frame_set[frame_set.focal_plane(), frame_set.field_angle()]
254 check_transform(
255 self, fp_to_fa, fp_xy_array, fa_xy_array, frame_set.focal_plane(), frame_set.field_angle()
256 )
257 # Test a composition.
258 pixel_to_fa_indirect = pixel_to_fp.then(fp_to_fa)
259 check_transform(
260 self,
261 pixel_to_fa_indirect,
262 pixel_xy_array,
263 fa_xy_array,
264 frame_set.detector(16),
265 frame_set.field_angle(),
266 )
267 pixel_to_fp_d, fp_to_fa_d = pixel_to_fa_indirect.decompose()
268 check_transform(
269 self, pixel_to_fp_d, pixel_xy_array, fp_xy_array, frame_set.detector(16), frame_set.focal_plane()
270 )
271 check_transform(
272 self, fp_to_fa_d, fp_xy_array, fa_xy_array, frame_set.focal_plane(), frame_set.field_angle()
273 )
274 fa_to_fp_d, fp_to_pixel_d = pixel_to_fa_indirect.inverted().decompose()
275 check_transform(
276 self, fa_to_fp_d, fa_xy_array, fp_xy_array, frame_set.field_angle(), frame_set.focal_plane()
277 )
278 check_transform(
279 self, fp_to_pixel_d, fp_xy_array, pixel_xy_array, frame_set.focal_plane(), frame_set.detector(16)
280 )
282 @unittest.skipUnless(DATA_DIR is not None, "TESTDATA_IMAGES_DIR is not in the environment.")
283 def test_detector_wcs(self) -> None:
284 """Test the Transform/SkyProjection representation of a detector
285 WCS.
286 """
287 try:
288 from lsst.afw.image import ExposureFitsReader
289 except ImportError:
290 raise unittest.SkipTest("'lsst.afw.image' could not be imported.") from None
291 assert DATA_DIR is not None, "Guaranteed by decorator."
292 filename = os.path.join(DATA_DIR, "dp2", "legacy", "visit_image.fits")
293 reader = ExposureFitsReader(filename)
294 legacy_wcs = reader.readWcs()
295 wcs_bbox = Box.from_legacy(reader.readDetector().getBBox())
296 subimage_bbox = Box.from_legacy(reader.readBBox())
297 detector_frame = DetectorFrame(**DP2_VISIT_DETECTOR_DATA_ID, bbox=wcs_bbox)
298 sky_projection = SkyProjection.from_legacy(legacy_wcs, detector_frame)
299 assert sky_projection.fits_approximation is not None
300 compare_sky_projection_to_legacy_wcs(self, sky_projection, legacy_wcs, detector_frame, subimage_bbox)
301 # When we convert from a legacy SkyWcs, the internal AST Mapping needs
302 # to really be an AST FrameSet in order to be able to convert back.
303 self.assertIn("Begin FrameSet", sky_projection.show())
304 compare_sky_projection_to_legacy_wcs(
305 self, sky_projection, sky_projection.to_legacy(), detector_frame, subimage_bbox
306 )
307 self.assertIn("Begin FrameSet", sky_projection.fits_approximation.show())
308 compare_sky_projection_to_legacy_wcs(
309 self,
310 sky_projection.fits_approximation,
311 sky_projection.fits_approximation.to_legacy(),
312 detector_frame,
313 subimage_bbox,
314 is_fits=True,
315 )
316 with RoundtripJson(self, sky_projection, "SkyProjection") as roundtrip:
317 pass
318 compare_sky_projection_to_legacy_wcs(
319 self, roundtrip.result, legacy_wcs, detector_frame, subimage_bbox
320 )
321 # The AST FrameSet-ness needs to propagate through serialization.
322 self.assertIn("Begin FrameSet", roundtrip.result.show())
323 compare_sky_projection_to_legacy_wcs(
324 self, sky_projection, roundtrip.result.to_legacy(), detector_frame, subimage_bbox
325 )
326 with RoundtripJson(self, sky_projection.fits_approximation, "SkyProjection") as roundtrip:
327 pass
328 compare_sky_projection_to_legacy_wcs(
329 self,
330 roundtrip.result,
331 legacy_wcs.getFitsApproximation(),
332 detector_frame,
333 subimage_bbox,
334 is_fits=True,
335 )
336 self.assertIn("Begin FrameSet", roundtrip.result.show())
337 compare_sky_projection_to_legacy_wcs(
338 self,
339 sky_projection.fits_approximation,
340 roundtrip.result.to_legacy(),
341 detector_frame,
342 subimage_bbox,
343 is_fits=True,
344 )
347@dataclasses.dataclass
348class FrameSetTestHolder:
349 """A top-level object that holds a CameraFrameSet and a transform
350 extracted from it, for testing archive pointers and frame set references.
351 """
353 frames: CameraFrameSet
354 pixels_to_fp: Transform[DetectorFrame, FocalPlaneFrame]
356 def serialize[P: pydantic.BaseModel](self, archive: OutputArchive[P]) -> FrameSetTestHolderModel[P]:
357 frames_model = archive.serialize_frame_set(
358 "frames", self.frames, self.frames.serialize, key=id(self.frames)
359 )
360 pixels_to_fp_model = archive.serialize_direct(
361 "pixels_to_fp", functools.partial(self.pixels_to_fp.serialize, use_frame_sets=True)
362 )
363 return FrameSetTestHolderModel[P](frames=frames_model, pixels_to_fp=pixels_to_fp_model)
365 @staticmethod
366 def _get_archive_tree_type[P: pydantic.BaseModel](
367 pointer_type: type[P],
368 ) -> type[FrameSetTestHolderModel[P]]:
369 return FrameSetTestHolderModel[pointer_type] # type: ignore
372class FrameSetTestHolderModel[P: pydantic.BaseModel](ArchiveTree):
373 """The serialization model for FrameSetTestHolder."""
375 SCHEMA_NAME: ClassVar[str] = "_test_frame_set_holder"
376 SCHEMA_VERSION: ClassVar[str] = "1.0.0"
377 MIN_READ_VERSION: ClassVar[int] = 1
378 PUBLIC_TYPE: ClassVar[type] = FrameSetTestHolder
380 frames: CameraFrameSetSerializationModel | P
381 pixels_to_fp: TransformSerializationModel[P]
383 def deserialize(self, archive: InputArchive[Any]) -> FrameSetTestHolder:
384 assert not isinstance(self.frames, CameraFrameSetSerializationModel), "Archive pointer expected."
385 frames = archive.deserialize_pointer(
386 self.frames, CameraFrameSetSerializationModel, CameraFrameSetSerializationModel.deserialize
387 )
388 pixels_to_fp = self.pixels_to_fp.deserialize(archive)
389 return FrameSetTestHolder(frames, pixels_to_fp)
392if __name__ == "__main__":
393 unittest.main()