Coverage for python / lsst / meas / algorithms / brightStarStamps.py: 27%
168 statements
« prev ^ index » next coverage.py v7.13.5, created at 2026-04-23 08:32 +0000
1# This file is part of meas_algorithms.
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/>.
22"""Collection of small images (postage stamps) centered on bright stars."""
24from __future__ import annotations
26__all__ = ["BrightStarStamp", "BrightStarStamps"]
28from collections.abc import Sequence
29from dataclasses import dataclass
31import numpy as np
33from lsst.afw.detection import Psf
34from lsst.afw.fits import Fits, readMetadata
35from lsst.afw.geom import SkyWcs
36from lsst.afw.image import ImageFitsReader, MaskedImageF, MaskFitsReader
37from lsst.afw.table.io import InputArchive, OutputArchive
38from lsst.daf.base import PropertyList
39from lsst.geom import Angle, Point2D, degrees
40from lsst.meas.algorithms.stamps import AbstractStamp
41from lsst.utils.introspection import get_full_type_name
44@dataclass
45class BrightStarStamp(AbstractStamp):
46 """A single postage stamp centered on a bright star.
48 Attributes
49 ----------
50 stamp_im : `~lsst.afw.image.MaskedImageF`
51 The pixel data for this stamp.
52 psf : `~lsst.afw.detection.Psf`, optional
53 The point-spread function for this star.
54 wcs : `~lsst.afw.geom.SkyWcs`, optional
55 World coordinate system associated with the stamp.
56 visit : `int`, optional
57 Visit number of the observation.
58 detector : `int`, optional
59 Detector ID within the visit.
60 ref_id : `int`, optional
61 Reference catalog ID of the star.
62 ref_mag : `float`, optional
63 Reference catalog magnitude of the star.
64 position : `~lsst.geom.Point2D`, optional
65 Center position of the star on the detector in pixel coordinates.
66 focal_plane_radius : `float`, optional
67 Radial distance from the focal plane center in tangent-plane pixels.
68 focal_plane_angle : `~lsst.geom.Angle`, optional
69 Azimuthal angle on the focal plane (counterclockwise from +X).
70 scale : `float`, optional
71 Flux scaling factor applied to the PSF model.
72 scale_err : `float`, optional
73 Error in the flux scale.
74 pedestal : `float`, optional
75 Background pedestal level.
76 pedestal_err : `float`, optional
77 Error on the pedestal.
78 pedestal_scale_cov : `float`, optional
79 Covariance between pedestal and scale.
80 gradient_x : `float`, optional
81 Background gradient in the X direction.
82 gradient_y : `float`, optional
83 Background gradient in the Y direction.
84 curvature_x: `float`, optional
85 Background curvature in the X direction.
86 curvature_y: `float`, optional
87 Background curvature in the Y direction.
88 curvature_xy: `float`, optional
89 The xy component of the second order fit.
90 global_reduced_chi_squared : `float`, optional
91 Reduced chi-squared for the global model fit.
92 global_degrees_of_freedom : `int`, optional
93 Degrees of freedom for the global model fit.
94 psf_reduced_chi_squared : `float`, optional
95 Reduced chi-squared for the PSF fit.
96 psf_degrees_of_freedom : `int`, optional
97 Degrees of freedom for the PSF fit.
98 psf_masked_flux_fraction : `float`, optional
99 Fraction of flux masked in the PSF.
101 Notes
102 -----
103 This class is designed to be used with `BrightStarStamps`, which manages
104 collections of these stamps and handles reading/writing them to FITS files.
105 The `factory` class method provides a standard interface to construct
106 instances from image data and metadata, while the `_getMetadata` method
107 extracts metadata for storage in FITS headers.
108 """
110 stamp_im: MaskedImageF
111 psf: Psf | None
112 wcs: SkyWcs | None
113 visit: int | None
114 detector: int | None
115 ref_id: int | None
116 ref_mag: float | None
117 position: Point2D | None
118 focal_plane_radius: float | None
119 focal_plane_angle: Angle | None
120 scale: float | None
121 scale_err: float | None
122 pedestal: float | None
123 pedestal_err: float | None
124 pedestal_scale_cov: float | None
125 gradient_x: float | None
126 gradient_y: float | None
127 curvature_x: float | None
128 curvature_y: float | None
129 curvature_xy: float | None
130 global_reduced_chi_squared: float | None
131 global_degrees_of_freedom: int | None
132 psf_reduced_chi_squared: float | None
133 psf_degrees_of_freedom: int | None
134 psf_masked_flux_fraction: float | None
136 # Mapping of metadata keys to attribute names
137 _metadata_attribute_map = {
138 "VISIT": "visit",
139 "DETECTOR": "detector",
140 "REF_ID": "ref_id",
141 "REF_MAG": "ref_mag",
142 "POSITION_X": "position.x",
143 "POSITION_Y": "position.y",
144 "FOCAL_PLANE_RADIUS": "focal_plane_radius",
145 "FOCAL_PLANE_ANGLE_DEGREES": "focal_plane_angle",
146 "SCALE": "scale",
147 "SCALE_ERR": "scale_err",
148 "PEDESTAL": "pedestal",
149 "PEDESTAL_ERR": "pedestal_err",
150 "PEDESTAL_SCALE_COV": "pedestal_scale_cov",
151 "GRADIENT_X": "gradient_x",
152 "GRADIENT_Y": "gradient_y",
153 "CURVATURE_X": "curvature_x",
154 "CURVATURE_Y": "curvature_y",
155 "CURVATURE_XY": "curvature_xy",
156 "GLOBAL_REDUCED_CHI_SQUARED": "global_reduced_chi_squared",
157 "GLOBAL_DEGREES_OF_FREEDOM": "global_degrees_of_freedom",
158 "PSF_REDUCED_CHI_SQUARED": "psf_reduced_chi_squared",
159 "PSF_DEGREES_OF_FREEDOM": "psf_degrees_of_freedom",
160 "PSF_MASKED_FLUX_FRACTION": "psf_masked_flux_fraction",
161 }
163 def _getMetadata(self) -> PropertyList:
164 """Extract metadata from the stamp's attributes.
166 This method constructs a `PropertyList` containing metadata
167 extracted from the stamp's attributes. It is used when writing the
168 stamp to a FITS file to store relevant metadata in the FITS headers.
170 Returns
171 -------
172 metadata : `PropertyList`
173 A `PropertyList` containing the metadata, or `None` if no
174 metadata attributes are defined.
175 """
176 metadata = PropertyList()
177 for metadata_key, attribute_name in self._metadata_attribute_map.items():
178 if "." in attribute_name:
179 top_attr, sub_attr = attribute_name.split(".")
180 value = getattr(getattr(self, top_attr), sub_attr)
181 elif metadata_key == "FOCAL_PLANE_ANGLE_DEGREES":
182 value = getattr(self, attribute_name).asDegrees()
183 else:
184 value = getattr(self, attribute_name)
185 metadata[metadata_key] = value
186 return metadata
188 @property
189 def metadata(self) -> PropertyList:
190 """Return the stamp's metadata as a PropertyList."""
191 return self._getMetadata()
193 @classmethod
194 def factory(
195 cls,
196 stamp_im: MaskedImageF,
197 psf: Psf | None,
198 wcs: SkyWcs | None,
199 metadata: PropertyList,
200 ) -> BrightStarStamp:
201 """Construct a `BrightStarStamp` from image data and metadata.
203 This method provides a standard interface to create a `BrightStarStamp`
204 from its image data, PSF, WCS, and associated metadata.
205 It is used by the `BrightStarStamps.readFits` method to construct
206 individual bright star stamps from FITS files.
208 Parameters
209 ----------
210 stamp_im : `~lsst.afw.image.MaskedImageF`
211 Masked image for the stamp.
212 psf : `~lsst.afw.detection.Psf`, optional
213 Point-spread function for the stamp.
214 wcs : `~lsst.afw.geom.SkyWcs`, optional
215 World coordinate system for the stamp.
216 metadata : `PropertyList`
217 Metadata associated with the stamp, containing keys for all
218 required attributes.
220 Returns
221 -------
222 brightStarStamp : `BrightStarStamp`
223 The constructed `BrightStarStamp` instance.
224 """
225 kwargs = {}
227 for metadata_key, attribute_name in cls._metadata_attribute_map.items():
228 if "." in attribute_name: # for nested attributes like position.x
229 top_attr, sub_attr = attribute_name.split(".")
230 if top_attr not in kwargs: # avoid overwriting position
231 if top_attr == "position": # make an initial Point2D
232 kwargs[top_attr] = Point2D(0, 0)
233 setattr(kwargs[top_attr], sub_attr, metadata[metadata_key])
234 elif attribute_name == "focal_plane_angle":
235 kwargs[attribute_name] = Angle(metadata[metadata_key], degrees)
236 else:
237 kwargs[attribute_name] = metadata[metadata_key]
239 return cls(stamp_im=stamp_im, psf=psf, wcs=wcs, **kwargs)
242class BrightStarStamps(Sequence[BrightStarStamp]):
243 """A collection of bright star stamps.
245 Parameters
246 ----------
247 brightStarStamps : `Iterable` [`BrightStarStamp`]
248 Collection of `BrightStarStamp` instances.
249 metadata : `~lsst.daf.base.PropertyList`, optional
250 Global metadata associated with the collection.
251 """
253 def __init__(
254 self,
255 brightStarStamps: Sequence[BrightStarStamp],
256 metadata: PropertyList | None = None,
257 ):
258 self._stamps = list(brightStarStamps)
259 self._metadata = PropertyList() if metadata is None else metadata.deepCopy()
260 self.by_ref_id = {stamp.ref_id: stamp for stamp in self}
262 def __len__(self):
263 return len(self._stamps)
265 def __getitem__(self, index):
266 if isinstance(index, slice):
267 return BrightStarStamps(self._stamps[index], metadata=self._metadata)
268 return self._stamps[index]
270 def __iter__(self):
271 return iter(self._stamps)
273 @property
274 def metadata(self):
275 """Return the collection's global metadata as a PropertyList."""
276 return self._metadata
278 @classmethod
279 def readFits(cls, filename: str) -> BrightStarStamps:
280 """Make a `BrightStarStamps` object from a FITS file.
282 Parameters
283 ----------
284 filename : `str`
285 Name of the FITS file to read.
287 Returns
288 -------
289 brightStarStamps : `BrightStarStamps`
290 The constructed `BrightStarStamps` instance.
291 """
292 return cls.readFitsWithOptions(filename, None)
294 @classmethod
295 def readFitsWithOptions(cls, filename: str, options: PropertyList | None) -> BrightStarStamps:
296 """Make a `BrightStarStamps` object from a FITS file, with options.
298 Parameters
299 ----------
300 filename : `str`
301 Name of the FITS file to read.
302 options : `~lsst.daf.base.PropertyList`, optional
303 Options for reading the FITS file. Not currently used.
305 Returns
306 -------
307 brightStarStamps : `BrightStarStamps`
308 The constructed `BrightStarStamps` instance.
309 """
310 with Fits(filename, "r") as fits_file:
311 stamp_planes = {}
312 stamp_psf_ids = {}
313 stamp_wcs_ids = {}
314 stamp_metadata = {}
315 archive = None
317 for hdu_num in range(1, fits_file.countHdus()): # Skip primary HDU
318 metadata = readMetadata(filename, hdu=hdu_num)
319 extname = metadata["EXTNAME"]
320 stamp_id: int | None = metadata.get("EXTVER", None)
322 # Skip non-image BINTABLEs (except ARCHIVE_INDEX)
323 if metadata["XTENSION"] == "BINTABLE" and not metadata.get("ZIMAGE", False):
324 if extname != "ARCHIVE_INDEX":
325 continue
327 # Handle the archive index separately
328 if extname == "ARCHIVE_INDEX":
329 fits_file.setHdu(hdu_num)
330 archive = InputArchive.readFits(fits_file)
331 continue
332 elif metadata.get("EXTTYPE") == "ARCHIVE_DATA":
333 continue
335 # Select reader and dtype
336 if extname == "IMAGE":
337 reader = ImageFitsReader(filename, hdu=hdu_num)
338 dtype = np.dtype(MaskedImageF.dtype)
339 stamp_psf_ids[stamp_id] = metadata.pop("PSF", None)
340 stamp_wcs_ids[stamp_id] = metadata.pop("WCS", None)
341 stamp_metadata[stamp_id] = metadata
342 elif extname == "MASK":
343 reader = MaskFitsReader(filename, hdu=hdu_num)
344 dtype = None
345 elif extname == "VARIANCE":
346 reader = ImageFitsReader(filename, hdu=hdu_num)
347 dtype = np.dtype("float32")
348 else:
349 raise ValueError(f"Unknown extension type: {extname}")
351 if stamp_id is not None:
352 stamp_planes.setdefault(stamp_id, {})[extname.lower()] = reader.read(dtype=dtype)
354 primary_metadata = readMetadata(filename, hdu=0)
355 num_stamps = primary_metadata["N_STAMPS"]
357 if len(stamp_planes) != num_stamps:
358 raise ValueError(
359 f"Number of stamps read ({len(stamp_planes)}) does not agree with the "
360 f"number of stamps recorded in the primary HDU metadata ({num_stamps})."
361 )
362 if archive is None:
363 raise ValueError("No archive index was found in the FITS file; cannot read PSF or WCS.")
365 brightStarStamps = []
366 for stamp_id in range(1, num_stamps + 1): # Need to increment by one as EXTVER starts at 1
367 stamp = MaskedImageF(**stamp_planes[stamp_id])
368 psf = archive.get(stamp_psf_ids[stamp_id])
369 wcs = archive.get(stamp_wcs_ids[stamp_id])
370 brightStarStamps.append(BrightStarStamp.factory(stamp, psf, wcs, stamp_metadata[stamp_id]))
372 return cls(brightStarStamps, primary_metadata)
374 def writeFits(self, filename: str):
375 """Write this `BrightStarStamps` object to a FITS file.
377 Parameters
378 ----------
379 filename : `str`
380 Name of the FITS file to write.
381 """
382 metadata = self._metadata.deepCopy()
384 # Store metadata in the primary HDU
385 metadata["N_STAMPS"] = len(self._stamps)
386 metadata["VERSION"] = 2 # Record version number in case of future code changes
387 metadata["STAMPCLS"] = get_full_type_name(self)
389 # Create and write to the FITS file within a context manager
390 with Fits(filename, "w") as fits_file:
391 fits_file.createEmpty()
393 # Store Persistables in an OutputArchive
394 output_archive = OutputArchive()
395 stamp_psf_ids = []
396 stamp_wcs_ids = []
397 for stamp in self._stamps:
398 stamp_psf_ids.append(output_archive.put(stamp.psf))
399 stamp_wcs_ids.append(output_archive.put(stamp.wcs))
401 # Write to the FITS file
402 fits_file.writeMetadata(metadata)
403 del metadata
404 output_archive.writeFits(fits_file)
406 # Add all pixel data to extension HDUs; note: EXTVER should be 1-based
407 for stamp_id, (stamp, stamp_psf_id, stamp_wcs_id) in enumerate(
408 zip(self._stamps, stamp_psf_ids, stamp_wcs_ids),
409 start=1,
410 ):
411 metadata = PropertyList()
412 metadata.update({"EXTVER": stamp_id, "EXTNAME": "IMAGE"})
413 if stamp_metadata := stamp._getMetadata():
414 metadata.update(stamp_metadata)
415 metadata["PSF"] = stamp_psf_id
416 metadata["WCS"] = stamp_wcs_id
417 stamp.stamp_im.getImage().writeFits(filename, metadata=metadata, mode="a")
419 metadata = PropertyList()
420 metadata.update({"EXTVER": stamp_id, "EXTNAME": "MASK"})
421 stamp.stamp_im.getMask().writeFits(filename, metadata=metadata, mode="a")
423 metadata = PropertyList()
424 metadata.update({"EXTVER": stamp_id, "EXTNAME": "VARIANCE"})
425 stamp.stamp_im.getVariance().writeFits(filename, metadata=metadata, mode="a")