Coverage for tests/test_extended_psf.py: 99%
346 statements
« prev ^ index » next coverage.py v7.15.2, created at 2026-07-25 09:10 +0000
1# This file is part of pipe_tasks.
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# This program is free software: you can redistribute it and/or modify
10# it under the terms of the GNU General Public License as published by
11# the Free Software Foundation, either version 3 of the License, or
12# (at your option) any later version.
13#
14# This program is distributed in the hope that it will be useful,
15# but WITHOUT ANY WARRANTY; without even the implied warranty of
16# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17# GNU General Public License for more details.
18#
19# You should have received a copy of the GNU General Public License
20# along with this program. If not, see <https://www.gnu.org/licenses/>.
22import logging
23import unittest
24from typing import cast
26import astropy.units as u
27import lsst.images.tests
28import lsst.utils.tests
29import numpy as np
30from astropy.table import Table
31from lsst.afw.cameraGeom import FOCAL_PLANE, PIXELS
32from lsst.afw.cameraGeom.testUtils import CameraWrapper
33from lsst.afw.geom import makeCdMatrix, makeSkyWcs
34from lsst.afw.image import (
35 ExposureF,
36 ImageD,
37 ImageF,
38 MaskedImageF,
39 VisitInfo,
40 makePhotoCalibFromCalibZeroPoint,
41)
42from lsst.afw.math import (
43 REDUCE_INTERP_ORDER,
44 ApproximateControl,
45 BackgroundControl,
46 BackgroundList,
47 FixedKernel,
48 Interpolate,
49 makeBackground,
50)
51from lsst.geom import Box2I, Extent2I, Point2D, Point2I, SpherePoint, arcseconds, degrees
52from lsst.images import Box, Image, Mask, MaskedImage, MaskSchema
53from lsst.meas.algorithms import KernelPsf
54from lsst.pipe.tasks.extended_psf import (
55 ExtendedPsfCandidate,
56 ExtendedPsfCandidateInfo,
57 ExtendedPsfCandidates,
58 ExtendedPsfCutoutConfig,
59 ExtendedPsfCutoutTask,
60 ExtendedPsfFit,
61 ExtendedPsfImage,
62 ExtendedPsfImageInfo,
63 ExtendedPsfMoffatFit,
64 ExtendedPsfStackConfig,
65 ExtendedPsfStackTask,
66 ExtendedPsfSubtractConfig,
67 ExtendedPsfSubtractTask,
68)
71class ExtendedPsfCandidatesTestCase(lsst.utils.tests.TestCase):
72 """Test ExtendedPsfCandidate and ExtendedPsfCandidates data structures."""
74 @classmethod
75 def setUpClass(cls):
76 rng = np.random.Generator(np.random.MT19937(seed=5))
77 cutout_size = (25, 25)
79 # Generate simulated stars
80 extended_psf_candidates = []
81 cls.star_infos = []
82 mask_schema = MaskSchema([])
84 for i in range(3):
85 candidate_masked_image = MaskedImage(
86 image=Image(rng.random(cutout_size).astype(np.float32)),
87 mask=Mask(0, schema=mask_schema, shape=cutout_size),
88 variance=Image(rng.random(cutout_size).astype(np.float32)),
89 )
91 star_info = ExtendedPsfCandidateInfo(
92 visit=100 + i,
93 detector=200 + i,
94 ref_id=1000 + i,
95 ref_mag=10.0 + i,
96 position_x=float(rng.random()),
97 position_y=float(rng.random()),
98 focal_plane_radius=float(rng.random()) * u.mm,
99 focal_plane_angle=float(rng.random()) * u.rad,
100 )
101 cls.star_infos.append(star_info)
103 # Build a normalized 2-D Gaussian kernel image directly.
104 yy, xx = np.mgrid[: cutout_size[0], : cutout_size[1]]
105 cy = (cutout_size[0] - 1) / 2.0
106 cx = (cutout_size[1] - 1) / 2.0
107 sigma = 1.5
108 kernel_array = np.exp(-((yy - cy) ** 2 + (xx - cx) ** 2) / (2.0 * sigma**2))
109 kernel_array /= np.sum(kernel_array)
110 psf_kernel_image = Image(kernel_array.astype(np.float64))
112 extended_psf_candidates.append(
113 ExtendedPsfCandidate(
114 image=candidate_masked_image.image,
115 mask=candidate_masked_image.mask,
116 variance=candidate_masked_image.variance,
117 psf_kernel_image=psf_kernel_image,
118 star_info=star_info,
119 )
120 )
122 cls.global_metadata = {"CANDIDATES_TEST_KEY": "CANDIDATES_TEST_VALUE"}
123 cls.extended_psf_candidates = ExtendedPsfCandidates(extended_psf_candidates, cls.global_metadata)
125 @classmethod
126 def tearDownClass(cls):
127 del cls.extended_psf_candidates
128 del cls.star_infos
129 del cls.global_metadata
131 def test_fits_roundtrip(self):
132 """Test that ExtendedPsfCandidates can be serialized/deserialized."""
134 with lsst.images.tests.RoundtripFits(
135 self.extended_psf_candidates, storage_class="ExtendedPsfCandidates"
136 ) as roundtrip:
137 pass
138 extended_psf_candidates = roundtrip.result
140 global_metadata = extended_psf_candidates.metadata
141 self.assertEqual(self.global_metadata["CANDIDATES_TEST_KEY"], global_metadata["CANDIDATES_TEST_KEY"])
142 self.assertEqual(len(self.extended_psf_candidates), len(extended_psf_candidates))
143 for input_info, input_candidate, output_candidate in zip(
144 self.star_infos, self.extended_psf_candidates, extended_psf_candidates
145 ):
146 self.assertEqual(input_candidate.metadata, {})
147 self.assertEqual(output_candidate.metadata, {})
149 output_info = output_candidate.star_info
150 self.assertEqual(input_info.visit, output_info.visit)
151 self.assertEqual(input_info.detector, output_info.detector)
152 self.assertEqual(input_info.ref_id, output_info.ref_id)
153 self.assertAlmostEqual(input_info.ref_mag, output_info.ref_mag, places=10)
154 self.assertAlmostEqual(input_info.position_x, output_info.position_x, places=10)
155 self.assertAlmostEqual(input_info.position_y, output_info.position_y, places=10)
157 self.assertIsNotNone(input_info.focal_plane_radius)
158 self.assertIsNotNone(output_info.focal_plane_radius)
159 input_radius = cast(u.Quantity, input_info.focal_plane_radius)
160 output_radius = cast(u.Quantity, output_info.focal_plane_radius)
161 self.assertAlmostEqual(input_radius.to_value(u.mm), output_radius.to_value(u.mm), places=10)
163 self.assertIsNotNone(input_info.focal_plane_angle)
164 self.assertIsNotNone(output_info.focal_plane_angle)
165 input_angle = cast(u.Quantity, input_info.focal_plane_angle)
166 output_angle = cast(u.Quantity, output_info.focal_plane_angle)
167 self.assertAlmostEqual(input_angle.to_value(u.rad), output_angle.to_value(u.rad), places=10)
169 np.testing.assert_allclose(
170 input_candidate.image.array,
171 output_candidate.image.array,
172 rtol=0.0,
173 atol=1e-10,
174 )
175 np.testing.assert_array_equal(
176 input_candidate.mask.array,
177 output_candidate.mask.array,
178 )
179 np.testing.assert_allclose(
180 input_candidate.variance.array, output_candidate.variance.array, rtol=0.0, atol=1e-10
181 )
182 np.testing.assert_allclose(
183 input_candidate.psf_kernel_image.array,
184 output_candidate.psf_kernel_image.array,
185 rtol=0.0,
186 atol=1e-10,
187 )
189 def test_ref_id_map(self):
190 """Test that ref_id_map maps reference catalog IDs to candidates."""
191 ref_id_map = self.extended_psf_candidates.ref_id_map
192 self.assertEqual(len(ref_id_map), len(self.star_infos))
193 for i, star_info in enumerate(self.star_infos):
194 self.assertIn(star_info.ref_id, ref_id_map)
195 self.assertIs(ref_id_map[star_info.ref_id], self.extended_psf_candidates[i])
197 def test_slice_preserves_subcollection_and_metadata(self):
198 """Test that slicing candidates returns a sub-collection."""
199 sub = self.extended_psf_candidates[1:3] # length=2
200 self.assertIsInstance(sub, ExtendedPsfCandidates)
201 self.assertEqual(len(sub), 2)
202 self.assertEqual(sub.metadata, self.global_metadata)
203 self.assertIs(sub[0], self.extended_psf_candidates[1])
204 self.assertIs(sub[1], self.extended_psf_candidates[2])
206 def test_spatial_box_slicing(self):
207 """Test that candidates support spatial slicing with a Box."""
208 candidate = self.extended_psf_candidates[0]
209 subbox = Box.from_shape((10, 10), start=(5, 5))
210 sub = candidate[subbox]
211 self.assertIsInstance(sub, ExtendedPsfCandidate)
212 self.assertEqual(sub.bbox, subbox)
213 self.assertIs(sub.star_info, candidate.star_info)
214 np.testing.assert_array_equal(sub.image.array, candidate.image[subbox].array)
215 np.testing.assert_array_equal(sub.mask.array, candidate.mask[subbox].array)
217 def test_copy(self):
218 """Test that copy() produces an independent deep copy."""
219 candidate = self.extended_psf_candidates[0]
220 copied = candidate.copy()
221 self.assertIsInstance(copied, ExtendedPsfCandidate)
222 np.testing.assert_array_equal(copied.image.array, candidate.image.array)
223 np.testing.assert_array_equal(copied.variance.array, candidate.variance.array)
224 # Modifying the copy must not affect the original.
225 copied.image.array[:] = 0.0
226 self.assertFalse(np.all(candidate.image.array == 0.0))
229class ExtendedPsfImageTestCase(lsst.utils.tests.TestCase):
230 """Test ExtendedPsfImage data model, properties, and operations."""
232 def setUp(self):
233 image = Image(
234 np.arange(30, dtype=np.float32).reshape(5, 6),
235 yx0=(10, -3),
236 unit=u.nJy,
237 )
238 variance = Image(
239 np.full((5, 6), 2.5, dtype=np.float32),
240 bbox=image.bbox,
241 unit=u.nJy**2,
242 )
243 info = ExtendedPsfImageInfo(
244 n_stars=17,
245 )
246 fit = ExtendedPsfMoffatFit(
247 chi2=12.3,
248 reduced_chi2=1.23,
249 dof=10,
250 amplitude=0.8,
251 x_0=-0.2,
252 y_0=0.4,
253 gamma=2.3,
254 alpha=4.5,
255 )
256 self.extended_psf_image = ExtendedPsfImage(
257 image=image,
258 variance=variance,
259 info=info,
260 fit=fit,
261 metadata={"EPSF_TEST_KEY": "EPSF_TEST VALUE"},
262 )
264 def tearDown(self):
265 del self.extended_psf_image
267 def test_fits_roundtrip(self):
268 """Test that ExtendedPsfImage can be serialized/deserialized."""
269 with lsst.images.tests.RoundtripFits(
270 self.extended_psf_image,
271 storage_class="ExtendedPsfImage",
272 ) as roundtrip:
273 subbox = Box.from_shape((3, 3), start=(11, -1))
274 subimage = roundtrip.get(bbox=subbox)
275 expected_subimage = self.extended_psf_image[subbox]
277 np.testing.assert_allclose(
278 subimage.image.array,
279 expected_subimage.image.array,
280 rtol=0.0,
281 atol=1e-8,
282 )
283 np.testing.assert_allclose(
284 subimage.variance.array,
285 expected_subimage.variance.array,
286 rtol=0.0,
287 atol=1e-8,
288 )
289 self.assertEqual(subimage.info, self.extended_psf_image.info)
290 self.assertEqual(subimage.fit, self.extended_psf_image.fit)
292 roundtripped = roundtrip.result
293 np.testing.assert_allclose(
294 roundtripped.image.array,
295 self.extended_psf_image.image.array,
296 rtol=0.0,
297 atol=1e-8,
298 )
299 np.testing.assert_allclose(
300 roundtripped.variance.array,
301 self.extended_psf_image.variance.array,
302 rtol=0.0,
303 atol=1e-8,
304 )
305 self.assertEqual(roundtripped.info, self.extended_psf_image.info)
306 self.assertEqual(roundtripped.fit, self.extended_psf_image.fit)
307 self.assertEqual(roundtripped.metadata["EPSF_TEST_KEY"], "EPSF_TEST VALUE")
309 def test_unit_sky_projection_bbox_properties(self):
310 """Test ExtendedPsfImage properties: unit, sky_projection, and bbox."""
311 self.assertEqual(self.extended_psf_image.unit, u.nJy)
312 self.assertIsNone(self.extended_psf_image.sky_projection)
313 self.assertEqual(self.extended_psf_image.bbox, self.extended_psf_image.image.bbox)
315 def test_copy_independence(self):
316 """Test that copy() produces an independent deep copy."""
317 copied = self.extended_psf_image.copy()
318 np.testing.assert_array_equal(copied.image.array, self.extended_psf_image.image.array)
319 np.testing.assert_array_equal(copied.variance.array, self.extended_psf_image.variance.array)
320 self.assertEqual(copied.info, self.extended_psf_image.info)
321 self.assertEqual(copied.fit, self.extended_psf_image.fit)
322 # Modifying the copy must not affect the original.
323 copied.image.array[:] = 0.0
324 self.assertFalse(np.all(self.extended_psf_image.image.array == 0.0))
326 def test_setitem_subregion_assignment(self):
327 """Test __setitem__ writes image and variance into a subregion."""
328 target = ExtendedPsfImage(Image(np.zeros((5, 6), dtype=np.float32), yx0=(10, -3), unit=u.nJy))
329 subbox = Box.from_shape((3, 3), start=(11, -1))
330 target[subbox] = self.extended_psf_image[subbox]
331 np.testing.assert_array_equal(
332 target[subbox].image.array,
333 self.extended_psf_image[subbox].image.array,
334 )
335 np.testing.assert_array_equal(
336 target[subbox].variance.array,
337 self.extended_psf_image[subbox].variance.array,
338 )
340 def test_constructor_invalid_inputs(self):
341 """Test that ValueError is raised for invalid constructor inputs."""
342 # Mismatched bboxes
343 image = Image(np.ones((5, 6), dtype=np.float32))
344 variance = Image(np.ones((4, 6), dtype=np.float32))
345 with self.assertRaises(ValueError):
346 ExtendedPsfImage(image=image, variance=variance)
348 # Variance has wrong units (nJy instead of nJy**2).
349 image = Image(np.ones((5, 6), dtype=np.float32), unit=u.nJy)
350 variance_wrong_units = Image(np.ones((5, 6), dtype=np.float32), unit=u.nJy)
351 with self.assertRaises(ValueError):
352 ExtendedPsfImage(image=image, variance=variance_wrong_units)
354 # Image has no units but variance does.
355 image_no_units = Image(np.ones((5, 6), dtype=np.float32))
356 variance_with_units = Image(np.ones((5, 6), dtype=np.float32), unit=u.nJy**2)
357 with self.assertRaises(ValueError):
358 ExtendedPsfImage(image=image_no_units, variance=variance_with_units)
361class ExtendedPsfFitTestCase(lsst.utils.tests.TestCase):
362 """Test ExtendedPsfFit."""
364 def test_construction(self):
365 """Test construction with and without optional fields."""
366 # With all fields
367 fit = ExtendedPsfFit(chi2=10.5, dof=10, reduced_chi2=1.05)
368 self.assertEqual(fit.chi2, 10.5)
369 self.assertEqual(fit.dof, 10)
370 self.assertEqual(fit.reduced_chi2, 1.05)
372 # With only the required field
373 fit_minimal = ExtendedPsfFit(chi2=10.5)
374 self.assertEqual(fit_minimal.chi2, 10.5)
375 self.assertIsNone(fit_minimal.dof)
376 self.assertIsNone(fit_minimal.reduced_chi2)
379class ExtendedPsfSubtractionTestCase(lsst.utils.tests.TestCase):
380 """Integration tests for all extended PSF subtraction pipeline tasks."""
382 @classmethod
383 def setUpClass(cls):
384 rng = np.random.default_rng(42)
386 # Background coefficients
387 sigma = 60.0 # noise
388 pedestal = 3210.1
389 coef_x = 1e-2
390 coef_y = 2e-2
391 coef_x2 = 1e-5
392 coef_xy = 2e-5
393 coef_y2 = 3e-5
395 # Make an input exposure
396 wcs = makeSkyWcs(
397 crpix=Point2D(0, 0),
398 crval=SpherePoint(0, 0, degrees),
399 cdMatrix=makeCdMatrix(scale=0.2 * arcseconds, flipX=True),
400 )
401 cls.exposure = ExposureF(1001, 1001, wcs)
402 cls.exposure.setPhotoCalib(makePhotoCalibFromCalibZeroPoint(10 ** (0.4 * 30), 1.0))
403 ny, nx = cls.exposure.image.array.shape
404 grid_y, grid_x = np.mgrid[(-ny + 1) // 2 : ny // 2 + 1, (-nx + 1) // 2 : nx // 2 + 1]
405 cls.exposure.image.array[:] += rng.normal(scale=sigma, size=cls.exposure.image.array.shape)
406 cls.background_array = np.zeros_like(cls.exposure.image.array)
407 cls.background_array += pedestal
408 cls.background_array += coef_x * grid_x
409 cls.background_array += coef_y * grid_y
410 cls.background_array += coef_x2 * grid_x * grid_x
411 cls.background_array += coef_xy * grid_x * grid_y
412 cls.background_array += coef_y2 * grid_y * grid_y
413 cls.exposure.image.array += cls.background_array
414 cls.exposure.info.setVisitInfo(VisitInfo(id=12345))
415 camera = CameraWrapper().camera
416 detector = camera[10]
417 cls.exposure.setDetector(detector)
418 for mask_plane in [
419 "BAD",
420 "CR",
421 "CROSSTALK",
422 "EDGE",
423 "NO_DATA",
424 "SAT",
425 "SUSPECT",
426 "UNMASKEDNAN",
427 "NEIGHBOR",
428 ]:
429 _ = cls.exposure.mask.addMaskPlane(mask_plane)
430 cls.exposure.variance.array.fill(1.0)
432 # Make a background list
433 background_image = ImageF(cls.exposure.getBBox())
434 background_image.array[:] = cls.background_array
435 background_control = BackgroundControl(Interpolate.NATURAL_SPLINE)
436 background_control.setNxSample(background_image.getWidth())
437 background_control.setNySample(background_image.getHeight())
438 background_model = makeBackground(background_image, background_control)
439 cls.background = BackgroundList()
440 cls.background.append(
441 (
442 background_model,
443 Interpolate.NATURAL_SPLINE,
444 REDUCE_INTERP_ORDER,
445 ApproximateControl.UNKNOWN,
446 0,
447 0,
448 False,
449 )
450 )
452 # Make a table of extended PSF candidate stars
453 corners = cls.exposure.wcs.pixelToSky([Point2D(x) for x in cls.exposure.getBBox().getCorners()])
454 corner_ras = [corner.getRa().asDegrees() for corner in corners]
455 corner_decs = [corner.getDec().asDegrees() for corner in corners]
456 num_stars = 10
457 ras = rng.uniform(np.min(corner_ras), np.max(corner_ras), num_stars)
458 decs = rng.uniform(np.min(corner_decs), np.max(corner_decs), num_stars)
459 sky_coords = [SpherePoint(ra, dec, degrees) for ra, dec in zip(ras, decs)]
460 pixel_coords = cls.exposure.wcs.skyToPixel(sky_coords)
461 pixel_x = [coord.getX() for coord in pixel_coords]
462 pixel_y = [coord.getY() for coord in pixel_coords]
463 mags = rng.uniform(10, 20, num_stars)
464 fluxes = [cls.exposure.photoCalib.magnitudeToInstFlux(mag) for mag in mags]
465 mm_coords = detector.transform(pixel_coords, PIXELS, FOCAL_PLANE)
466 mm_coords_x = np.array([mm_coord.x for mm_coord in mm_coords])
467 mm_coords_y = np.array([mm_coord.y for mm_coord in mm_coords])
468 radius_mm = np.sqrt(mm_coords_x**2 + mm_coords_y**2)
469 theta_radians = np.arctan2(mm_coords_y, mm_coords_x)
470 cls.extended_psf_candidate_table = Table(
471 {
472 "id": np.arange(num_stars),
473 "coord_ra": ras,
474 "coord_dec": decs,
475 "phot_g_mean_flux": fluxes,
476 "mag": mags,
477 "pixel_x": pixel_x,
478 "pixel_y": pixel_y,
479 "radius_mm": radius_mm,
480 "angle_radians": theta_radians,
481 }
482 )
484 # Make a synthetic star
485 cutout_radius = 25
486 grid_y, grid_x = np.mgrid[-cutout_radius : cutout_radius + 1, -cutout_radius : cutout_radius + 1]
487 dist_from_center = np.sqrt(grid_x**2 + grid_y**2)
488 sigma = 1.5
489 psf_array = np.exp(-(dist_from_center**2) / (2 * sigma**2))
490 psf_array /= np.sum(psf_array)
491 fixed_kernel = FixedKernel(ImageD(psf_array))
492 kernel_psf = KernelPsf(fixed_kernel)
493 cls.exposure.setPsf(kernel_psf)
494 psf = kernel_psf.computeKernelImage(kernel_psf.getAveragePosition())
496 # Add synthetic stars to the exposure
497 footprints = ImageF(cls.exposure.getBBox())
498 cls.injected_stars_array = np.zeros_like(cls.exposure.image.array)
499 for candidate_id, candidate in enumerate(cls.extended_psf_candidate_table):
500 bbox_star = Box2I(Point2I(candidate["pixel_x"], candidate["pixel_y"]), Extent2I(1, 1)).dilatedBy(
501 cutout_radius
502 )
503 bbox_star_clipped = bbox_star.clippedTo(cls.exposure.getBBox())
504 candidate_im = MaskedImageF(bbox_star)
505 candidate_im.image.array = candidate["phot_g_mean_flux"] * psf.getArray()
506 candidate_im = candidate_im[bbox_star_clipped]
507 detection_threshold = cls.exposure.getPhotoCalib().magnitudeToInstFlux(25)
508 detected = candidate_im.image.array > detection_threshold
509 footprints[bbox_star_clipped].array = detected * (candidate_id + 1)
510 _ = candidate_im.mask.addMaskPlane("DETECTED")
511 candidate_im.mask.array[detected] |= candidate_im.mask.getPlaneBitMask("DETECTED")
512 candidate_im.variance.array.fill(1.0)
513 y_slice = slice(bbox_star_clipped.getMinY(), bbox_star_clipped.getMaxY() + 1)
514 x_slice = slice(bbox_star_clipped.getMinX(), bbox_star_clipped.getMaxX() + 1)
515 cls.injected_stars_array[y_slice, x_slice] += candidate_im.image.array
516 cls.exposure.maskedImage[bbox_star_clipped] += candidate_im
517 cls.footprints = footprints.array
519 # Run the cutout task
520 extendedPsfCutoutConfig = ExtendedPsfCutoutConfig()
521 extendedPsfCutoutTask = ExtendedPsfCutoutTask(config=extendedPsfCutoutConfig)
522 cls.extended_psf_candidates = extendedPsfCutoutTask._get_extended_psf_candidates(
523 input_exposure=cls.exposure,
524 input_background=None,
525 footprints=cls.footprints,
526 extended_psf_candidate_table=cls.extended_psf_candidate_table,
527 )
529 # Run the stack task
530 extendedPsfStackConfig = ExtendedPsfStackConfig()
531 extendedPsfStackTask = ExtendedPsfStackTask(config=extendedPsfStackConfig)
532 stack_result = extendedPsfStackTask.run(extended_psf_candidates=cls.extended_psf_candidates)
533 cls.extended_psf = stack_result.extended_psf if stack_result is not None else None
535 @classmethod
536 def tearDownClass(cls):
537 del cls.exposure
538 del cls.background
539 del cls.background_array
540 del cls.injected_stars_array
541 del cls.extended_psf_candidate_table
542 del cls.footprints
543 del cls.extended_psf_candidates
544 del cls.extended_psf
546 def test_cutout_task_candidate_extraction(self):
547 """Test ExtendedPsfCutoutTask candidate extraction."""
548 assert self.extended_psf_candidates is not None
549 self.assertAlmostEqual(
550 float(self.extended_psf_candidates.metadata["FOCAL_PLANE_RADIUS_MM_MIN"]), 5.22408977, 7
551 )
552 self.assertAlmostEqual(
553 float(self.extended_psf_candidates.metadata["FOCAL_PLANE_RADIUS_MM_MAX"]), 14.6045832, 7
554 )
555 self.assertEqual(len(self.extended_psf_candidates), len(self.extended_psf_candidate_table))
556 self.assertEqual(self.extended_psf_candidates[0].star_info.visit, 12345)
557 self.assertEqual(self.extended_psf_candidates[0].star_info.detector, 10)
559 for candidate, candidate_entry in zip(
560 self.extended_psf_candidates, self.extended_psf_candidate_table
561 ):
562 self.assertEqual(candidate.star_info.ref_id, candidate_entry["id"])
563 self.assertEqual(candidate.star_info.ref_mag, candidate_entry["mag"])
564 self.assertEqual(candidate.star_info.position_x, candidate_entry["pixel_x"])
565 self.assertEqual(candidate.star_info.position_y, candidate_entry["pixel_y"])
566 self.assertIsInstance(candidate.psf_kernel_image, Image)
567 self.assertEqual(candidate.psf_kernel_image.array.ndim, 2)
568 self.assertGreater(candidate.psf_kernel_image.array.shape[0], 0)
569 self.assertGreater(candidate.psf_kernel_image.array.shape[1], 0)
570 self.assertTrue(np.isfinite(candidate.psf_kernel_image.array).all())
571 self.assertGreater(candidate.psf_kernel_image.array.sum(), 0.0)
572 focal_plane_radius = cast(u.Quantity, candidate.star_info.focal_plane_radius)
573 focal_plane_angle = cast(u.Quantity, candidate.star_info.focal_plane_angle)
574 self.assertEqual(focal_plane_radius.to_value(u.mm), candidate_entry["radius_mm"])
575 self.assertEqual(focal_plane_angle.to_value(u.rad), candidate_entry["angle_radians"])
577 def test_stack_task_moffat_fitting(self):
578 """Test Moffat fitting."""
579 assert self.extended_psf is not None
580 self.assertAlmostEqual(np.sum(self.extended_psf.image.array), 0.8233417, places=2)
581 self.assertAlmostEqual(np.sum(self.extended_psf.variance.array), 0.0075618913, places=4)
582 fit = self.extended_psf.fit
583 self.assertAlmostEqual(fit.chi2, 107652.97393353, delta=5)
584 self.assertEqual(fit.dof, 62996)
585 self.assertAlmostEqual(fit.reduced_chi2, 1.7088858647141, places=2)
586 self.assertAlmostEqual(fit.amplitude, 0.078900464260488, places=2)
587 self.assertAlmostEqual(fit.x_0, -0.68834523633912, places=2)
588 self.assertAlmostEqual(fit.y_0, -0.069005412739451, places=2)
589 self.assertAlmostEqual(fit.gamma, 8.0966823485900, places=2)
590 self.assertAlmostEqual(fit.alpha, 16.048683662812, places=2)
592 def test_stack_task_no_candidates(self):
593 """Test that None returned when no candidates pass the radius check."""
594 config = ExtendedPsfStackConfig()
595 config.min_focal_plane_radius = 1e6 # Excludes all candidates.
596 task = ExtendedPsfStackTask(config=config)
597 with self.assertLogs(level=logging.INFO):
598 result = task.run(extended_psf_candidates=self.extended_psf_candidates)
599 self.assertIsNone(result)
601 def test_subtract_task(self):
602 """Test subtraction task on all synthetic stars.
604 This test validates subtraction against known injected-star flux
605 while exercising default background restore and re-estimation.
606 """
608 assert self.extended_psf is not None
610 config = ExtendedPsfSubtractConfig()
611 config.bad_mask_planes = []
612 task = ExtendedPsfSubtractTask(config=config)
614 preliminary_visit_image = self.exposure.clone()
615 preliminary_visit_image.image.array -= self.background_array
616 restored_input_image = self.exposure.image.array.copy()
618 star_table = self.extended_psf_candidate_table
620 with unittest.mock.patch.object(task, "_get_subtraction_star_table", return_value=star_table):
621 result = task.run(
622 preliminary_visit_image=preliminary_visit_image,
623 preliminary_visit_image_background=self.background,
624 extended_psf=self.extended_psf,
625 ref_obj_loader=object(),
626 )
628 output_exposure = result.preliminary_visit_image_extended_psf_subtracted
629 output_background = result.preliminary_visit_image_extended_psf_subtracted_background
630 self.assertIsNotNone(output_background)
632 restored_output_image = output_exposure.image.array + output_background.getImage().array
633 removed_image = restored_input_image - restored_output_image
634 expected_flux = float(np.sum(self.injected_stars_array))
635 removed_flux = float(np.sum(removed_image))
636 self.assertFloatsAlmostEqual(removed_flux, expected_flux, rtol=0.25)
638 metadata = output_exposure.getMetadata()
639 n_stars = len(star_table)
640 self.assertEqual(metadata.getScalar("EPSFSUB_ATTEMPTED"), n_stars)
641 self.assertEqual(metadata.getScalar("EPSFSUB_SUBTRACTED"), n_stars)
644class MemoryTestCase(lsst.utils.tests.MemoryTestCase):
645 pass
648def setup_module(module):
649 lsst.utils.tests.init()
652if __name__ == "__main__": 652 ↛ 653line 652 didn't jump to line 653 because the condition on line 652 was never true
653 lsst.utils.tests.init()
654 unittest.main()