Coverage for python/lsst/images/cells/_coadd.py: 58%

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. 

11 

12from __future__ import annotations 

13 

14__all__ = ("CellCoadd", "CellCoaddSerializationModel") 

15 

16import functools 

17from collections.abc import Mapping, Sequence 

18from types import EllipsisType 

19from typing import TYPE_CHECKING, Any, ClassVar, cast 

20 

21import astropy.io.fits 

22import astropy.units 

23import astropy.wcs 

24import numpy as np 

25import pydantic 

26 

27from .._backgrounds import BackgroundMap, BackgroundMapSerializationModel 

28from .._cell_grid import CellGrid, CellGridBounds, PatchDefinition 

29from .._geom import YX, Box 

30from .._image import Image, ImageSerializationModel 

31from .._mask import ( 

32 Mask, 

33 MaskPlane, 

34 MaskSchema, 

35 MaskSerializationModel, 

36 get_legacy_deep_coadd_mask_planes, 

37 get_legacy_non_cell_coadd_mask_planes, 

38) 

39from .._masked_image import MaskedImage, MaskedImageSerializationModel 

40from .._transforms import SkyProjection, SkyProjectionSerializationModel, TractFrame 

41from ..describe import DescribeOptions, FieldRole, Report, ReportField 

42from ..fields import BaseField, field_from_legacy_photo_calib 

43from ..serialization import InputArchive, InvalidParameterError, OutputArchive 

44from ._aperture_corrections import ( 

45 CellApertureCorrectionMapSerializationModel, 

46 CellField, 

47 approximate_aperture_corrections_from_legacy, 

48) 

49from ._provenance import CoaddProvenance, CoaddProvenanceSerializationModel 

50from ._psf import CellPointSpreadFunction, CellPointSpreadFunctionSerializationModel 

51 

52if TYPE_CHECKING: 

53 try: 

54 from lsst.afw.image import Exposure as LegacyExposure 

55 from lsst.cell_coadds import MultipleCellCoadd as LegacyMultipleCellCoadd 

56 from lsst.skymap import TractInfo 

57 except ImportError: 

58 type LegacyExposure = Any # type: ignore[no-redef] 

59 type LegacyMultipleCellCoadd = Any # type: ignore[no-redef] 

60 type TractInfo = Any # type: ignore[no-redef] 

61 

62 

63class CellCoadd(MaskedImage): 

64 """A coadd comprised of cells on a regular grid. 

65 

66 Parameters 

67 ---------- 

68 image 

69 The main image plane. If this has a `.SkyProjection`, it will be used 

70 for all planes unless a ``sky_projection`` is passed separately. 

71 mask 

72 A bitmask image that annotates the main image plane. Must have the 

73 same bounding box as ``image`` if provided. Any attached 

74 ``sky_projection`` is replaced (possibly by `None`). 

75 variance 

76 The per-pixel uncertainty of the main image as an image of variance 

77 values. Must have the same bounding box as ``image`` if provided, and 

78 its units must be the square of ``image.unit`` or `None`. 

79 Values default to ``1.0``. Any attached ``sky_projection`` is replaced 

80 (possibly by `None`). 

81 mask_fractions 

82 A mapping from an input-image mask plane name to an image of the 

83 weights sums of that plane. 

84 noise_realizations 

85 A sequence of images with Monte Carlo realizations of the noise in 

86 the coadd. 

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 projection is already attached to ``image``. 

93 band 

94 Name of the band. 

95 psf 

96 Effective point-spread function for the coadd. The missing cells 

97 reported by ``psf.bounds`` are assumed to apply to all image data for 

98 that cell as well (i.e. there is a PSF for a cell if and only if 

99 there is image data for that cell). 

100 aperture_corrections 

101 Aperture corrections for different photometry algorithms. 

102 patch 

103 Identifiers and geometry of the full patch, if the image is confined 

104 to a single patch. When present, the cell grid of the PSF and 

105 provenance (if provideD) must be the full patch grid, even if its 

106 bounds select a subset of that area. 

107 provenance 

108 Information about the images that went into the coadd. 

109 backgrounds 

110 Background models associated with this image. 

111 """ 

112 

113 def __init__( 

114 self, 

115 image: Image, 

116 *, 

117 mask: Mask | None = None, 

118 variance: Image | None = None, 

119 mask_fractions: Mapping[str, Image] | None = None, 

120 noise_realizations: Sequence[Image] = (), 

121 mask_schema: MaskSchema | None = None, 

122 sky_projection: SkyProjection[TractFrame] | None = None, 

123 band: str | None = None, 

124 psf: CellPointSpreadFunction, 

125 aperture_corrections: Mapping[str, CellField] | None = None, 

126 patch: PatchDefinition | None = None, 

127 provenance: CoaddProvenance | None = None, 

128 backgrounds: BackgroundMap | None = None, 

129 ) -> None: 

130 super().__init__( 

131 image, 

132 mask=mask, 

133 variance=variance, 

134 mask_schema=mask_schema, 

135 sky_projection=sky_projection, 

136 ) 

137 if self.image.unit is None: 137 ↛ 138line 137 didn't jump to line 138 because the condition on line 137 was never true

138 raise TypeError("The image component of a CellCoadd must have units.") 

139 if self.image.sky_projection is None: 139 ↛ 140line 139 didn't jump to line 140 because the condition on line 139 was never true

140 raise TypeError("The sky_projection component of a CellCoadd cannot be None.") 

141 if not isinstance(self.image.sky_projection.pixel_frame, TractFrame): 141 ↛ 142line 141 didn't jump to line 142 because the condition on line 141 was never true

142 raise TypeError("The sky_projection's pixel frame must be a TractFrame for CellCoadd.") 

143 self._mask_fractions = dict(mask_fractions) if mask_fractions is not None else {} 

144 self._noise_realizations = list(noise_realizations) 

145 self._band = band 

146 if psf.bounds.bbox != self.bbox: 

147 psf = psf[self.bbox.intersection(psf.bounds.bbox)] 

148 self._psf = psf 

149 self._aperture_corrections = dict(aperture_corrections) if aperture_corrections is not None else {} 

150 for ap_corr_name, ap_corr_field in self._aperture_corrections.items(): 

151 if ap_corr_field.bounds.grid != self.grid: 151 ↛ 152line 151 didn't jump to line 152 because the condition on line 151 was never true

152 raise ValueError( 

153 f"Grids for cell PSF and {ap_corr_name} aperture corrections are not consistent." 

154 ) 

155 self._patch = patch 

156 self._provenance = provenance 

157 if self._provenance and not self._patch: 157 ↛ 158line 157 didn't jump to line 158 because the condition on line 157 was never true

158 raise TypeError("A CellCoadd cannot carry provenance without a patch definition.") 

159 self._backgrounds = backgrounds if backgrounds is not None else BackgroundMap() 

160 

161 @property 

162 def skymap(self) -> str: 

163 """Name of the skymap (`str`).""" 

164 return self.sky_projection.pixel_frame.skymap 

165 

166 @property 

167 def tract(self) -> int: 

168 """ID of the tract (`int`).""" 

169 return self.sky_projection.pixel_frame.tract 

170 

171 @property 

172 def patch(self) -> PatchDefinition: 

173 """Identifiers and geometry of the full patch, if the image is confined 

174 to a single patch (`PatchDefinition`). 

175 """ 

176 if self._patch is None: 

177 raise AttributeError("Coadd has no patch information.") 

178 return self._patch 

179 

180 @property 

181 def band(self) -> str | None: 

182 """Name of the band (`str` or `None`).""" 

183 return self._band 

184 

185 @property 

186 def mask_fractions(self) -> Mapping[str, Image]: 

187 """A mapping from an input-image mask plane name to an image of the 

188 weights sums of that plane 

189 (`~collections.abc.Mapping` [`str`, `.Image`]). 

190 """ 

191 return self._mask_fractions 

192 

193 @property 

194 def noise_realizations(self) -> Sequence[Image]: 

195 """A sequence of images with Monte Carlo realizations of the noise in 

196 the coadd (`~collections.abc.Sequence` [`.Image`]). 

197 """ 

198 return self._noise_realizations 

199 

200 @property 

201 def unit(self) -> astropy.units.UnitBase: 

202 """The units of the image plane (`astropy.units.Unit`).""" 

203 return cast(astropy.units.UnitBase, super().unit) 

204 

205 @property 

206 def sky_projection(self) -> SkyProjection[TractFrame]: 

207 """The projection that maps the pixel grid to the sky 

208 (`.SkyProjection` [`.TractFrame`]). 

209 """ 

210 return cast(SkyProjection[TractFrame], super().sky_projection) 

211 

212 @property 

213 def psf(self) -> CellPointSpreadFunction: 

214 """Effective point-spread function for the coadd 

215 (`CellPointSpreadFunction`). 

216 """ 

217 return self._psf 

218 

219 @property 

220 def aperture_corrections(self) -> Mapping[str, CellField]: 

221 """Aperture corrections for different photometry algorithms 

222 (`dict` [`str`, `CellField`]). 

223 """ 

224 return self._aperture_corrections 

225 

226 @property 

227 def bounds(self) -> CellGridBounds: 

228 """The grid of cells that overlap this coadd and a set of missing 

229 cells (`CellGridBounds`). 

230 """ 

231 return self._psf.bounds 

232 

233 @property 

234 def grid(self) -> CellGrid: 

235 """The grid of cells that overlap this coadd (`CellGrid`).""" 

236 return self._psf.bounds.grid 

237 

238 @property 

239 def provenance(self) -> CoaddProvenance: 

240 """Information about the images that went into the coadd 

241 (`CoaddProvenance` or `None`). 

242 """ 

243 if self._provenance is None: 

244 raise AttributeError("Coadd has no provenance information.") 

245 return self._provenance 

246 

247 @property 

248 def backgrounds(self) -> BackgroundMap: 

249 """A mapping of backgrounds associated with this image 

250 (`.BackgroundMap`). 

251 """ 

252 return self._backgrounds 

253 

254 def __getitem__(self, bbox: Box | EllipsisType) -> CellCoadd: 

255 bbox, _ = self._handle_getitem_args(bbox) 

256 psf = self.psf[bbox] 

257 return self._transfer_metadata( 

258 CellCoadd( 

259 self.image[bbox], 

260 mask=self.mask[bbox], 

261 variance=self.variance[bbox], 

262 sky_projection=self.sky_projection, 

263 mask_fractions={k: v[bbox] for k, v in self._mask_fractions.items()}, 

264 noise_realizations=[v[bbox] for v in self._noise_realizations], 

265 band=self.band, 

266 psf=psf, 

267 patch=self._patch, 

268 provenance=( 

269 self._provenance.subset(psf.bounds.cell_indices()) 

270 if self._provenance is not None 

271 else None 

272 ), 

273 backgrounds=self._backgrounds, 

274 aperture_corrections=self._aperture_corrections.copy(), 

275 ), 

276 bbox=bbox, 

277 ) 

278 

279 def _describe(self, options: DescribeOptions = DescribeOptions(), /) -> Report: 

280 """Return a `Report` describing this cell coadd. 

281 

282 Parameters 

283 ---------- 

284 options : `DescribeOptions`, optional 

285 Rendering options; forwarded to all children. 

286 """ 

287 fields = [ 

288 ReportField(label="skymap", value=self.skymap, role=FieldRole.DERIVED), 

289 ReportField(label="tract", value=self.tract, role=FieldRole.DERIVED), 

290 ReportField(label="patch", value=self._patch, role=FieldRole.DERIVED), 

291 ReportField(label="band", value=self.band, role=FieldRole.DERIVED), 

292 ReportField(label="bbox", value=self.bbox, repr_value=repr(self.bbox), role=FieldRole.DERIVED), 

293 ] 

294 summary = f"CellCoadd({self.bbox!s}, tract={self.tract})" 

295 if options.brief: 

296 return Report(type_name="CellCoadd", summary=summary, fields=fields) 

297 # Components with no report of their own get a line each, so that 

298 # nothing the coadd carries is absent from its description. 

299 if self._mask_fractions: 

300 fields.append( 

301 ReportField( 

302 label="mask_fractions", value=", ".join(self._mask_fractions), role=FieldRole.DERIVED 

303 ) 

304 ) 

305 if self._noise_realizations: 

306 n_noise = len(self._noise_realizations) 

307 fields.append( 

308 ReportField( 

309 label="noise_realizations", 

310 value=( 

311 f"{n_noise} image{'s' if n_noise != 1 else ''} " 

312 f"(dtype {self._noise_realizations[0].array.dtype})" 

313 ), 

314 role=FieldRole.DERIVED, 

315 ) 

316 ) 

317 if self._aperture_corrections: 

318 n_ap_corr = len(self._aperture_corrections) 

319 # Never print all the keys even in detailed mode since there 

320 # are usually many of them and they swamp the report. 

321 ap_corr = f"{n_ap_corr} field{'s' if n_ap_corr != 1 else ''}" 

322 fields.append(ReportField(label="aperture_corrections", value=ap_corr, role=FieldRole.DERIVED)) 

323 if self._provenance is None: 

324 # A coadd without provenance is a meaningful state: it is what 

325 # reading with provenance=False gives, and it is what makes 

326 # to_legacy_cell_coadd fail. 

327 fields.append(ReportField(label="provenance", value="none", role=FieldRole.DERIVED)) 

328 child = options.for_child() 

329 plane = options.for_child("sky_projection", "bbox") 

330 children: dict[str, Report] = { 

331 "image": self.image._describe(plane), 

332 "mask": self.mask._describe(plane), 

333 "variance": self.variance._describe(plane), 

334 "sky_projection": self.sky_projection._describe(child, bbox=self.bbox), 

335 "psf": self.psf._describe(child), 

336 } 

337 if self._provenance is not None: 

338 # The provenance property raises when there is none. Its cells 

339 # come from this coadd, so pass them down for the coverage ratio. 

340 children["provenance"] = self._provenance._describe(child, bounds=self.bounds) 

341 children["backgrounds"] = self.backgrounds._describe(child) 

342 return Report( 

343 type_name="CellCoadd", 

344 summary=summary, 

345 fields=fields, 

346 children=children, 

347 ) 

348 

349 def copy(self) -> CellCoadd: 

350 """Deep-copy the coadd.""" 

351 return self._transfer_metadata( 

352 CellCoadd( 

353 image=self._image.copy(), 

354 mask=self._mask.copy(), 

355 variance=self._variance.copy(), 

356 sky_projection=self.sky_projection, 

357 mask_fractions={k: v.copy() for k, v in self._mask_fractions.items()}, 

358 noise_realizations=[v.copy() for v in self._noise_realizations], 

359 band=self.band, 

360 psf=self.psf, 

361 patch=self._patch, 

362 provenance=self._provenance, 

363 backgrounds=self._backgrounds.copy(), 

364 aperture_corrections=self._aperture_corrections.copy(), 

365 ), 

366 copy=True, 

367 ) 

368 

369 def apply_background(self, name: str | None) -> None: 

370 """Subtract the background with the given name, modifying the image 

371 in place. 

372 

373 If ``name`` is `None`, restore the original background. 

374 

375 Parameters 

376 ---------- 

377 name 

378 Name of the background to subtract, or `None` to restore the 

379 original background. 

380 """ 

381 current_bg = self.backgrounds.subtracted 

382 if current_bg is not None and name == current_bg.name: 382 ↛ 383line 382 didn't jump to line 383 because the condition on line 382 was never true

383 return 

384 

385 adjustment: astropy.units.Quantity | None = None 

386 if current_bg is not None: 

387 adjustment = current_bg.field.render(self.bbox, dtype=self.image.array.dtype).quantity 

388 if name is not None: 

389 new_bg = self.backgrounds[name] 

390 new_values = new_bg.field.render(self.bbox, dtype=self.image.array.dtype).quantity 

391 adjustment = -new_values if adjustment is None else adjustment - new_values 

392 

393 # Resolve and render every model before changing either pixels or 

394 # state. Computing one combined adjustment also ensures a failed unit 

395 # conversion cannot leave a previously subtracted background restored 

396 # while the map still says it is subtracted. 

397 if adjustment is not None: 397 ↛ 399line 397 didn't jump to line 399 because the condition on line 397 was always true

398 self.image.quantity += adjustment 

399 self._backgrounds._subtracted = name 

400 

401 def serialize(self, archive: OutputArchive[Any]) -> CellCoaddSerializationModel: 

402 """Serialize the image to an output archive. 

403 

404 Parameters 

405 ---------- 

406 archive 

407 Archive to write to. 

408 """ 

409 serialized_image = archive.serialize_direct( 

410 "image", 

411 functools.partial(self.image.serialize, save_projection=False, tile_shape=self.grid.cell_shape), 

412 ) 

413 serialized_mask = archive.serialize_direct( 

414 "mask", 

415 functools.partial(self.mask.serialize, save_projection=False, tile_shape=self.grid.cell_shape), 

416 ) 

417 serialized_variance = archive.serialize_direct( 

418 "variance", 

419 functools.partial( 

420 self.variance.serialize, save_projection=False, tile_shape=self.grid.cell_shape 

421 ), 

422 ) 

423 serialized_projection = archive.serialize_direct("sky_projection", self.sky_projection.serialize) 

424 serialized_mask_fractions = { 

425 k: archive.serialize_direct( 

426 f"mask_fractions/{k}", 

427 functools.partial( 

428 v.serialize, 

429 save_projection=False, 

430 tile_shape=self.grid.cell_shape, 

431 options_name="mask_fractions", 

432 ), 

433 ) 

434 for k, v in self.mask_fractions.items() 

435 } 

436 serialized_noise_realizations = [ 

437 archive.serialize_direct( 

438 f"noise_realizations/{n}", 

439 functools.partial( 

440 v.serialize, save_projection=False, tile_shape=self.grid.cell_shape, options_name="image" 

441 ), 

442 ) 

443 for n, v in enumerate(self.noise_realizations) 

444 ] 

445 serialized_psf = archive.serialize_direct("psf", self.psf.serialize) 

446 serialized_aperture_corrections = archive.serialize_direct( 

447 "aperture_corrections", 

448 functools.partial( 

449 CellApertureCorrectionMapSerializationModel.serialize, self.aperture_corrections 

450 ), 

451 ) 

452 serialized_provenance = ( 

453 archive.serialize_direct("provenance", self._provenance.serialize) 

454 if self._provenance is not None 

455 else None 

456 ) 

457 serialized_backgrounds = archive.serialize_direct("background", self._backgrounds.serialize) 

458 return CellCoaddSerializationModel( 

459 image=serialized_image, 

460 mask=serialized_mask, 

461 variance=serialized_variance, 

462 sky_projection=serialized_projection, 

463 mask_fractions=serialized_mask_fractions, 

464 noise_realizations=serialized_noise_realizations, 

465 band=self._band, 

466 psf=serialized_psf, 

467 aperture_corrections=serialized_aperture_corrections, 

468 patch=self._patch, 

469 provenance=serialized_provenance, 

470 backgrounds=serialized_backgrounds, 

471 metadata=self._metadata, 

472 ) 

473 

474 @staticmethod 

475 def _get_archive_tree_type[P: pydantic.BaseModel]( 

476 pointer_type: type[P], 

477 ) -> type[CellCoaddSerializationModel[P]]: 

478 """Return the serialization model type for this object for an archive 

479 type that uses the given pointer type. 

480 """ 

481 return CellCoaddSerializationModel[pointer_type] # type: ignore 

482 

483 @staticmethod 

484 def from_legacy_cell_coadd( 

485 legacy: LegacyMultipleCellCoadd, 

486 *, 

487 plane_map: Mapping[str, MaskPlane] | None = None, 

488 tract_info: TractInfo, 

489 bbox: Box | None = None, 

490 ) -> CellCoadd: 

491 """Convert from a `lsst.cell_coadds.MultipleCellCoadd` instance. 

492 

493 Parameters 

494 ---------- 

495 legacy 

496 A `lsst.cell_coadds.MultipleCellCoadd` instance to convert. 

497 plane_map 

498 A mapping from legacy mask plane name to the new plane name and 

499 description. 

500 tract_info 

501 Information about the full tract. 

502 bbox 

503 Bounding box of the image. The default is to include just the 

504 bounding box of the valid cells, which may not cover a full patch. 

505 """ 

506 from lsst.geom import Box2I 

507 

508 if plane_map is None: 

509 plane_map = get_legacy_deep_coadd_mask_planes() 

510 if bbox is None: 

511 legacy_bbox = Box2I() 

512 for single_cell in legacy.cells.values(): 

513 legacy_bbox.include(single_cell.inner.bbox) 

514 else: 

515 legacy_bbox = bbox.to_legacy() 

516 legacy_stitched = legacy.stitch(legacy_bbox) 

517 unit = astropy.units.Unit(legacy.units.value) 

518 tract_bbox = Box.from_legacy(tract_info.getBBox()) 

519 sky_projection = SkyProjection.from_legacy( 

520 legacy.wcs, 

521 TractFrame( 

522 skymap=legacy.identifiers.skymap, 

523 tract=legacy.identifiers.tract, 

524 bbox=tract_bbox, 

525 ), 

526 pixel_bounds=tract_bbox, 

527 ) 

528 band = legacy.identifiers.band 

529 image = Image.from_legacy(legacy_stitched.image, unit=unit) 

530 mask = Mask.from_legacy(legacy_stitched.mask, plane_map=plane_map) 

531 variance = Image.from_legacy(legacy_stitched.variance, unit=unit**2) 

532 noise_realizations = [ 

533 Image.from_legacy(noise_image) for noise_image in legacy_stitched.noise_realizations 

534 ] 

535 mask_fractions = ( 

536 {"rejected": Image.from_legacy(legacy_stitched.mask_fractions)} 

537 if legacy_stitched.mask_fractions is not None 

538 else {} 

539 ) 

540 psf = CellPointSpreadFunction.from_legacy(legacy_stitched.psf, image.bbox) 

541 aperture_corrections = { 

542 ap_corr_name: CellField.from_legacy_aperture_correction(legacy_ap_corr, psf.bounds) 

543 for ap_corr_name, legacy_ap_corr in legacy_stitched.ap_corr_map.items() 

544 } 

545 patch_info = tract_info[legacy.identifiers.patch] 

546 patch = PatchDefinition( 

547 id=patch_info.getSequentialIndex(), 

548 index=YX(y=legacy.identifiers.patch.y, x=legacy.identifiers.patch.x), 

549 inner_bbox=Box.from_legacy(patch_info.getInnerBBox()), 

550 cells=CellGrid.from_legacy(legacy.grid), 

551 ) 

552 provenance = CoaddProvenance.from_legacy(legacy) 

553 return CellCoadd( 

554 image=image, 

555 mask=mask, 

556 variance=variance, 

557 mask_fractions=mask_fractions, 

558 noise_realizations=noise_realizations, 

559 sky_projection=sky_projection, 

560 band=band, 

561 psf=psf, 

562 aperture_corrections=aperture_corrections, 

563 patch=patch, 

564 provenance=provenance, 

565 ) 

566 

567 def to_legacy_cell_coadd( 

568 self, copy: bool | None = None, plane_map: Mapping[str, MaskPlane] | None = None 

569 ) -> LegacyMultipleCellCoadd: 

570 """Convert to a `lsst.cell_coadds.MultipleCellCoadd` instance. 

571 

572 Parameters 

573 ---------- 

574 copy 

575 If `True`, always copy the image and variance pixel data. 

576 If `False`, return a view, and raise `TypeError` if the pixel data 

577 is read-only (this is not supported by afw). If `None`, only copy 

578 if the pixel data is read-only. Mask pixel data is always copied. 

579 plane_map 

580 A mapping from legacy mask plane name to the new plane name and 

581 description. 

582 """ 

583 from frozendict import frozendict 

584 

585 from lsst.cell_coadds import CellIdentifiers as LegacyCellIdentifiers 

586 from lsst.cell_coadds import CoaddUnits as LegacyCoaddUnits 

587 from lsst.cell_coadds import CommonComponents as LegacyCommonComponents 

588 from lsst.cell_coadds import MultipleCellCoadd as LegacyMultipleCellCoadd 

589 from lsst.cell_coadds import OwnedImagePlanes as LegacyOwnedImagePlanes 

590 from lsst.cell_coadds import PatchIdentifiers as LegacyPatchIdentifiers 

591 from lsst.cell_coadds import SingleCellCoadd as LegacySingleCellCoadd 

592 from lsst.skymap import Index2D as LegacyIndex2D 

593 

594 if plane_map is None: 

595 plane_map = get_legacy_deep_coadd_mask_planes() 

596 if self.unit != astropy.units.nJy: 

597 raise ValueError("CellCoadd.to_legacy requires nJy pixel units.") 

598 if self.bbox != self.bounds.bbox: 

599 raise ValueError("MultipleCellCoadd requires its bounding box to lie on the cell grid.") 

600 legacy_grid = self.grid.to_legacy() 

601 visit_polygons = self.provenance.to_legacy_polygon_map() 

602 legacy_common = LegacyCommonComponents( 

603 units=LegacyCoaddUnits.nJy, 

604 wcs=self.sky_projection.to_legacy(), 

605 band=self.band, 

606 identifiers=LegacyPatchIdentifiers( 

607 self.skymap, 

608 self.tract, 

609 LegacyIndex2D(x=self.patch.index.x, y=self.patch.index.y), 

610 band=self.band, 

611 ), 

612 visit_polygons=visit_polygons, 

613 ) 

614 legacy_inputs = self.provenance.to_legacy_cell_coadd_inputs(visit_polygons.keys()) 

615 cells: list[LegacySingleCellCoadd] = [] 

616 for cell_index in self.bounds.cell_indices(): 

617 cell_bbox = self.grid.bbox_of(cell_index) 

618 # Legacy type only has room for one mask_fractions plane. 

619 legacy_mask_fractions = ( 

620 next(iter(self.mask_fractions.values()))[cell_bbox].to_legacy(copy=copy) 

621 if self.mask_fractions 

622 else None 

623 ) 

624 legacy_planes = LegacyOwnedImagePlanes( 

625 image=self.image[cell_bbox].to_legacy(copy=copy), 

626 mask=self.mask[cell_bbox].to_legacy(plane_map), 

627 variance=self.variance[cell_bbox].to_legacy(copy=copy), 

628 mask_fractions=legacy_mask_fractions, 

629 noise_realizations=[n[cell_bbox].to_legacy(copy=copy) for n in self.noise_realizations], 

630 ) 

631 legacy_aperture_correction_map = frozendict( 

632 { 

633 name: field.value_in_cell(cell_index) 

634 if cell_index not in field.bounds.missing 

635 else np.nan 

636 for name, field in self.aperture_corrections.items() 

637 } 

638 ) 

639 cells.append( 

640 LegacySingleCellCoadd( 

641 legacy_planes, 

642 psf=self.psf[cell_index].to_legacy(copy=copy), 

643 inner_bbox=cell_bbox.to_legacy(), 

644 common=legacy_common, 

645 inputs=legacy_inputs[cell_index.to_legacy()], 

646 identifiers=LegacyCellIdentifiers( 

647 self.skymap, 

648 self.tract, 

649 legacy_common.identifiers.patch, 

650 band=self.band, 

651 cell=cell_index.to_legacy(), 

652 ), 

653 aperture_correction_map=legacy_aperture_correction_map, 

654 ) 

655 ) 

656 return LegacyMultipleCellCoadd( 

657 cells, 

658 legacy_grid, 

659 outer_cell_size=self.grid.cell_shape.to_legacy_int_extent(), 

660 psf_image_size=self.psf.kernel_bbox.shape.to_legacy_int_extent(), 

661 common=legacy_common, 

662 inner_bbox=self.bbox.to_legacy(), 

663 ) 

664 

665 @staticmethod 

666 def from_legacy_lazy_coadd( 

667 legacy: LegacyExposure, 

668 *, 

669 instrument: str, 

670 plane_map: Mapping[str, MaskPlane] | None = None, 

671 tract_info: TractInfo, 

672 skymap: str | None = None, 

673 patch: int | None = None, 

674 band: str | None = None, 

675 ) -> CellCoadd: 

676 """Convert from a `lsst.afw.image.Exposure` instance with a lazy 

677 PSF and/or aperture corrections. 

678 

679 This approximates the given image's PSF and aperture corrections by 

680 evaluating them at cell centers and assuming they are valid for the 

681 entirety of the cell. 

682 

683 Parameters 

684 ---------- 

685 legacy 

686 A `lsst.afw.image.Exposure` instance to convert. 

687 instrument 

688 Name of the instrument that produced all inputs to this coadd. 

689 plane_map 

690 A mapping from legacy mask plane name to the new plane name and 

691 description. 

692 tract_info 

693 Information about the full tract. 

694 skymap 

695 Name of the skymap. 

696 patch 

697 ID of the patch this coadd belongs to. 

698 band 

699 Name of the bandpass filter. 

700 """ 

701 if plane_map is None: 

702 plane_map = get_legacy_non_cell_coadd_mask_planes() 

703 legacy_md = legacy.getMetadata() 

704 if skymap is None: 

705 skymap = legacy_md.get("LSST BUTLER DATAID SKYMAP") 

706 if skymap is None: 

707 raise ValueError("No skymap name provided and no skymap name found in header.") 

708 if patch is None: 

709 patch = legacy_md.get("LSST BUTLER DATAID PATCH") 

710 if patch is None: 

711 raise ValueError("No patch ID provided and no patch ID found in header.") 

712 if band is None: 

713 band = legacy_md.get("LSST BUTLER DATAID BAND") 

714 if band is None: 

715 raise ValueError("No band name provided and no band name found in header.") 

716 tract_bbox = Box.from_legacy(tract_info.getBBox()) 

717 sky_projection = SkyProjection.from_legacy( 

718 legacy.wcs, 

719 TractFrame( 

720 skymap=skymap, 

721 tract=tract_info.getId(), 

722 bbox=tract_bbox, 

723 ), 

724 pixel_bounds=tract_bbox, 

725 ) 

726 image = Image.from_legacy(legacy.image, unit=astropy.units.nJy) 

727 if (legacy_photo_calib := legacy.getPhotoCalib()) is not None and ( 

728 field_from_legacy_photo_calib(legacy_photo_calib, image.bbox) is not None 

729 ): 

730 raise ValueError( 

731 "Legacy coadd pixels must already be calibrated to nJy (an identity PhotoCalib)." 

732 ) 

733 mask = Mask.from_legacy(legacy.mask, plane_map=plane_map) 

734 variance = Image.from_legacy(legacy.variance, unit=astropy.units.nJy**2) 

735 patch_definition = PatchDefinition.from_legacy_patch_info(tract_info[patch]) 

736 psf_bounds = CellGridBounds(grid=patch_definition.cells, bbox=patch_definition.cells.bbox)[ 

737 Box.from_legacy(legacy.getBBox()) 

738 ] 

739 psf = CellPointSpreadFunction.approximate_from_legacy(legacy.psf, psf_bounds) 

740 aperture_corrections = approximate_aperture_corrections_from_legacy( 

741 legacy.info.getApCorrMap(), psf.bounds 

742 ) 

743 provenance = CoaddProvenance.from_legacy_lazy( 

744 legacy.info.getCoaddInputs(), 

745 bounds=psf.bounds, 

746 sky_projection=sky_projection, 

747 instrument=instrument, 

748 ) 

749 return CellCoadd( 

750 image=image, 

751 mask=mask, 

752 variance=variance, 

753 mask_fractions={}, 

754 noise_realizations=[], 

755 sky_projection=sky_projection, 

756 band=band, 

757 psf=psf, 

758 aperture_corrections=aperture_corrections, 

759 patch=patch_definition, 

760 provenance=provenance, 

761 ) 

762 

763 def to_legacy( 

764 self, copy: bool | None = None, plane_map: Mapping[str, MaskPlane] | None = None 

765 ) -> LegacyExposure: 

766 """Convert to a `lsst.afw.image.Exposure` instance. 

767 

768 Parameters 

769 ---------- 

770 copy 

771 If `True`, always copy the image and variance pixel data. 

772 If `False`, return a view, and raise `TypeError` if the pixel data 

773 is read-only (this is not supported by afw). If `None`, only copy 

774 if the pixel data is read-only. Mask pixel data is always copied. 

775 plane_map 

776 A mapping from legacy mask plane name to the new plane name and 

777 description. 

778 

779 Returns 

780 ------- 

781 `lsst.afw.image.Exposure` 

782 A legacy representation of the coadd. This will have its ``wcs``, 

783 ``psf``, ``filter``, ``photoCalib``, and ``metadata`` components 

784 set. The ``apCorrMap`` component is not set, because there is no 

785 true `lsst.afw.math.BoundedField` representation for cell-coadd 

786 aperture corrections, and the ``coaddInputs`` component is not set 

787 because that data structure cannot fully capture cell-coadd 

788 provenance. 

789 

790 Notes 

791 ----- 

792 This method requires the `provenance` attribute to have been populated 

793 at construction. 

794 """ 

795 from lsst.afw.image import Exposure as LegacyExposure 

796 from lsst.afw.image import FilterLabel as LegacyFilterLabel 

797 

798 if plane_map is None: 

799 plane_map = get_legacy_deep_coadd_mask_planes() 

800 legacy_masked_image = super().to_legacy(copy=copy, plane_map=plane_map) 

801 result = LegacyExposure(legacy_masked_image, dtype=self.image.array.dtype) 

802 result_info = result.info 

803 result_info.setWcs(self.sky_projection.to_legacy()) 

804 result_info.setPsf(self.psf.to_legacy()) 

805 result_info.setFilter(LegacyFilterLabel.fromBand(self.band))