Coverage for python/lsst/images/cells/_provenance.py: 26%
146 statements
« prev ^ index » next coverage.py v7.14.1, created at 2026-06-08 01:51 -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
14__all__ = ("CoaddProvenance", "CoaddProvenanceSerializationModel")
16from collections.abc import Iterable
17from typing import TYPE_CHECKING, Any, ClassVar
19import astropy.table
20import astropy.units as u
21import numpy as np
22import pydantic
24from .._cell_grid import CellIJ
25from .._polygon import Polygon
26from ..serialization import ArchiveTree, InputArchive, InvalidParameterError, OutputArchive, TableModel
28if TYPE_CHECKING:
29 try:
30 from lsst.cell_coadds import CoaddInputs as LegacyCoaddInputs
31 from lsst.cell_coadds import MultipleCellCoadd, ObservationIdentifiers
32 from lsst.skymap import Index2D
33 except ImportError:
34 type Index2D = Any # type: ignore[no-redef]
35 type LegacyCoaddInputs = Any # type: ignore[no-redef]
36 type MultipleCellCoadd = Any # type: ignore[no-redef]
37 type ObservationIdentifiers = Any # type: ignore[no-redef]
40class CoaddProvenance:
41 """A pair of tables that record the inputs to a cell-based coadd.
43 Parameters
44 ----------
45 inputs
46 A table of {visit, detector} combinations that contribute to any cell
47 in the coadd.
48 contributions
49 A table of {visit, detector, cell} combinations that describe how an
50 observation contributed to a cell.
52 Notes
53 -----
54 This object can represent the provenance of a whole patch, a single cell,
55 or anything in between. In the single-cell case, the ``inputs`` and
56 ``contributions`` tables have the same number of rows (but may not be
57 ordered the same way!).
58 """
60 def __init__(self, inputs: astropy.table.Table, contributions: astropy.table.Table):
61 self._inputs = inputs
62 self._contributions = contributions
64 _INPUT_TABLE_COLUMNS: ClassVar[list[tuple[str, type, str]]] = [
65 ("instrument", np.object_, "Name of the instrument."),
66 ("visit", np.uint64, "ID of the visit."),
67 ("detector", np.uint16, "ID of the detector."),
68 ("physical_filter", np.object_, "Full name of the bandpass filter."),
69 ("day_obs", np.uint32, "Observation night as a YYYYMMDD integer."),
70 (
71 "polygon",
72 np.object_,
73 (
74 "Polygon that approximates the overlap of the observation and the coadd patch, "
75 "in coadd coordinates."
76 ),
77 ),
78 ]
80 _CONTRIBUTION_TABLE_COLUMNS: ClassVar[list[tuple[str, type, str, u.UnitBase | None]]] = [
81 ("cell_i", np.uint16, "Y-axis index of the cell within the patch.", None),
82 ("cell_j", np.uint16, "X-axis index of the cell within the patch.", None),
83 ("instrument", np.object_, "Name of the instrument.", None),
84 ("visit", np.uint64, "ID of the visit.", None),
85 ("detector", np.uint16, "ID of the detector.", None),
86 ("overlaps_center", np.bool_, "Whether a this observation overlaps the center of the cell.", None),
87 ("overlap_fraction", np.float64, "Fraction of the cell that is covered by the overlap region.", None),
88 ("weight", np.float64, "Weight to be used for this input in this cell.", None),
89 ("psf_shape_xx", np.float64, "Second order moments of the PSF.", u.pix**2),
90 ("psf_shape_yy", np.float64, "Second order moments of the PSF.", u.pix**2),
91 ("psf_shape_xy", np.float64, "Second order moments of the PSF.", u.pix**2),
92 (
93 "psf_shape_flag",
94 np.bool_,
95 "Flag indicating whether the PSF shape measurement was successful.",
96 None,
97 ),
98 ]
100 @classmethod
101 def make_empty_input_table(cls, n_rows: int) -> astropy.table.Table:
102 """Make an empty `inputs` table with a set number of rows."""
103 return astropy.table.Table(
104 [
105 astropy.table.Column(name=name, length=n_rows, dtype=dtype, description=description)
106 for name, dtype, description in cls._INPUT_TABLE_COLUMNS
107 ]
108 )
110 @classmethod
111 def make_empty_contribution_table(cls, n_rows: int) -> astropy.table.Table:
112 """Make an empty `contributions` table with a set number of rows."""
113 return astropy.table.Table(
114 [
115 astropy.table.Column(
116 name=name, length=n_rows, dtype=dtype, description=description, unit=unit
117 )
118 for name, dtype, description, unit in cls._CONTRIBUTION_TABLE_COLUMNS
119 ]
120 )
122 @property
123 def inputs(self) -> astropy.table.Table:
124 """A table of {visit, detector} combinations that contribute to any
125 cell in the coadd.
126 """
127 return self._inputs
129 @property
130 def contributions(self) -> astropy.table.Table:
131 """A table of {visit, detector, cell} combinations that describe how an
132 observation contributed to a cell.
133 """
134 return self._contributions
136 def __getitem__(self, cell: CellIJ) -> CoaddProvenance:
137 return self.subset([cell])
139 def subset(self, cells: Iterable[CellIJ]) -> CoaddProvenance:
140 """Return a new provenance object with just the given cells."""
141 cells_to_keep = astropy.table.Table(
142 rows=[(index.i, index.j) for index in cells],
143 names=["cell_i", "cell_j"],
144 dtype=[np.uint16, np.uint16],
145 )
146 contributions = astropy.table.join(self._contributions, cells_to_keep)
147 assert contributions.columns.keys() == {name for name, _, _, _ in self._CONTRIBUTION_TABLE_COLUMNS}
148 inputs = astropy.table.join(contributions["instrument", "visit", "detector"], self._inputs)
149 assert inputs.columns.keys() == {name for name, _, _ in self._INPUT_TABLE_COLUMNS}
150 return CoaddProvenance(inputs=inputs, contributions=contributions)
152 def serialize(self, archive: OutputArchive[Any]) -> CoaddProvenanceSerializationModel:
153 """Serialize the provenance to an output archive.
155 Parameters
156 ----------
157 archive
158 Archive to write to.
159 """
160 inputs = self._inputs.copy(copy_data=False)
161 contributions = self._contributions.copy(copy_data=False)
162 instrument = CoaddProvenanceSerializationModel._fix_str_for_serialization(
163 "instrument", inputs, contributions
164 )
165 physical_filter = CoaddProvenanceSerializationModel._fix_str_for_serialization(
166 "physical_filter", inputs
167 )
168 CoaddProvenanceSerializationModel._fix_polygon_for_serialization(inputs)
169 inputs_model = archive.add_table(inputs, name="inputs")
170 contributions_model = archive.add_table(contributions, name="contributions")
171 return CoaddProvenanceSerializationModel(
172 instrument=instrument,
173 physical_filter=physical_filter,
174 inputs=inputs_model,
175 contributions=contributions_model,
176 )
178 @staticmethod
179 def from_legacy(legacy_cell_coadd: MultipleCellCoadd) -> CoaddProvenance:
180 """Extract provenance from a legacy
181 `lsst.cell_coadds.MultipleCellCoadd` object.
182 """
183 inputs = CoaddProvenance.make_empty_input_table(len(legacy_cell_coadd.common.visit_polygons))
184 for n, (legacy_identifiers, legacy_polygon) in enumerate(
185 legacy_cell_coadd.common.visit_polygons.items()
186 ):
187 inputs["instrument"][n] = legacy_identifiers.instrument
188 inputs["visit"][n] = legacy_identifiers.visit
189 inputs["detector"][n] = legacy_identifiers.detector
190 inputs["physical_filter"][n] = legacy_identifiers.physical_filter
191 inputs["day_obs"][n] = legacy_identifiers.day_obs
192 inputs["polygon"][n] = Polygon.from_legacy(legacy_polygon)
193 n_contributions = 0
194 for legacy_cell in legacy_cell_coadd.cells.values():
195 n_contributions += len(legacy_cell.inputs)
196 contributions = CoaddProvenance.make_empty_contribution_table(n_contributions)
197 n = 0
198 for legacy_cell in legacy_cell_coadd.cells.values():
199 for legacy_identifiers, legacy_inputs in legacy_cell.inputs.items():
200 contributions["cell_i"][n] = legacy_cell.identifiers.cell.y
201 contributions["cell_j"][n] = legacy_cell.identifiers.cell.x
202 contributions["instrument"][n] = legacy_identifiers.instrument
203 contributions["visit"][n] = legacy_identifiers.visit
204 contributions["detector"][n] = legacy_identifiers.detector
205 contributions["overlaps_center"][n] = legacy_inputs.overlaps_center
206 contributions["overlap_fraction"][n] = legacy_inputs.overlap_fraction
207 contributions["weight"][n] = legacy_inputs.weight
208 contributions["psf_shape_xx"][n] = legacy_inputs.psf_shape.getIxx()
209 contributions["psf_shape_yy"][n] = legacy_inputs.psf_shape.getIyy()
210 contributions["psf_shape_xy"][n] = legacy_inputs.psf_shape.getIxy()
211 contributions["psf_shape_flag"][n] = legacy_inputs.psf_shape_flag
212 n += 1
213 return CoaddProvenance(inputs=inputs, contributions=contributions)
216class CoaddProvenanceSerializationModel(ArchiveTree):
217 """A Pydantic model used to represent a serialized `CoaddProvenance`.
219 Notes
220 -----
221 We can't rewrite the Astropy tables directly into the archive (e.g. as
222 FITS binary tables for a FITS archive), because:
224 - `str` columns are a huge pain in both Numpy and FITS;
225 - the polygon columns need to be rewritten as array-valued columns.
227 To deal with the string columns (``instrument`` and ``physical_filter``)
228 we do dictionary compression: we map each distinct value of those columns
229 to an integer, and then we save that mapping to the model while saving
230 an integer version of that column in the table. But if there is actually
231 only one value in that column (the most common case by far) we just drop
232 the column and store that value directly in the model.
233 """
235 SCHEMA_NAME: ClassVar[str] = "coadd_provenance"
236 SCHEMA_VERSION: ClassVar[str] = "1.0.0"
237 MIN_READ_VERSION: ClassVar[int] = 1
238 PUBLIC_TYPE: ClassVar[type] = CoaddProvenance
240 instrument: str | dict[str, int] = pydantic.Field(
241 description=(
242 "Instrument name for all inputs to this coadd, or a mapping from "
243 "instrument name to the integer used in its place in the tables."
244 )
245 )
246 physical_filter: str | dict[str, int] = pydantic.Field(
247 description="Physical filter name for all inputs to this coadd."
248 )
249 inputs: TableModel = pydantic.Field(description="Table of all inputs to the coadd.")
250 contributions: TableModel = pydantic.Field(description="Table of per-cell contributions to the coadd.")
252 def deserialize(self, archive: InputArchive[Any], **kwargs: Any) -> CoaddProvenance:
253 """Deserialize a provenance from an input archive.
255 Parameters
256 ----------
257 archive
258 Archive to read from.
260 Notes
261 -----
262 While `CoaddProvenance.subset` can be used to filter provenance
263 information down to just certain cells, there is no advantage to be
264 had from doing this during deserialization (the table data is not
265 ordered by cell, and hence there's read-slicing we can do).
266 """
267 if kwargs:
268 raise InvalidParameterError(f"Unrecognized parameters for CoaddProvenance: {set(kwargs.keys())}.")
269 inputs = archive.get_table(self.inputs)
270 contributions = archive.get_table(self.contributions)
271 CoaddProvenanceSerializationModel._fix_str_for_deserialization(
272 "instrument", self.instrument, inputs, contributions
273 )
274 CoaddProvenanceSerializationModel._fix_str_for_deserialization(
275 "physical_filter", self.physical_filter, inputs
276 )
277 CoaddProvenanceSerializationModel._fix_polygon_for_deserialization(inputs)
278 for name, _, description in CoaddProvenance._INPUT_TABLE_COLUMNS:
279 inputs.columns[name].description = description
280 for name, _, description, unit in CoaddProvenance._CONTRIBUTION_TABLE_COLUMNS:
281 contributions.columns[name].description = description
282 contributions.columns[name].unit = unit
283 return CoaddProvenance(inputs=inputs, contributions=contributions)
285 @staticmethod
286 def _fix_str_for_serialization(column: str, *tables: astropy.table.Table) -> str | dict[str, int]:
287 """Rewrite a string column as an integer column or drop it.
289 Parameters
290 ----------
291 column
292 Name of the column to rewrite.
293 *tables
294 One or more astropy tables to rewrite. The first table is assumed
295 to have all values for this column that might appear in any other
296 tables.
298 Returns
299 -------
300 `str` | `dict` [`str`, `int`]
301 If there is only one unique value for this column in the first
302 table, that value (and the column will have been dropped from
303 all givne tables). If the tables are empty, the column is
304 dropped and an empty `dict` is returned. In all other cases the
305 given column is replaced with an integer column in all given
306 tables and the mapping from strings to integers is returned.
307 """
308 result: str | dict[str, int] = {name: n for n, name in enumerate(sorted(set(tables[0][column])))}
309 match len(result):
310 case 0:
311 pass
312 case 1:
313 (result,) = result.keys() # type: ignore[union-attr]
314 case _:
315 for table in tables:
316 table.columns[column] = astropy.table.Column(
317 data=[result[k] for k in table.columns[column]],
318 name=column,
319 dtype=np.uint8,
320 description=f"Integer mapped to {column} name.",
321 )
322 return result
323 # If we didn't remap to an integer (case 0 and 1 above), delete the
324 # column.
325 for table in tables:
326 del table.columns[column]
327 return result
329 @staticmethod
330 def _fix_str_for_deserialization(
331 column: str, value: str | dict[str, int], *tables: astropy.table.Table
332 ) -> None:
333 """Rewrite an integer column back to a string one.
335 Parameters
336 ----------
337 column
338 Name of the column to rewrite.
339 value
340 Value or mapping of values returned by
341 `_fix_str_for_serialization`.
342 tables
343 Tables to rewrite this column in.
344 """
345 match value:
346 case str():
347 for table in tables:
348 table.columns[column] = astropy.table.Column([value] * len(table), dtype=object)
349 case dict():
350 mapping = {v: k for k, v in value.items()}
351 for table in tables:
352 table.columns[column] = astropy.table.Column(
353 [mapping[k] for k in table[column]], dtype=object
354 )
356 @staticmethod
357 def _fix_polygon_for_serialization(inputs: astropy.table.Table) -> None:
358 """Rewrite a polygon `object` column as a pair of array-valued columns
359 and an array-size column.
361 Parameters
362 ----------
363 inputs
364 A copy of the in-memory coadd inputs table to modify in-place into
365 its serialization form.
366 """
367 max_n_vertices = max(p.n_vertices for p in inputs["polygon"])
368 inputs["n_vertices"] = astropy.table.Column(
369 [p.n_vertices for p in inputs["polygon"]],
370 name="n_vertices",
371 dtype=np.uint8,
372 description="Number of polygon vertices.",
373 )
374 inputs["x_vertices"] = astropy.table.Column(
375 name="x_vertices",
376 dtype=np.float64,
377 length=len(inputs),
378 shape=(max_n_vertices,),
379 description="X coordinates of polygon vertices, in tract coordinates.",
380 )
381 inputs["x_vertices"][:, :] = np.nan
382 inputs["y_vertices"] = astropy.table.Column(
383 name="y_vertices",
384 dtype=np.float64,
385 length=len(inputs),
386 shape=(max_n_vertices,),
387 description="Y coordinates of polygon vertices, in tract coordinates.",
388 )
389 inputs["y_vertices"][:, :] = np.nan
390 for i, polygon in enumerate(inputs["polygon"]):
391 inputs["n_vertices"][i] = polygon.n_vertices
392 inputs["x_vertices"][i][: polygon.n_vertices] = polygon.x_vertices
393 inputs["y_vertices"][i][: polygon.n_vertices] = polygon.y_vertices
394 del inputs["polygon"]
396 @staticmethod
397 def _fix_polygon_for_deserialization(inputs: astropy.table.Table) -> None:
398 """Rewrite a a pair of array-valued columns and an array-size column
399 into a polygon `object` column.
401 Parameters
402 ----------
403 inputs
404 The serialized version of the coadd inputs table, to be modified
405 in-place into its in-memory form.
406 """
407 polygons = [
408 Polygon(x_vertices=x_vertices[:n_vertices], y_vertices=y_vertices[:n_vertices])
409 for n_vertices, x_vertices, y_vertices in zip(
410 inputs["n_vertices"], inputs["x_vertices"], inputs["y_vertices"]
411 )
412 ]
413 del inputs["n_vertices"]
414 del inputs["x_vertices"]
415 del inputs["y_vertices"]
416 inputs["polygon"] = astropy.table.Column(polygons, name="polygon", dtype=np.object_)