Coverage for python / lsst / images / tests / _checks.py: 11%

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. 

11 

12from __future__ import annotations 

13 

14__all__ = ( 

15 "arrays_to_legacy_points", 

16 "assert_close", 

17 "assert_equal_allow_nan", 

18 "assert_images_equal", 

19 "assert_masked_images_equal", 

20 "assert_masks_equal", 

21 "assert_projections_equal", 

22 "assert_psfs_equal", 

23 "check_astropy_wcs_interface", 

24 "check_projection", 

25 "check_transform", 

26 "compare_amplifier_to_legacy", 

27 "compare_aperture_corrections_to_legacy", 

28 "compare_cell_coadd_to_legacy", 

29 "compare_detector_to_legacy", 

30 "compare_field_to_legacy", 

31 "compare_image_to_legacy", 

32 "compare_mask_to_legacy", 

33 "compare_masked_image_to_legacy", 

34 "compare_observation_summary_stats_to_legacy", 

35 "compare_projection_to_legacy_wcs", 

36 "compare_psf_to_legacy", 

37 "compare_visit_image_to_legacy", 

38 "legacy_coords_to_astropy", 

39 "legacy_points_to_xy_array", 

40) 

41 

42import dataclasses 

43import math 

44import unittest 

45from collections.abc import Mapping 

46from typing import TYPE_CHECKING, Any, Literal, cast 

47 

48import astropy.units as u 

49import astropy.wcs.wcsapi 

50import numpy as np 

51from astropy.coordinates import SkyCoord 

52 

53from .._geom import XY, YX, BoundsError, Box 

54from .._image import Image 

55from .._mask import Mask, MaskPlane 

56from .._masked_image import MaskedImage 

57from .._observation_summary_stats import ObservationSummaryStats 

58from .._transforms import DetectorFrame, Frame, Projection, SkyFrame, TractFrame, Transform 

59from .._visit_image import VisitImage 

60from ..aperture_corrections import ApertureCorrectionMap 

61from ..cameras import Amplifier, Detector, DetectorType, ReadoutCorner 

62from ..cells import CellCoadd, CellIJ, CoaddProvenance 

63from ..fields import BaseField 

64from ..psfs import PointSpreadFunction 

65 

66if TYPE_CHECKING: 

67 try: 

68 from lsst.cell_coadds import MultipleCellCoadd 

69 except ImportError: 

70 type MultipleCellCoadd = Any # type: ignore[no-redef] 

71 

72 

73def assert_close( 

74 tc: unittest.TestCase, 

75 a: np.ndarray | u.Quantity | float, 

76 b: np.ndarray | u.Quantity | float, 

77 **kwargs: Any, 

78) -> None: 

79 """Test that two arrays, floats, or quantities are almost equal. 

80 

81 Parameters 

82 ---------- 

83 tc 

84 Test case object with assert methods to use. 

85 a 

86 Array, float, or quantity to compare. 

87 b 

88 Array, float, or quantity to compare. 

89 **kwargs 

90 Forwarded to `astropy.units.allclose`. 

91 """ 

92 tc.assertTrue(u.allclose(a, b, **kwargs), msg=f"{a} != {b}") 

93 

94 

95def assert_equal_allow_nan(tc: unittest.TestCase, a: float, b: float) -> None: 

96 """Test that two floating point values are equal, with nan == nan.""" 

97 try: 

98 tc.assertEqual(a, b) 

99 except AssertionError: 

100 if not (math.isnan(a) and math.isnan(b)): 

101 raise 

102 

103 

104def assert_images_equal( 

105 tc: unittest.TestCase, 

106 a: Image, 

107 b: Image, 

108 *, 

109 rtol: float = 0.0, 

110 atol: float = 0.0, 

111 expect_view: bool | Literal["array"] | None = None, 

112) -> None: 

113 """Assert that two images are equal or nearly equal.""" 

114 tc.assertEqual(a.bbox, b.bbox) 

115 tc.assertEqual(a.unit, b.unit) 

116 assert_projections_equal(tc, a.projection, b.projection) 

117 if expect_view is not None: 

118 tc.assertEqual(np.may_share_memory(a.array, b.array), bool(expect_view)) 

119 if expect_view == "array": 

120 tc.assertEqual(a.metadata, b.metadata) 

121 else: 

122 tc.assertEqual(a.metadata is b.metadata, expect_view) 

123 if not expect_view: 

124 assert_close(tc, a.array, b.array, atol=atol, rtol=rtol) 

125 tc.assertEqual(a.metadata, b.metadata) 

126 

127 

128def assert_masks_equal(tc: unittest.TestCase, a: Mask, b: Mask) -> None: 

129 """Assert that two masks are equal or nearly equal.""" 

130 tc.assertEqual(a.bbox, b.bbox) 

131 tc.assertEqual(a.schema, b.schema) 

132 tc.assertEqual(a.metadata, b.metadata) 

133 assert_projections_equal(tc, a.projection, b.projection) 

134 np.testing.assert_array_equal(a.array, b.array) 

135 

136 

137def assert_masked_images_equal( 

138 tc: unittest.TestCase, 

139 a: MaskedImage, 

140 b: MaskedImage, 

141 *, 

142 rtol: float = 0.0, 

143 atol: float = 0.0, 

144 expect_view: bool | None = None, 

145) -> None: 

146 """Assert that two masked images are equal or nearly equal.""" 

147 tc.assertEqual(a.metadata, b.metadata) 

148 assert_projections_equal(tc, a.projection, b.projection) 

149 assert_images_equal(tc, a.image, b.image, rtol=rtol, atol=atol, expect_view=expect_view) 

150 assert_masks_equal(tc, a.mask, b.mask) 

151 assert_images_equal(tc, a.variance, b.variance, rtol=rtol, atol=atol, expect_view=expect_view) 

152 

153 

154def assert_psfs_equal( 

155 tc: unittest.TestCase, 

156 psf1: PointSpreadFunction, 

157 psf2: PointSpreadFunction, 

158 points: YX[np.ndarray] | XY[np.ndarray] | None = None, 

159) -> int: 

160 """Compare two PSF objets. 

161 

162 Parameters 

163 ---------- 

164 tc 

165 Test case object with assert methods to use. 

166 psf1 

167 Point-spread function to test. 

168 psf2 

169 The other point-spread function to test. 

170 points 

171 Points to evaluate the PSFs at. If not provided, the intersection of 

172 the PSF bounding boxes are used to generate points on a grid. 

173 

174 Returns 

175 ------- 

176 `int` 

177 The number of points actually tested. 

178 """ 

179 if points is None: 

180 points = psf1.bounds.bbox.intersection(psf1.bounds.bbox).meshgrid(3).map(np.ravel) 

181 

182 tc.assertEqual(psf1.kernel_bbox, psf2.kernel_bbox) 

183 

184 n_points_tested: int = 0 

185 for x, y in zip(points.x, points.y): 

186 if not psf1.bounds.contains(x=x, y=y): 

187 with tc.assertRaises(BoundsError): 

188 psf2.compute_kernel_image(x=x, y=y) 

189 continue 

190 tc.assertEqual(psf1.compute_kernel_image(x=x, y=y), psf2.compute_kernel_image(x=x, y=y)) 

191 tc.assertEqual(psf1.compute_stellar_bbox(x=x, y=y), psf2.compute_stellar_bbox(x=x, y=y)) 

192 tc.assertEqual(psf1.compute_stellar_image(x=x, y=y), psf2.compute_stellar_image(x=x, y=y)) 

193 n_points_tested += 1 

194 return n_points_tested 

195 

196 

197def compare_image_to_legacy( 

198 tc: unittest.TestCase, image: Image, legacy_image: Any, expect_view: bool | None = None 

199) -> None: 

200 """Compare an `.Image` object to a legacy `lsst.afw.image.Image` object.""" 

201 tc.assertEqual(image.bbox, Box.from_legacy(legacy_image.getBBox())) 

202 if expect_view is not None: 

203 tc.assertEqual(np.may_share_memory(image.array, legacy_image.array), expect_view) 

204 if not expect_view: 

205 np.testing.assert_array_equal(image.array, legacy_image.array) 

206 

207 

208def compare_mask_to_legacy( 

209 tc: unittest.TestCase, mask: Mask, legacy_mask: Any, plane_map: Mapping[str, MaskPlane] | None = None 

210) -> None: 

211 """Compare a `.Mask` object to a legacy `lsst.afw.image.Mask` object.""" 

212 tc.assertEqual(mask.bbox, Box.from_legacy(legacy_mask.getBBox())) 

213 if plane_map is None: 

214 plane_map = {plane.name: plane for plane in mask.schema if plane is not None} 

215 for old_name, new_plane in plane_map.items(): 

216 np.testing.assert_array_equal( 

217 (legacy_mask.array & legacy_mask.getPlaneBitMask(old_name)).astype(bool), 

218 mask.get(new_plane.name), 

219 ) 

220 

221 

222def compare_masked_image_to_legacy( 

223 tc: unittest.TestCase, 

224 masked_image: MaskedImage, 

225 legacy_masked_image: Any, 

226 *, 

227 plane_map: Mapping[str, MaskPlane] | None = None, 

228 expect_view: bool | None = None, 

229 alternates: Mapping[str, Any] | None = None, 

230) -> None: 

231 """Compare a `.MaskedImage` object to a legacy `lsst.afw.image.MaskedImage` 

232 object. 

233 

234 Parameters 

235 ---------- 

236 tc 

237 Test case to use for asserts. 

238 masked_image 

239 New image to test. 

240 legacy_masked_image 

241 Legacy image to test against. 

242 plane_map 

243 Mapping between new and legacy mask planes. 

244 expect_view 

245 Whether to test that the image and variance arrays do or do not share 

246 memory. 

247 alternates 

248 A mapping of other versions of one or more (new) components to also 

249 check against the legacy versions of those components. 

250 """ 

251 compare_image_to_legacy(tc, masked_image.image, legacy_masked_image.getImage(), expect_view=expect_view) 

252 compare_mask_to_legacy(tc, masked_image.mask, legacy_masked_image.getMask(), plane_map=plane_map) 

253 compare_image_to_legacy( 

254 tc, masked_image.variance, legacy_masked_image.getVariance(), expect_view=expect_view 

255 ) 

256 if alternates: 

257 if image := alternates.get("image"): 

258 compare_image_to_legacy(tc, image, legacy_masked_image.getImage(), expect_view=expect_view) 

259 if mask := alternates.get("mask"): 

260 compare_mask_to_legacy(tc, mask, legacy_masked_image.getMask(), plane_map=plane_map) 

261 if variance := alternates.get("variance"): 

262 compare_image_to_legacy(tc, variance, legacy_masked_image.getVariance(), expect_view=expect_view) 

263 

264 

265def compare_visit_image_to_legacy( 

266 tc: unittest.TestCase, 

267 visit_image: VisitImage, 

268 legacy_exposure: Any, 

269 *, 

270 plane_map: Mapping[str, MaskPlane] | None = None, 

271 expect_view: bool | None = None, 

272 instrument: str, 

273 visit: int, 

274 detector: int, 

275 alternates: Mapping[str, Any] | None = None, 

276) -> None: 

277 """Compare a `.VisitImage` object to a legacy `lsst.afw.image.Exposure` 

278 object. 

279 

280 Parameters 

281 ---------- 

282 tc 

283 Test case to use for asserts. 

284 visit_image 

285 New image to test. 

286 legacy_exposure 

287 Legacy image to test against. 

288 plane_map 

289 Mapping between new and legacy mask planes. 

290 expect_view 

291 Whether to test that the image and variance arrays do or do not share 

292 memory. 

293 instrument 

294 Expected instrument name. 

295 visit 

296 Expected visit ID. 

297 detector 

298 Expected detector ID. 

299 alternates 

300 A mapping of other versions of one or more (new) components to also 

301 check against the legacy versions of those components. 

302 """ 

303 compare_masked_image_to_legacy( 

304 tc, 

305 visit_image, 

306 legacy_exposure.getMaskedImage(), 

307 plane_map=plane_map, 

308 expect_view=expect_view, 

309 alternates=alternates, 

310 ) 

311 detector_bbox = Box.from_legacy(legacy_exposure.getDetector().getBBox()) 

312 compare_projection_to_legacy_wcs( 

313 tc, 

314 visit_image.projection, 

315 legacy_exposure.getWcs(), 

316 DetectorFrame(instrument=instrument, visit=visit, detector=detector, bbox=detector_bbox), 

317 visit_image.bbox, 

318 ) 

319 tc.assertIs(visit_image.projection, visit_image.mask.projection) 

320 tc.assertIs(visit_image.projection, visit_image.variance.projection) 

321 compare_psf_to_legacy(tc, visit_image.psf, legacy_exposure.getPsf()) 

322 compare_observation_summary_stats_to_legacy( 

323 tc, visit_image.summary_stats, legacy_exposure.info.getSummaryStats() 

324 ) 

325 compare_detector_to_legacy(tc, visit_image.detector, legacy_exposure.getDetector(), is_raw_assembled=True) 

326 # Make a tiny box for Field comparisons that need to make arrays; that can 

327 # get expensive otherwisre. 

328 tiny_bbox = detector_bbox.local[2:4, 3:6] 

329 compare_aperture_corrections_to_legacy( 

330 tc, visit_image.aperture_corrections, legacy_exposure.info.getApCorrMap(), tiny_bbox 

331 ) 

332 if alternates: 

333 if projection := alternates.get("projection"): 

334 compare_projection_to_legacy_wcs( 

335 tc, 

336 projection, 

337 legacy_exposure.getWcs(), 

338 DetectorFrame(instrument=instrument, visit=visit, detector=detector, bbox=detector_bbox), 

339 visit_image.bbox, 

340 ) 

341 if psf := alternates.get("psf"): 

342 compare_psf_to_legacy(tc, psf, legacy_exposure.getPsf()) 

343 if summary_stats := alternates.get("summary_stats"): 

344 compare_observation_summary_stats_to_legacy( 

345 tc, summary_stats, legacy_exposure.info.getSummaryStats() 

346 ) 

347 if detector_obj := alternates.get("detector"): 

348 compare_detector_to_legacy(tc, detector_obj, legacy_exposure.getDetector(), is_raw_assembled=True) 

349 if obs_info := alternates.get("obs_info"): 

350 visitInfo = legacy_exposure.visitInfo 

351 tc.assertEqual(obs_info.instrument, visitInfo.getInstrumentLabel()) 

352 if aperture_corrections := alternates.get("aperture_corrections"): 

353 compare_aperture_corrections_to_legacy( 

354 tc, aperture_corrections, legacy_exposure.info.getApCorrMap(), tiny_bbox 

355 ) 

356 

357 

358def compare_cell_coadd_to_legacy( 

359 tc: unittest.TestCase, 

360 cell_coadd: CellCoadd, 

361 legacy_cell_coadd: MultipleCellCoadd, 

362 *, 

363 tract_bbox: Box, 

364 plane_map: Mapping[str, MaskPlane] | None = None, 

365 alternates: Mapping[str, Any] | None = None, 

366 psf_points: XY[np.ndarray] | YX[np.ndarray] | None = None, 

367) -> None: 

368 """Compare a `.cells.CellCoadd` object to a legacy 

369 `lsst.cell_coadds.MultipleCellCoadd` object. 

370 

371 Parameters 

372 ---------- 

373 tc 

374 Test case to use for asserts. 

375 cell_coadd 

376 New coadd to test. 

377 legacy_cell_coadd 

378 Legacy coadd to test against. 

379 tract_bbox 

380 Bounding box of the full tract. 

381 psf_points 

382 Points to use to compare the PSFs. 

383 plane_map 

384 Mapping between new and legacy mask planes. 

385 alternates 

386 A mapping of other versions of one or more (new) components to also 

387 check against the legacy versions of those components. 

388 """ 

389 legacy_stitched = legacy_cell_coadd.stitch(cell_coadd.bbox.to_legacy()) 

390 compare_image_to_legacy(tc, cell_coadd.image, legacy_stitched.image, expect_view=False) 

391 compare_mask_to_legacy(tc, cell_coadd.mask, legacy_stitched.mask, plane_map=plane_map) 

392 compare_image_to_legacy(tc, cell_coadd.variance, legacy_stitched.variance, expect_view=False) 

393 if legacy_stitched.mask_fractions is not None: 

394 compare_image_to_legacy( 

395 tc, cell_coadd.mask_fractions["rejected"], legacy_stitched.mask_fractions, expect_view=False 

396 ) 

397 for n in range(legacy_stitched.n_noise_realizations): 

398 compare_image_to_legacy( 

399 tc, cell_coadd.noise_realizations[n], legacy_stitched.noise_realizations[n], expect_view=False 

400 ) 

401 tc.assertEqual(cell_coadd.skymap, legacy_stitched.identifiers.skymap) 

402 tc.assertEqual(cell_coadd.tract, legacy_stitched.identifiers.tract) 

403 tc.assertEqual(cell_coadd.patch.index.x, legacy_stitched.identifiers.patch.x) 

404 tc.assertEqual(cell_coadd.patch.index.y, legacy_stitched.identifiers.patch.y) 

405 tc.assertEqual(cell_coadd.band, legacy_stitched.identifiers.band) 

406 tc.assertTrue(tract_bbox.contains(cell_coadd.patch.outer_bbox)) 

407 tc.assertTrue(cell_coadd.patch.outer_bbox.contains(cell_coadd.patch.inner_bbox)) 

408 tc.assertTrue(cell_coadd.patch.outer_bbox.contains(cell_coadd.bbox)) 

409 tc.assertEqual(cell_coadd.unit, u.Unit(legacy_cell_coadd.common.units.value)) 

410 tc.assertTrue(cell_coadd.bounds.bbox.contains(cell_coadd.bbox)) 

411 tc.assertTrue(cell_coadd.grid.bbox.contains(cell_coadd.bbox)) 

412 compare_projection_to_legacy_wcs( 

413 tc, 

414 cell_coadd.projection, 

415 legacy_cell_coadd.wcs, 

416 TractFrame( 

417 skymap=legacy_cell_coadd.identifiers.skymap, 

418 tract=legacy_cell_coadd.identifiers.tract, 

419 bbox=tract_bbox, 

420 ), 

421 cell_coadd.bbox, 

422 is_fits=True, 

423 ) 

424 tc.assertIs(cell_coadd.projection, cell_coadd.mask.projection) 

425 tc.assertIs(cell_coadd.projection, cell_coadd.variance.projection) 

426 compare_psf_to_legacy( 

427 tc, cell_coadd.psf, legacy_stitched.psf, expect_legacy_raise_on_out_of_bounds=True, points=psf_points 

428 ) 

429 compare_cell_coadd_provenance_to_legacy(tc, cell_coadd.provenance, legacy_cell_coadd) 

430 if alternates: 

431 if projection := alternates.get("projection"): 

432 compare_projection_to_legacy_wcs( 

433 tc, 

434 projection, 

435 legacy_stitched.wcs, 

436 TractFrame( 

437 skymap=legacy_cell_coadd.identifiers.skymap, 

438 tract=legacy_cell_coadd.identifiers.tract, 

439 bbox=tract_bbox, 

440 ), 

441 cell_coadd.bbox, 

442 is_fits=True, 

443 ) 

444 if psf := alternates.get("psf"): 

445 compare_psf_to_legacy(tc, psf, legacy_stitched.psf, points=psf_points) 

446 if provenance := alternates.get("provenance"): 

447 compare_cell_coadd_provenance_to_legacy(tc, provenance, legacy_cell_coadd) 

448 

449 

450def compare_cell_coadd_provenance_to_legacy( 

451 tc: unittest.TestCase, provenance: CoaddProvenance, legacy_cell_coadd: MultipleCellCoadd 

452) -> None: 

453 """Compare a `.cells.CoaddProvenance` object to a legacy 

454 `lsst.cell_coadds.MultipleCellCoadd` object. 

455 

456 Parameters 

457 ---------- 

458 tc 

459 Test case to use for asserts. 

460 provenance 

461 New provenance object to test. 

462 legacy_cell_coadd 

463 Legacy coadd to test against. 

464 """ 

465 from lsst.cell_coadds import ObservationIdentifiers 

466 

467 for legacy_cell in legacy_cell_coadd.cells.values(): 

468 cell_index = CellIJ.from_legacy(legacy_cell.identifiers.cell) 

469 prov = provenance[cell_index] 

470 legacy_table = astropy.table.Table( 

471 rows=[ 

472 [ 

473 ids.instrument, 

474 ids.visit, 

475 ids.detector, 

476 ids.day_obs, 

477 ids.physical_filter, 

478 legacy_input.overlaps_center, 

479 legacy_input.overlap_fraction, 

480 legacy_input.weight, 

481 legacy_input.psf_shape.getIxx(), 

482 legacy_input.psf_shape.getIyy(), 

483 legacy_input.psf_shape.getIxy(), 

484 legacy_input.psf_shape_flag, 

485 ] 

486 for ids, legacy_input in legacy_cell.inputs.items() 

487 ], 

488 dtype=[ 

489 np.object_, 

490 np.uint64, 

491 np.uint16, 

492 np.uint32, 

493 np.object_, 

494 np.bool_, 

495 np.float64, 

496 np.float64, 

497 np.float64, 

498 np.float64, 

499 np.float64, 

500 np.bool_, 

501 ], 

502 names=[ 

503 "instrument", 

504 "visit", 

505 "detector", 

506 "day_obs", 

507 "physical_filter", 

508 "overlaps_center", 

509 "overlap_fraction", 

510 "weight", 

511 "psf_shape_xx", 

512 "psf_shape_yy", 

513 "psf_shape_xy", 

514 "psf_shape_flag", 

515 ], 

516 ) 

517 # For a single cell all 'inputs' are also 'contributions'. 

518 tc.assertEqual(len(legacy_cell.inputs), len(prov.inputs)) 

519 tc.assertEqual(len(legacy_cell.inputs), len(prov.contributions)) 

520 prov.inputs.sort(["instrument", "visit", "detector"]) 

521 prov.contributions.sort(["instrument", "visit", "detector"]) 

522 legacy_table.sort(["instrument", "visit", "detector"]) 

523 np.testing.assert_array_equal(prov.inputs["instrument"], prov.contributions["instrument"]) 

524 np.testing.assert_array_equal(prov.inputs["visit"], prov.contributions["visit"]) 

525 np.testing.assert_array_equal(prov.inputs["detector"], prov.contributions["detector"]) 

526 np.testing.assert_array_equal(prov.inputs["instrument"], legacy_table["instrument"]) 

527 np.testing.assert_array_equal(prov.inputs["visit"], legacy_table["visit"]) 

528 np.testing.assert_array_equal(prov.inputs["detector"], legacy_table["detector"]) 

529 np.testing.assert_array_equal(prov.inputs["physical_filter"], legacy_table["physical_filter"]) 

530 np.testing.assert_array_equal(prov.inputs["day_obs"], legacy_table["day_obs"]) 

531 np.testing.assert_array_equal(prov.contributions["overlaps_center"], legacy_table["overlaps_center"]) 

532 np.testing.assert_array_equal( 

533 prov.contributions["overlap_fraction"], legacy_table["overlap_fraction"] 

534 ) 

535 np.testing.assert_array_equal(prov.contributions["weight"], legacy_table["weight"]) 

536 np.testing.assert_array_equal(prov.contributions["psf_shape_xx"], legacy_table["psf_shape_xx"]) 

537 np.testing.assert_array_equal(prov.contributions["psf_shape_yy"], legacy_table["psf_shape_yy"]) 

538 np.testing.assert_array_equal(prov.contributions["psf_shape_xy"], legacy_table["psf_shape_xy"]) 

539 np.testing.assert_array_equal(prov.contributions["psf_shape_flag"], legacy_table["psf_shape_flag"]) 

540 for row in prov.inputs: 

541 polygon_key = ObservationIdentifiers(**{k: row[k] for k in row.keys() if k != "polygon"}) 

542 legacy_polygon = legacy_cell_coadd.common.visit_polygons[polygon_key] 

543 tc.assertEqual(legacy_polygon, row["polygon"].to_legacy()) 

544 

545 

546def compare_psf_to_legacy( 

547 tc: unittest.TestCase, 

548 psf: PointSpreadFunction, 

549 legacy_psf: Any, 

550 points: YX[np.ndarray] | XY[np.ndarray] | None = None, 

551 expect_legacy_raise_on_out_of_bounds: bool = False, 

552) -> int: 

553 """Compare a PSF model object to its legacy interface. 

554 

555 Parameters 

556 ---------- 

557 tc 

558 Test case object with assert methods to use. 

559 psf 

560 Point-spread function to test. 

561 legacy_psf 

562 Legacy `lsst.afw.detection.Psf` instance to compare with. 

563 points 

564 Points to evaluate the PSFs at. If not provided, the intersection of 

565 the PSF bounding boxes are used to generate points on a grid. 

566 expect_legacy_raise_on_out_of_bounds 

567 If `True`, expect ``legacy_psf`` to raise 

568 `lsst.afw.detection.InvalidPsfError` when evaluated at a position 

569 considered out-of-bounds by ``psf``. 

570 

571 Returns 

572 ------- 

573 `int` 

574 The number of points actually tested. 

575 """ 

576 from lsst.afw.detection import InvalidPsfError 

577 

578 if points is None: 

579 points = psf.bounds.bbox.meshgrid(n=3).map(np.ravel) 

580 legacy_points = arrays_to_legacy_points(points.x, points.y) 

581 n_points_tested: int = 0 

582 for p in legacy_points: 

583 if not psf.bounds.contains(x=p.x, y=p.y): 

584 if expect_legacy_raise_on_out_of_bounds: 

585 with tc.assertRaises(InvalidPsfError): 

586 legacy_psf.computeKernelImage(p) 

587 continue 

588 tc.assertEqual(psf.kernel_bbox, Box.from_legacy(legacy_psf.computeKernelBBox(p))) 

589 tc.assertEqual( 

590 psf.compute_kernel_image(x=p.x, y=p.y), Image.from_legacy(legacy_psf.computeKernelImage(p)) 

591 ) 

592 tc.assertEqual( 

593 psf.compute_stellar_bbox(x=p.x, y=p.y), Box.from_legacy(legacy_psf.computeImageBBox(p)) 

594 ) 

595 tc.assertEqual(psf.compute_stellar_image(x=p.x, y=p.y), Image.from_legacy(legacy_psf.computeImage(p))) 

596 n_points_tested += 1 

597 return n_points_tested 

598 

599 

600def compare_field_to_legacy( 

601 tc: unittest.TestCase, 

602 field: BaseField, 

603 legacy_field: Any, 

604 subimage_bbox: Box, 

605) -> None: 

606 """Test a Field object by comparing it to an equivalent 

607 `lsst.afw.math.BoundedField`. 

608 

609 Parameters 

610 ---------- 

611 tc 

612 Test case object with assert methods to use. 

613 field 

614 Field to test. 

615 legacy_field : ``lsst.afw.math.BoundedField`` 

616 Equivalent legacy bounded field. 

617 subimage_bbox 

618 Bounding box for full-image tests. 

619 """ 

620 tc.assertEqual(field.bounds.bbox, Box.from_legacy(legacy_field.getBBox())) 

621 # Pixel coordinates to test the numpy array interface with. 

622 pixel_xy = field.bounds.bbox.meshgrid(n=5).map(np.ravel) 

623 assert_close(tc, field(x=pixel_xy.x, y=pixel_xy.y), legacy_field.evaluate(pixel_xy.x, pixel_xy.y)) 

624 legacy_image_1 = Image(0, bbox=subimage_bbox, dtype=np.float64).to_legacy() 

625 legacy_field.addToImage(legacy_image_1, overlapOnly=True) 

626 assert_images_equal(tc, field.render(subimage_bbox), Image.from_legacy(legacy_image_1), rtol=1e-13) 

627 

628 

629def compare_aperture_corrections_to_legacy( 

630 tc: unittest.TestCase, 

631 aperture_corrections: ApertureCorrectionMap, 

632 legacy_ap_corr_map: Any, 

633 subimage_bbox: Box, 

634) -> None: 

635 """Test an aperture correction `dict` by comparing it to an equivalent 

636 `lsst.afw.image.ApCorrMap`. 

637 

638 Parameters 

639 ---------- 

640 tc 

641 Test case object with assert methods to use. 

642 aperture_corrections 

643 Dictionary to test. 

644 legacy_ap_corr_map : ``lsst.afw.image.ApCorrMap`` 

645 Equivalent legacy aperture correction map. 

646 subimage_bbox 

647 Bounding box for full-image tests. 

648 """ 

649 tc.assertEqual(aperture_corrections.keys(), set(legacy_ap_corr_map.keys())) 

650 for name, field in aperture_corrections.items(): 

651 compare_field_to_legacy(tc, field, legacy_ap_corr_map[name], subimage_bbox) 

652 

653 

654def compare_observation_summary_stats_to_legacy( 

655 tc: unittest.TestCase, 

656 summary_stats: ObservationSummaryStats, 

657 legacy_summary_stats: Any, 

658) -> None: 

659 """Test an ObservationSummaryStats object by comparing it to an equivalent 

660 `lsst.afw.image.ExposureSummaryStats`. 

661 

662 Parameters 

663 ---------- 

664 tc 

665 Test case object with assert methods to use. 

666 summary_stats 

667 Struct to test. 

668 legacy : ``lsst.afw.image.ExposureSummaryStats`` 

669 Equivalent legacy struct. 

670 """ 

671 for field in dataclasses.fields(legacy_summary_stats): 

672 a = getattr(legacy_summary_stats, field.name) 

673 b = getattr(summary_stats, field.name) 

674 if isinstance(b, tuple): 

675 for ai, bi in zip(a, b): 

676 tc.assertTrue(ai == bi or (math.isnan(ai) and math.isnan(bi)), f"{field.name}: {a} != {b}") 

677 else: 

678 tc.assertTrue(a == b or (math.isnan(a) and math.isnan(b)), f"{field.name}: {a} != {b}") 

679 

680 

681def compare_projection_to_legacy_wcs[F: Frame]( 

682 tc: unittest.TestCase, 

683 projection: Projection[F], 

684 legacy_wcs: Any, 

685 pixel_frame: F, 

686 subimage_bbox: Box, 

687 is_fits: bool = False, 

688) -> None: 

689 """Test a Projection object by comparing it to an equivalent 

690 `lsst.afw.geom.SkyWcs`. 

691 

692 Parameters 

693 ---------- 

694 tc 

695 Test case object with assert methods to use. 

696 projection 

697 Projection to test. 

698 legacy_wcs : ``lsst.afw.geom.SkyWcs`` 

699 Equivalent legacy WCS. 

700 pixel_frame 

701 Expected pixel frame for the projection. 

702 subimage_bbox 

703 Bounding box of points to generate for tests. 

704 is_fits 

705 Whether this projection is expected to be exactly representable as a 

706 FIT WCS. If `False` it is assumed to have a FITS approximation 

707 attached instead. 

708 """ 

709 # Pixel coordinates to test on over the subimage region of interest: 

710 pixel_xy = subimage_bbox.meshgrid(step=50).map(np.ravel) 

711 # Array indices of those pixel values (subtract off bbox starts): 

712 subimage_array_xy = XY(x=pixel_xy.x - subimage_bbox.x.start, y=pixel_xy.y - subimage_bbox.y.start) 

713 sky_coords = legacy_coords_to_astropy( 

714 legacy_wcs.pixelToSky(arrays_to_legacy_points(pixel_xy.x, pixel_xy.y)) 

715 ) 

716 # Test transforming with the Projection itself, which also tests its 

717 # nested Transform and an Astropy High-Level WCS view with no origin 

718 # change. 

719 check_projection(tc, projection, pixel_xy, sky_coords, pixel_frame) 

720 # Also test the Astropy High-Level WCS view with an origin change to 

721 # array indices. 

722 check_astropy_wcs_interface( 

723 tc, projection.as_astropy(subimage_bbox), subimage_array_xy, sky_coords, pixel_atol=1e-5 

724 ) 

725 if is_fits: 

726 fits_wcs = projection.as_fits_wcs(subimage_bbox, allow_approximation=True) 

727 assert fits_wcs is not None 

728 check_astropy_wcs_interface(tc, fits_wcs, subimage_array_xy, sky_coords, pixel_atol=1e-5) 

729 # Use that FITS approximation to check that we can make a 

730 # Projection from a FITS WCS, too. 

731 fits_projection = Projection.from_fits_wcs(fits_wcs, pixel_frame) 

732 check_projection( 

733 tc, 

734 fits_projection, 

735 subimage_array_xy, 

736 sky_coords, 

737 pixel_frame, 

738 pixel_atol=1e-5, 

739 ) 

740 # We want Projections we create from a FITS WCS to be backed by an 

741 # AST FrameSet so we can convert them into legacy 

742 # `lsst.afw.geom.SkyWcs` objects if desired. 

743 tc.assertIn("Begin FrameSet", fits_projection.show()) 

744 else: 

745 tc.assertIsNone(projection.as_fits_wcs(subimage_bbox, allow_approximation=False)) 

746 # The legacy SkyWcs should instead have a FITS approximation 

747 # attached; run the same tests on that. 

748 assert projection.fits_approximation is not None 

749 compare_projection_to_legacy_wcs( 

750 tc, 

751 projection.fits_approximation, 

752 legacy_wcs.getFitsApproximation(), 

753 pixel_frame, 

754 subimage_bbox, 

755 is_fits=True, 

756 ) 

757 

758 

759def check_transform[I: Frame, O: Frame]( 

760 tc: unittest.TestCase, 

761 transform: Transform[I, O], 

762 input_xy: XY[np.ndarray], 

763 output_xy: XY[np.ndarray],