Coverage for tests / test_stitched_psf.py: 20%
123 statements
« prev ^ index » next coverage.py v7.13.5, created at 2026-04-23 08:47 +0000
1# This file is part of cell_coadds.
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 pickle
23import unittest
25import numpy as np
27import lsst.geom as geom
28import lsst.meas.base.tests
29import lsst.utils.tests
30from lsst.afw.detection import GaussianPsf, InvalidPsfError
31from lsst.afw.image import ImageD
32from lsst.cell_coadds import GridContainer, StitchedPsf, UniformGrid
33from lsst.meas.algorithms import SingleGaussianPsf
34from lsst.skymap import Index2D
37class StitchedPsfTestCase(lsst.utils.tests.TestCase):
38 """Test the methods of StitchedPsf class."""
40 @classmethod
41 def setUpClass(cls) -> None: # noqa: D102
42 super().setUpClass()
44 shape = Index2D(x=3, y=2)
46 cls.psf_size = 15
47 cls.psf_sigmas = {
48 Index2D(x=0, y=0): 1.2,
49 Index2D(x=0, y=1): 0.7,
50 Index2D(x=1, y=0): 0.9,
51 Index2D(x=1, y=1): 1.1,
52 Index2D(x=2, y=0): 1.3,
53 Index2D(x=2, y=1): 0.8,
54 }
56 gc = GridContainer[ImageD](shape)
57 for key, sigma in cls.psf_sigmas.items():
58 # It does not matter where we compute the kernel image.
59 # Pick any point.
60 gc[key] = GaussianPsf(cls.psf_size, cls.psf_size, sigma).computeKernelImage(geom.Point2D(1, 1))
62 grid = UniformGrid(cell_size=geom.Extent2I(cls.psf_size, cls.psf_size), shape=shape)
64 cls.psf = StitchedPsf(gc, grid)
66 cls.test_positions = (
67 (geom.Point2D(5, 4), Index2D(x=0, y=0)), # inner point in lower left
68 (geom.Point2D(6, 24), Index2D(x=0, y=1)), # inner point in upper left
69 (geom.Point2D(23.2, 7.8), Index2D(x=1, y=0)), # inner point in lower middle
70 (geom.Point2D(21, 22), Index2D(x=1, y=1)), # inner point in upper middle
71 (geom.Point2D(37, 9.4), Index2D(x=2, y=0)), # inner point in lower right
72 (geom.Point2D(42, 24), Index2D(x=2, y=1)), # inner point in upper right
73 # Some points that lie on the border
74 (geom.Point2D(31, 24), Index2D(x=2, y=1)), # inner point in upper right
75 (geom.Point2D(44, 0), Index2D(x=2, y=0)), # inner point in lower right
76 (geom.Point2D(17, 16), Index2D(x=1, y=1)), # inner point in upper middle
77 (geom.Point2D(15, 8), Index2D(x=1, y=0)), # inner point in lower middle
78 (geom.Point2D(0, 29), Index2D(x=0, y=1)), # inner point in upper left
79 (geom.Point2D(0, 0), Index2D(x=0, y=0)), # inner point in lower left
80 )
82 def test_resized(self):
83 """Test that the resized method works as it should."""
84 original_bbox = self.psf.computeBBox(self.psf.getAveragePosition())
85 # Resizing to original size should leave everything unchanged.
86 psf = self.psf.resized(self.psf_size, self.psf_size)
87 self.assertEqual(psf.getAveragePosition(), self.psf.getAveragePosition())
88 self.assertEqual(
89 psf.computeBBox(psf.getAveragePosition()), self.psf.computeBBox(self.psf.getAveragePosition())
90 )
91 self.assertImagesEqual(
92 psf.computeKernelImage(psf.getAveragePosition()),
93 self.psf.computeKernelImage(self.psf.getAveragePosition()),
94 )
95 self.assertEqual(psf, self.psf)
97 # Resize to a rectangular postage stamp
98 psf = self.psf.resized(25, 21)
99 self.assertEqual(self.psf.computeBBox(self.psf.getAveragePosition()), original_bbox)
100 self.assertEqual(psf.computeBBox(psf.getAveragePosition()).getDimensions(), geom.Extent2I(25, 21))
101 self.assertEqual(
102 psf.computeBBox(psf.getAveragePosition()).getCenter(),
103 self.psf.computeBBox(self.psf.getAveragePosition()).getCenter(),
104 )
105 self.assertNotEqual(psf, self.psf)
107 # Test that resizing to even dimensions throws an error
108 with self.assertRaises(ValueError):
109 psf = self.psf.resized(26, 22) # even, even
110 with self.assertRaises(ValueError):
111 psf = self.psf.resized(25, 22) # odd, even
112 with self.assertRaises(ValueError):
113 psf = self.psf.resized(26, 21) # even, odd
115 def test_clone(self):
116 """Test that the clone method works."""
117 psf = self.psf
118 cloned = psf.clone()
119 # Check that the cloned version is an exact replica of the original.
120 # self.assertEqual(psf, cloned) # Should there be an __eq__ method?
121 self.assertEqual(psf.getAveragePosition(), cloned.getAveragePosition())
122 self.assertEqual(
123 psf.computeBBox(psf.getAveragePosition()), cloned.computeBBox(cloned.getAveragePosition())
124 )
125 self.assertImagesEqual(
126 psf.computeKernelImage(psf.getAveragePosition()),
127 cloned.computeKernelImage(cloned.getAveragePosition()),
128 )
129 self.assertEqual(psf, cloned)
130 psf = psf.resized(41, 41)
131 # Now that one of them has been resized, they should not be equal.
132 self.assertNotEqual(
133 psf.computeBBox(psf.getAveragePosition()), cloned.computeBBox(cloned.getAveragePosition())
134 )
135 with self.assertRaises(TypeError):
136 self.assertImagesEqual(
137 psf.computeKernelImage(psf.getAveragePosition()),
138 cloned.computeKernelImage(cloned.getAveragePosition()),
139 )
140 self.assertNotEqual(psf, cloned)
142 def test_computeKernelImage(self):
143 """Test the computeKernelImage method for a StitchedPsf object."""
144 stitched_psf = self.psf
145 psf_bbox = geom.Box2I(
146 geom.Point2I(-(self.psf_size // 2), -(self.psf_size // 2)),
147 geom.Extent2I(self.psf_size, self.psf_size),
148 )
150 for position, cell_index in self.test_positions:
151 image1 = stitched_psf.computeKernelImage(position)
152 image2 = SingleGaussianPsf(
153 self.psf_size, self.psf_size, self.psf_sigmas[cell_index]
154 ).computeKernelImage(position)
155 # Small differences may exist due to differences in evaluating
156 # GaussianPsf vs. SingleGaussianPsf
157 self.assertImagesAlmostEqual(image1, image2, atol=2e-16)
158 self.assertEqual(image1.getBBox(), psf_bbox)
160 def test_computeBBox(self):
161 """Test the computeBBox method for a StitchedPsf object."""
162 psf = self.psf
163 psf_bbox = geom.Box2I(
164 geom.Point2I(-(self.psf_size // 2), -(self.psf_size // 2)),
165 geom.Extent2I(self.psf_size, self.psf_size),
166 )
168 for position, _ in self.test_positions:
169 bbox = psf.computeBBox(position)
170 self.assertEqual(bbox, psf_bbox)
172 def test_computeShape(self):
173 """Test the results from the computeShape method on a StitchedPsf
174 object matches the true input.
175 """
176 stitched_psf = self.psf
177 for position, cell_index in self.test_positions:
178 psf_shape = stitched_psf.computeShape(position) # check we can compute shape
179 self.assertIsNot(psf_shape.getIxx(), np.nan)
180 self.assertIsNot(psf_shape.getIyy(), np.nan)
181 self.assertIsNot(psf_shape.getIxy(), np.nan)
183 # Moments measured from pixellated images are significantly
184 # underestimated for small PSFs.
185 if self.psf_sigmas[cell_index] >= 1.0:
186 self.assertAlmostEqual(psf_shape.getIxx(), self.psf_sigmas[cell_index] ** 2, delta=1e-3)
187 self.assertAlmostEqual(psf_shape.getIyy(), self.psf_sigmas[cell_index] ** 2, delta=1e-3)
188 self.assertAlmostEqual(psf_shape.getIxy(), 0.0)
190 def test_computeApertureFlux(self):
191 """Test that the results from the computeApertureFlux method on a
192 StitchedPsf object returns the analytical results for a Gaussian PSF.
193 """
194 stitched_psf = self.psf
195 for position, cell_index in self.test_positions:
196 flux1sigma = stitched_psf.computeApertureFlux(self.psf_sigmas[cell_index], position=position)
197 self.assertAlmostEqual(flux1sigma, 0.39, delta=5e-2)
199 flux3sigma = stitched_psf.computeApertureFlux(
200 3.0 * self.psf_sigmas[cell_index], position=position
201 )
202 self.assertAlmostEqual(flux3sigma, 0.97, delta=2e-2)
204 def test_computeImage(self):
205 """Test the computeImage method for a StitchedPsf object produces
206 the same result as that on GaussianPsf for Gaussian PSFs.
207 """
208 stitched_psf = self.psf
209 psf_extent = geom.Extent2I(self.psf_size, self.psf_size)
211 for position, cell_index in self.test_positions:
212 image1 = stitched_psf.computeImage(position)
213 image2 = GaussianPsf(self.psf_size, self.psf_size, self.psf_sigmas[cell_index]).computeImage(
214 position
215 )
216 self.assertImagesEqual(image1, image2)
217 self.assertEqual(image1.getBBox().getDimensions(), psf_extent)
219 def test_computeImage_computeKernelImage(self):
220 """Test that computeImage called at integer points gives the same
221 result as calling computeKernelImage.
222 """
223 stitched_psf = self.psf
224 for position, _cell_index in self.test_positions:
225 pos = geom.Point2D(geom.Point2I(position)) # round to integer
226 image1 = stitched_psf.computeKernelImage(pos)
227 image2 = stitched_psf.computeImage(pos)
228 self.assertImagesEqual(image1, image2)
230 def test_pickle(self):
231 """Test that StitchedPsf objects can be pickled and unpickled."""
232 self.assertTrue(self.psf.isPersistable())
233 stream = pickle.dumps(self.psf)
234 psf = pickle.loads(stream)
235 self.assertEqual(psf.getAveragePosition(), self.psf.getAveragePosition())
236 self.assertEqual(
237 psf.computeBBox(psf.getAveragePosition()), self.psf.computeBBox(self.psf.getAveragePosition())
238 )
239 self.assertImagesEqual(
240 psf.computeKernelImage(psf.getAveragePosition()),
241 self.psf.computeKernelImage(self.psf.getAveragePosition()),
242 )
243 self.assertEqual(psf, self.psf)
245 def test_invalid_psf_error_no_inputs(self):
246 """Test that InvalidPsfError is raised when no inputs exist."""
247 stitched_psf = self.psf
249 # Test positions that are outside the grid coverage
250 # The grid covers indices (0,0), (0,1), (1,0), (1,1), (2,0), (2,1)
251 # Test position that maps to a non-existent grid cell
252 invalid_position = geom.Point2D(-10, -10) # Far outside the grid
253 err_msg = "No inputs exists at position."
254 # Test that computeBBox raises InvalidPsfError
255 with self.assertRaises(InvalidPsfError, msg=err_msg):
256 stitched_psf.computeBBox(invalid_position)
258 # Test that computeKernelImage raises InvalidPsfError
259 with self.assertRaises(InvalidPsfError, msg=err_msg):
260 stitched_psf.computeKernelImage(invalid_position)
263class TestMemory(lsst.utils.tests.MemoryTestCase):
264 """Test for memory/resource leaks."""
267def setup_module(module): # noqa: D103
268 lsst.utils.tests.init()
271if __name__ == "__main__": 271 ↛ 272line 271 didn't jump to line 272 because the condition on line 271 was never true
272 lsst.utils.tests.init()
273 unittest.main()