Coverage for python/lsst/images/_difference_image.py: 24%
207 statements
« prev ^ index » next coverage.py v7.16.2, created at 2026-10-14 19:06 +0000
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__ = ("DifferenceImage", "DifferenceImageSerializationModel", "DifferenceImageTemplateInfo")
16import logging
17import math
18import operator
19import uuid
20from collections.abc import Iterable, Mapping
21from types import EllipsisType
22from typing import TYPE_CHECKING, Any, ClassVar, Literal, cast
24import astropy.units
25import pydantic
26from astro_metadata_translator import ObservationInfo
28from ._backgrounds import BackgroundMap
29from ._geom import Bounds, Box, NoOverlapError
30from ._image import Image
31from ._mask import Mask, MaskPlane, MaskSchema, get_legacy_difference_image_mask_planes
32from ._observation_summary_stats import ObservationSummaryStats
33from ._polygon import Polygon
34from ._transforms import DetectorFrame, SkyProjection, TractFrame, Transform
35from ._visit_image import VisitImage, VisitImageSerializationModel
36from .aperture_corrections import (
37 ApertureCorrectionMap,
38)
39from .cameras import Detector
40from .convolution_kernels import ConvolutionKernel, ConvolutionKernelSerializationModel
41from .describe import DescribeOptions, FieldRole, Report, ReportField, ReportTable
42from .fields import Field
43from .psfs import (
44 PointSpreadFunction,
45)
46from .serialization import (
47 ArchiveReadError,
48 InputArchive,
49 InvalidParameterError,
50 MetadataValue,
51 OutputArchive,
52)
54if TYPE_CHECKING:
55 from lsst.daf.butler import DataId, DimensionRecord
57 try:
58 from lsst.afw.geom import SkyWcs as LegacySkyWcs
59 from lsst.afw.image import Exposure as LegacyExposure
60 from lsst.geom import Box2I as LegacyBox2I
61 from lsst.meas.algorithms import CoaddPsf as LegacyCoaddPsf
62 except ImportError:
63 type LegacyBox2I = Any # type: ignore[no-redef]
64 type LegacyExposure = Any # type: ignore[no-redef]
65 type LegacyCoaddPsf = Any # type: ignore[no-redef]
66 type LegacySkyWcs = Any # type: ignore[no-redef]
69class DifferenceImage(VisitImage):
70 """An image that is the PSF-matched difference of two other images.
72 Parameters
73 ----------
74 image
75 The main image plane. If this has a `SkyProjection`, it will be used
76 for all planes unless a ``sky_projection`` is passed separately.
77 mask
78 A bitmask image that annotates the main image plane. Must have the
79 same bounding box as ``image`` if provided. Any attached
80 ``sky_projection`` is replaced (possibly by `None`).
81 variance
82 The per-pixel uncertainty of the main image as an image of variance
83 values. Must have the same bounding box as ``image`` if provided, and
84 its units must be the square of ``image.unit`` or `None`.
85 Values default to ``1.0``. Any attached sky_projection is replaced
86 (possibly by `None`).
87 mask_schema
88 Schema for the mask plane. Must be provided if and only if ``mask`` is
89 not provided.
90 sky_projection
91 Projection that maps the pixel grid to the sky. Can only be `None` if
92 a ``sky_projection`` is already attached to ``image``.
93 bounds
94 The region where this image's pixels and other properties are valid.
95 If not provided, the bounding box of the image is used. Other
96 components (``psf``, ``sky_projection``, ``aperture_corrections``,
97 etc.) are assumed to have their own bounds which may or may not be the
98 same as the image bounds. If ``bounds`` extends beyond the image
99 bounding box, the intersection between ``bounds`` and the image
100 bounding box is used instead.
101 obs_info
102 General information about this visit in standardized form.
103 summary_stats
104 Summary statistics associated with this visit. Initialized to default
105 values if not provided.
106 photometric_scaling
107 Field that can be used to multiply a post-ISR image units to yield
108 calibrated image units. This may be a scaling that was already
109 applied (so dividing by it will recover the post-ISR units) or a
110 scaling that has not been applied, depending on ``image.unit``.
111 psf
112 Point-spread function model for this image, or an exception explaining
113 why it could not be read (to be raised if the PSF is requested later).
114 detector
115 Geometry and electronic information for the detector attached to this
116 image.
117 aperture_corrections : `dict` [`str`, `~fields.BaseField`]
118 Mapping from photometry algorithm name to the aperture correction for
119 that algorithm.
120 backgrounds
121 Background models associated with this image.
122 band
123 Name of the passband the image was observed with (this is a shorter,
124 less specific version of ``obs_info.physical_filter``).
125 kernel
126 The convolution kernel used to match the (warped) template to the
127 science image.
128 templates
129 Information about the template coadds that went into this difference
130 image.
131 metadata
132 Arbitrary flexible metadata to associate with the image.
134 Notes
135 -----
136 This class assumes that the difference has been performed on the pixel
137 grid of the 'science image' (i.e. a single observation, like `VisitImage`),
138 and most of the attributes of `DifferenceImage` correspond to the science
139 image. The 'template image' is assumed to be comprised of one or more
140 resampled coadd images stitched together.
142 The `DifferenceImage` class can also be used to represent the stitched
143 template itself; while this makes the naming a bit confusing, the type has
144 the right state to play this role.
145 """
147 def __init__(
148 self,
149 image: Image,
150 *,
151 mask: Mask | None = None,
152 variance: Image | None = None,
153 mask_schema: MaskSchema | None = None,
154 sky_projection: SkyProjection[DetectorFrame] | None = None,
155 bounds: Bounds | None = None,
156 obs_info: ObservationInfo | None = None,
157 summary_stats: ObservationSummaryStats | None = None,
158 photometric_scaling: Field | None = None,
159 psf: PointSpreadFunction | ArchiveReadError,
160 detector: Detector,
161 aperture_corrections: ApertureCorrectionMap | None = None,
162 backgrounds: BackgroundMap | None = None,
163 band: str,
164 kernel: ConvolutionKernel | None = None,
165 templates: Iterable[DifferenceImageTemplateInfo] | None = None,
166 metadata: dict[str, MetadataValue] | None = None,
167 ) -> None:
168 super().__init__(
169 image,
170 mask=mask,
171 variance=variance,
172 mask_schema=mask_schema,
173 sky_projection=sky_projection,
174 bounds=bounds,
175 obs_info=obs_info,
176 summary_stats=summary_stats,
177 photometric_scaling=photometric_scaling,
178 psf=psf,
179 detector=detector,
180 aperture_corrections=aperture_corrections,
181 backgrounds=backgrounds,
182 band=band,
183 metadata=metadata,
184 )
185 self._kernel = kernel
186 self._templates = list(templates) if templates is not None else None
188 @staticmethod
189 def _from_visit_image(
190 visit_image: VisitImage,
191 kernel: ConvolutionKernel | None,
192 templates: Iterable[DifferenceImageTemplateInfo] | None,
193 ) -> DifferenceImage:
194 return visit_image._transfer_metadata(
195 DifferenceImage(
196 visit_image.image,
197 mask=visit_image.mask,
198 variance=visit_image.variance,
199 sky_projection=visit_image.sky_projection,
200 bounds=visit_image.bounds,
201 obs_info=visit_image.obs_info,
202 summary_stats=visit_image.summary_stats,
203 photometric_scaling=visit_image.photometric_scaling,
204 psf=visit_image._psf, # get private attr to avoid triggering on ArchiveReadError early.
205 detector=visit_image.detector,
206 aperture_corrections=visit_image.aperture_corrections,
207 backgrounds=visit_image.backgrounds,
208 kernel=kernel,
209 templates=templates,
210 band=visit_image.band,
211 ),
212 )
214 @property
215 def kernel(self) -> ConvolutionKernel:
216 """The convolution kernel used to match the (warped) template
217 to the science image (`.convolution_kernels.ConvolutionKernel`).
218 """
219 if self._kernel is None:
220 raise AttributeError("This difference image does not have a kernel attached.")
221 return self._kernel
223 @kernel.setter
224 def kernel(self, kernel: ConvolutionKernel) -> None:
225 self._kernel = kernel
227 @kernel.deleter
228 def kernel(self) -> None:
229 self._kernel = None
231 @property
232 def templates(self) -> list[DifferenceImageTemplateInfo]:
233 """Information about the template coadds that went into this
234 difference image (`list` [`DifferenceImageTemplateInfo`]).
235 """
236 if self._templates is None: 236 ↛ 237line 236 didn't jump to line 237 because the condition on line 236 was never true
237 raise AttributeError("This difference image does not have any template information attached.")
238 return self._templates
240 @templates.setter
241 def templates(self, templates: Iterable[DifferenceImageTemplateInfo]) -> None:
242 self._templates = list(templates)
244 @templates.deleter
245 def templates(self) -> None:
246 self._templates = None
248 def __getitem__(self, bbox: Box | EllipsisType) -> DifferenceImage:
249 if bbox is ...:
250 return self
251 return self._from_visit_image(
252 super().__getitem__(bbox), kernel=self._kernel, templates=self._templates
253 )
255 def _describe(self, options: DescribeOptions = DescribeOptions(), /) -> Report:
256 """Return a `Report` describing this difference image.
258 Parameters
259 ----------
260 options : `DescribeOptions`, optional
261 Rendering options; forwarded to the base-class report.
262 """
263 report = super()._describe(options)
264 report.type_name = "DifferenceImage"
265 report.summary = f"DifferenceImage({self.image!s}, {list(self.mask.schema.names)})"
266 if options.brief:
267 return report
268 if self._templates: 268 ↛ 272line 268 didn't jump to line 272 because the condition on line 268 was always true
269 fields, table = DifferenceImageTemplateInfo._describe_templates(self._templates)
270 report.fields.extend(fields)
271 report.tables.append(table)
272 if self._kernel is not None: 272 ↛ 274line 272 didn't jump to line 274 because the condition on line 272 was always true
273 report.children["kernel"] = self._kernel._describe(options.for_child())
274 return report
276 def copy(self, *, copy_detector: bool = False) -> DifferenceImage:
277 """Deep-copy the difference image.
279 Parameters
280 ----------
281 copy_detector
282 Whether to deep-copy the `detector` attribute.
283 """
284 return self._from_visit_image(
285 super().copy(copy_detector=copy_detector), kernel=self._kernel, templates=self._templates
286 )
288 def convert_unit(
289 self,
290 unit: astropy.units.UnitBase = astropy.units.nJy,
291 copy: Literal["as-needed"] | bool = True,
292 copy_detector: bool = False,
293 ) -> DifferenceImage:
294 """Return an equivalent image with different pixel units.
296 Parameters
297 ----------
298 unit
299 The unit to transform to. This may be any of the following:
301 - any unit directly relatable to the current units via Astropy;
302 - any unit relatable to the product of the current units with the
303 `photometric_scaling` (i.e. if the current image is in
304 instrumental units but we know how to calibrate them)
305 - any unit relatable to the quotient of the current units with the
306 `photometric_scaling` (i.e. if the current image is in
307 calibrated units and we want to revert back to instrumental
308 units).
309 copy
310 Whether to copy the images and other components. If `True`, all
311 components that aren't controlled by some other argument will
312 always be deep-copied. If `False`, the operation will fail if the
313 image is not already in the right units. If ``as-needed``, only
314 the image and variance will be copied, and only if they are not
315 already in the right units.
316 copy_detector
317 Whether to deep-copy the `detector` attribute.
319 Returns
320 -------
321 `DifferenceImage`
322 An image with the given units.
323 """
324 return self._from_visit_image(
325 super().convert_unit(unit, copy=copy, copy_detector=copy_detector),
326 kernel=self._kernel,
327 templates=self._templates,
328 )
330 def serialize(self, archive: OutputArchive[Any]) -> DifferenceImageSerializationModel[Any]:
331 result = self._serialize_impl(DifferenceImageSerializationModel, archive)
332 if self._kernel is not None:
333 result.kernel = archive.serialize_direct("kernel", self._kernel.serialize)
334 else:
335 result.kernel = None
336 result.templates = self._templates
337 return result
339 @staticmethod
340 def _get_archive_tree_type[P: pydantic.BaseModel](
341 pointer_type: type[P],
342 ) -> type[DifferenceImageSerializationModel[P]]:
343 """Return the serialization model type for this object for an archive
344 type that uses the given pointer type.
345 """
346 return DifferenceImageSerializationModel[pointer_type] # type: ignore
348 @staticmethod
349 def from_legacy( # type: ignore[override]
350 legacy: LegacyExposure,
351 *,
352 exposure_record: DimensionRecord | None,
353 unit: astropy.units.UnitBase | None = None,
354 plane_map: Mapping[str, MaskPlane] | None = None,
355 ) -> DifferenceImage:
356 """Convert from an `lsst.afw.image.Exposure` instance.
358 Parameters
359 ----------
360 legacy
361 An `lsst.afw.image.Exposure` instance that will share image and
362 variance (but not mask) pixel data with the returned object.
363 exposure_record
364 The ``exposure`` dimension record for this observation, as a
365 source of additional required metadata. Pass `None` explicitly to
366 convert without some `~astro_metadata_translator.ObservationInfo`
367 fields (with a warning).
368 unit
369 Units of the image. If not provided, the ``BUNIT`` metadata
370 key will be used, if available.
371 plane_map
372 A mapping from legacy mask plane name to the new plane name and
373 description. If `None` (default),
374 `get_legacy_difference_image_mask_planes` is used.
375 """
376 if plane_map is None:
377 plane_map = get_legacy_difference_image_mask_planes()
378 return DifferenceImage._from_visit_image(
379 VisitImage.from_legacy(legacy, unit=unit, plane_map=plane_map, exposure_record=exposure_record),
380 kernel=None,
381 templates=None,
382 )
384 def to_legacy(
385 self, *, copy: bool | None = None, plane_map: Mapping[str, MaskPlane] | None = None
386 ) -> LegacyExposure:
387 """Convert to an `lsst.afw.image.Exposure` instance.
389 Parameters
390 ----------
391 copy
392 If `True`, always copy the image and variance pixel data.
393 If `False`, return a view, and raise `TypeError` if the pixel data
394 is read-only (this is not supported by afw). If `None`, only copy
395 if the pixel data is read-only. Mask pixel data is always copied.
396 plane_map
397 A mapping from legacy mask plane name to the new plane name and
398 description. If `None` (default),
399 `get_legacy_difference_image_mask_planes` is used.
400 """
401 if plane_map is None:
402 plane_map = get_legacy_difference_image_mask_planes()
403 return super().to_legacy(copy=copy, plane_map=plane_map)
405 @staticmethod
406 def read_legacy( # type: ignore[override]
407 filename: str,
408 *,
409 exposure_record: DimensionRecord | None,
410 preserve_quantization: bool = False,
411 plane_map: Mapping[str, MaskPlane] | None = None,
412 component: Literal[
413 "bbox",
414 "image",
415 "mask",
416 "variance",
417 "sky_projection",
418 "psf",
419 "detector",
420 "photometric_scaling",
421 "obs_info",
422 "summary_stats",
423 "aperture_corrections",
424 "bounds",
425 "unit",
426 ]
427 | None = None,
428 ) -> Any:
429 """Read a FITS file written by `lsst.afw.image.Exposure.writeFits`.
431 Parameters
432 ----------
433 filename
434 Full name of the file.
435 exposure_record
436 The ``exposure`` dimension record for this observation, as a
437 source of additional required metadata. Pass `None` explicitly to
438 convert without some `~astro_metadata_translator.ObservationInfo`
439 fields (with a warning).
440 preserve_quantization
441 If `True`, ensure that writing the masked image back out again will
442 exactly preserve quantization-compressed pixel values. This causes
443 the image and variance plane arrays to be marked as read-only and
444 stores the original binary table data for those planes in memory.
445 If the `MaskedImage` is copied, the precompressed pixel values are
446 not transferred to the copy.
447 plane_map
448 A mapping from legacy mask plane name to the new plane name and
449 description. If `None` (default),
450 `get_legacy_difference_image_mask_planes` is used.
451 component
452 A component to read instead of the full image.
453 """
454 if plane_map is None:
455 plane_map = get_legacy_difference_image_mask_planes()
456 result = VisitImage.read_legacy(
457 filename,
458 exposure_record=exposure_record,
459 preserve_quantization=preserve_quantization,
460 plane_map=plane_map,
461 component=component,
462 )
463 if component is None:
464 return DifferenceImage._from_visit_image(result, kernel=None, templates=None)
465 return result
468class DifferenceImageTemplateInfo(pydantic.BaseModel, ser_json_inf_nan="constants"):
469 """Information about how a template image contributed to a difference
470 image.
471 """
473 skymap: str = pydantic.Field(description="Name of the skymap that defines the tract/patch tiling.")
474 tract: int = pydantic.Field(description="ID of the tract (each tract is a different projection).")
475 patch: int = pydantic.Field(
476 description="ID of the patch (all patches within a tract share a projection)."
477 )
478 dataset_id: uuid.UUID = pydantic.Field(
479 description="Universally unique butler identifier for this template.",
480 )
481 dataset_run: str = pydantic.Field(description="Name of the butler RUN collection for this template.")
482 bounds: Polygon = pydantic.Field(
483 description=(
484 "The approximate intersection of the template and the science image, "
485 "in the science image's pixel coordinate system."
486 )
487 )
488 psf_shape_xx: float = pydantic.Field(description="Second moment of the effective PSF of the template.")
489 psf_shape_yy: float = pydantic.Field(description="Second moment of the effective PSF of the template.")
490 psf_shape_xy: float = pydantic.Field(description="Second moment of the effective PSF of the template.")
491 psf_shape_flag: bool = pydantic.Field(
492 description="Flag set if the second moments of the effective template PSF could not be computed."
493 )
495 @staticmethod
496 def from_legacy(
497 detector_frame: DetectorFrame,
498 legacy_template_psf: LegacyCoaddPsf,
499 legacy_template_metadata: Mapping[str, Any],
500 coadd_data_ids_by_uuid: Mapping[uuid.UUID, DataId],
501 coadd_dataset_type: str = "template_coadd",
502 log: logging.Logger | None = None,
503 ) -> list[DifferenceImageTemplateInfo]:
504 """Construct a list of template information structs from information
505 stored in a legacy stitched template image.
507 Parameters
508 ----------
509 detector_frame
510 Coordinate system and bounding box of the science image.
511 legacy_template_psf
512 The lazy-evaluation PSF model for the stitched template; used to
513 extract the tract and patch IDs of the coadds actually used and
514 their PSF models.
515 legacy_template_metadata
516 The FITS-style metadata of the stitched template; used to extract
517 butler UUIDs and RUN collection names for all *potential* input
518 coadds.
519 coadd_data_ids_by_uuid
520 A mapping from butler dataset ID to ``{tract, patch, band}`` data
521 ID for all coadds that may have contributed to the template. May
522 be a much larger superset of the needed datasets.
523 coadd_dataset_type
524 The name of the coadd template dataset type.
525 log
526 Logger to use for diagnostic messages.
528 Returns
529 -------
530 `list` [`DifferenceImageTemplateInfo`]
531 One struct for each input coadd that overlaps the science image,
532 sorted by tract and patch.
534 Raises
535 ------
536 KeyError
537 Raised if a PSF component's ``(tract, patch)`` has no matching
538 input coadd in the metadata, or if an input coadd's dataset ID is
539 missing from ``coadd_data_ids_by_uuid``.
540 RuntimeError
541 Raised if the input coadds come from more than one skymap.
543 See Also
544 --------
545 from_legacy_psf
546 """
547 n_inputs = legacy_template_metadata["LSST BUTLER N_INPUTS"]
548 butler_info: dict[tuple[int, int], tuple[uuid.UUID, str]] = {}
549 skymap: str | None = None
550 for n in range(n_inputs):
551 if legacy_template_metadata[f"LSST BUTLER INPUT {n} DATASETTYPE"] == coadd_dataset_type:
552 input_id = uuid.UUID(legacy_template_metadata[f"LSST BUTLER INPUT {n} ID"])
553 input_run = legacy_template_metadata[f"LSST BUTLER INPUT {n} RUN"]
554 input_data_id = coadd_data_ids_by_uuid[input_id]
555 if skymap is None:
556 skymap = cast(str, input_data_id["skymap"])
557 elif skymap != input_data_id["skymap"]:
558 raise RuntimeError("Cannot handle multiple skymaps in the inputs to a single template.")
559 butler_info[cast(int, input_data_id["tract"]), cast(int, input_data_id["patch"])] = (
560 input_id,
561 input_run,
562 )
563 return DifferenceImageTemplateInfo.from_legacy_psf(
564 detector_frame, legacy_template_psf, cast(str, skymap), butler_info, log=log
565 )
567 @staticmethod
568 def from_legacy_psf(
569 detector_frame: DetectorFrame,
570 legacy_template_psf: LegacyCoaddPsf,
571 skymap: str,
572 butler_info: Mapping[tuple[int, int], tuple[uuid.UUID, str]],
573 log: logging.Logger | None = None,
574 ) -> list[DifferenceImageTemplateInfo]:
575 """Construct a list of template information structs from a legacy
576 stitched template PSF and butler provenance held by the caller.
578 Parameters
579 ----------
580 detector_frame
581 Coordinate system and bounding box of the science image.
582 legacy_template_psf
583 The lazy-evaluation PSF model for the stitched template; used to
584 extract the tract and patch IDs of the coadds actually used and
585 their PSF models.
586 skymap
587 Name of the skymap that defines the tract and patch tiling.
588 butler_info
589 A mapping from ``(tract, patch)`` to the butler dataset ID and
590 RUN collection name of that coadd. May be a superset of the
591 coadds that contributed to the template.
592 log
593 Logger to use for diagnostic messages.
595 Returns
596 -------
597 `list` [`DifferenceImageTemplateInfo`]
598 One struct for each input coadd that overlaps the science image,
599 sorted by tract and patch.
601 Raises
602 ------
603 KeyError
604 Raised if a PSF component's ``(tract, patch)`` is missing from
605 ``butler_info``.
607 See Also
608 --------
609 from_legacy
611 Notes
612 -----
613 A task that builds the template has its input coadd references in
614 hand and can call this directly. `from_legacy` is for a caller that
615 has only the stored template, and reads the same provenance out of
616 its FITS metadata.
617 """
618 from lsst.afw.geom import makeWcsPairTransform
620 result: list[DifferenceImageTemplateInfo] = []
621 # A "component" of this PSF is an input {tract, patch} coadd.
622 for n in range(legacy_template_psf.getComponentCount()):
623 tract = legacy_template_psf.getTract(n)
624 patch = legacy_template_psf.getPatch(n)
625 dataset_id, dataset_run = butler_info[tract, patch]
626 patch_bbox = Box.from_legacy(legacy_template_psf.getBBox(n))
627 coadd_frame = TractFrame(
628 skymap=skymap,
629 tract=tract,
630 # This bbox is supposed to be the full tract bbox, but this
631 # frame is just a temporary and we don't have access to that.
632 # (If this ever becomes not-a-temporary, we could add a skymap
633 # argument).
634 bbox=patch_bbox,
635 )
636 detector_to_coadd = Transform.from_legacy(
637 makeWcsPairTransform(
638 # CoaddPsf method names did not anticipate being used for
639 # detector-level templates, so this is confusing:
640 legacy_template_psf.getCoaddWcs(), # this is the template_detector WCS!
641 legacy_template_psf.getWcs(n), # this is the template_coadd WCS!
642 ),
643 detector_frame,
644 coadd_frame,
645 )
646 coadd_to_detector = detector_to_coadd.inverted()
647 try:
648 # We transform the detector bbox to each coadd frame, do the
649 # intersection there, and then transform the intersection back
650 # to the detector frame, because we do not trust detector WCSs
651 # beyond the detector bounding box; they can be polynomials
652 # that extrapolate badly. Coadd WCSs in contrast are simple
653 # projections.
654 tmp_bounds = (
655 Polygon.from_box(detector_frame.bbox)
656 .transform(detector_to_coadd)
657 .intersection(patch_bbox)
658 ).transform(coadd_to_detector)
659 # Unfortunately doing the intersection in the coadd coordinate
660 # system means the transformed intersection might not quite be
661 # contained by the detector bounding box, due to floating-point
662 # round-off error. Intersect one more time to tidy it up.
663 bounds = tmp_bounds.intersection(detector_frame.bbox)
664 assert isinstance(bounds, Polygon), (
665 "The operations above should not change the region's fundamental topology."
666 )
667 except NoOverlapError:
668 if log is not None:
669 log.exception(
670 "No overlap with tract=%s, patch=%s; skipping provenance for that template.",
671 tract,
672 patch,
673 )
674 continue
675 try:
676 psf_shape = legacy_template_psf.computeShape(bounds.centroid.to_legacy_float_point())
677 except Exception:
678 if log is not None:
679 log.exception(
680 "Could not compute PSF shape for template coadd with tract=%s, patch=%s", tract, patch
681 )
682 else:
683 raise
684 psf_shape = None
685 result.append(
686 DifferenceImageTemplateInfo(
687 skymap=skymap,
688 tract=tract,
689 patch=patch,
690 dataset_id=dataset_id,
691 dataset_run=dataset_run,
692 bounds=bounds,
693 psf_shape_xx=psf_shape.getIxx() if psf_shape is not None else math.nan,
694 psf_shape_yy=psf_shape.getIyy() if psf_shape is not None else math.nan,
695 psf_shape_xy=psf_shape.getIxy() if psf_shape is not None else math.nan,
696 psf_shape_flag=psf_shape is None,
697 )
698 )
699 result.sort(key=lambda item: (item.tract, item.patch))
700 return result
702 _HOISTABLE_REPORT_VALUES: ClassVar[tuple[tuple[str, str, str], ...]] = (
703 ("skymap", "skymap", "Skymap"),
704 ("dataset_run", "template run", "Run"),
705 )
706 """Template attributes a report shows once above the templates table when
707 every template carries the same one, as ``(attribute, field label, column
708 header)``.
710 A value belongs here when it is wide enough to crowd the table and usually
711 uniform across a difference image's templates, and when it does not
712 identify the row: see the notes on `_describe_templates`. Adding an
713 attribute here is all it takes to cover it; a value that turns out to vary
714 falls back to its column on its own.
715 """
717 def _report_psf_sigma(self) -> str:
718 """Return the determinant radius of a template's effective PSF,
719 formatted for a report cell.
721 Returns
722 -------
723 sigma : `str`
724 Determinant radius in pixels, or ``"n/a"`` where the second moments
725 do not describe an ellipse.
727 Notes
728 -----
729 This is the same quantity as `ObservationSummaryStats.psfSigma`, so the
730 template PSFs can be compared directly against the science image's.
731 """
732 if self.psf_shape_flag:
733 return "n/a"
734 determinant = self.psf_shape_xx * self.psf_shape_yy - self.psf_shape_xy**2
735 if not math.isfinite(determinant) or determinant <= 0.0:
736 return "n/a"
737 return f"{determinant**0.25:.3f}"
739 @staticmethod
740 def _describe_templates(
741 templates: list[DifferenceImageTemplateInfo],
742 ) -> tuple[list[ReportField], ReportTable]:
743 """Return the report elements describing a difference image's
744 templates.
746 Parameters
747 ----------
748 templates
749 Templates to describe; must not be empty.
751 Returns
752 -------
753 fields : `list` [ `~lsst.images.ReportField` ]
754 The `_HOISTABLE_REPORT_VALUES` that every template shares.
755 table : `~lsst.images.ReportTable`
756 One row per template.
758 Notes
759 -----
760 Some of what a template carries is usually, but not dependably, the
761 same for every template of one difference image, and is long enough to
762 crowd out the columns that do differ. Each such value is shown once
763 above the table where every template shares it, and stays a column
764 where they do not, so nothing is lost when the usual case does not
765 hold. `_HOISTABLE_REPORT_VALUES` lists them.
767 What identifies a row is never hoisted, however uniform it happens to
768 be. A detector overlapping a single tract gives every template the
769 same tract, but tract and patch together name the coadd a row
770 describes, and splitting that pair between a field and a column would
771 leave each row unable to say what it is.
772 """
773 fields: list[ReportField] = []
774 columns: list[str] = []
775 getters: list[Any] = []
776 for attribute, label, column in DifferenceImageTemplateInfo._HOISTABLE_REPORT_VALUES:
777 getter = operator.attrgetter(attribute)
778 values = {getter(template) for template in templates}
779 if len(values) == 1:
780 fields.append(ReportField(label=label, value=values.pop(), role=FieldRole.DERIVED))
781 else:
782 columns.append(column)
783 getters.append(getter)
784 columns.extend(["Tract", "Patch", "PSF \N{GREEK SMALL LETTER SIGMA}"])
785 getters.extend([lambda t: t.tract, lambda t: t.patch, DifferenceImageTemplateInfo._report_psf_sigma])
786 if any(template.psf_shape_flag for template in templates):
787 columns.append("PSF flag")
788 getters.append(lambda t: "set" if t.psf_shape_flag else "")
789 columns.append("Dataset ID")
790 getters.append(lambda t: t.dataset_id)
791 return fields, ReportTable(
792 title="Templates",
793 columns=columns,
794 rows=[[getter(template) for getter in getters] for template in templates],
795 )
798class DifferenceImageSerializationModel[P: pydantic.BaseModel](VisitImageSerializationModel[P]):
799 """A Pydantic model used to represent a serialized `DifferenceImage`."""
801 SCHEMA_NAME: ClassVar[str] = "difference_image"
802 SCHEMA_VERSION: ClassVar[str] = "1.1.0"
803 MIN_READ_VERSION: ClassVar[int] = 1
804 PUBLIC_TYPE: ClassVar[type] = DifferenceImage
806 kernel: ConvolutionKernelSerializationModel | None = pydantic.Field(
807 description="The convolution kernel used to match the (warped) template to the science image."
808 )
809 templates: list[DifferenceImageTemplateInfo] | None = pydantic.Field(
810 description="Information about the template coadds that went into this difference image"
811 )
813 def deserialize(
814 self, archive: InputArchive[Any], *, bbox: Box | None = None, **kwargs: Any
815 ) -> DifferenceImage:
816 if kwargs: 816 ↛ 817line 816 didn't jump to line 817 because the condition on line 816 was never true
817 raise InvalidParameterError(f"Unrecognized parameters for DifferenceImage: {set(kwargs.keys())}.")
818 kernel = self.kernel.deserialize(archive) if self.kernel is not None else None
819 return DifferenceImage._from_visit_image(
820 super().deserialize(archive, bbox=bbox), kernel=kernel, templates=self.templates
821 )
823 def deserialize_component(self, component: str, archive: InputArchive[Any], **kwargs: Any) -> Any:
824 if kwargs and component not in ("image", "mask", "variance", "masked_image"):
825 raise InvalidParameterError(
826 f"Unsupported parameters for DifferenceImage component {component}: {set(kwargs.keys())}."
827 )
828 return super().deserialize_component(component, archive, **kwargs)