Coverage for python/lsst/obs/base/formatters/fitsExposure.py: 15%

368 statements  

« prev     ^ index     » next       coverage.py v7.14.1, created at 2026-05-29 08:25 +0000

1# This file is part of obs_base. 

2# 

3# Developed for the LSST Data Management System. 

4# This product includes software developed by the LSST Project 

5# (http://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 <http://www.gnu.org/licenses/>. 

21 

22__all__ = ( 

23 "FitsExposureFormatter", 

24 "FitsImageFormatter", 

25 "FitsImageFormatterBase", 

26 "FitsMaskFormatter", 

27 "FitsMaskedImageFormatter", 

28 "StandardFitsImageFormatterBase", 

29 "standardizeAmplifierParameters", 

30) 

31 

32import hashlib 

33import json 

34import logging 

35import threading 

36import uuid 

37import warnings 

38from abc import abstractmethod 

39from collections.abc import Mapping, Set 

40from io import BytesIO 

41from typing import TYPE_CHECKING, Any, ClassVar, NamedTuple, Protocol 

42 

43import astropy.io.fits 

44import numpy as np 

45 

46import lsst.geom 

47from lsst.afw.cameraGeom import AmplifierGeometryComparison, AmplifierIsolator 

48from lsst.afw.fits import CompressionOptions, MemFileManager 

49from lsst.afw.geom.wcsUtils import getImageXY0FromMetadata 

50from lsst.afw.image import ( 

51 ExposureFitsReader, 

52 ExposureInfo, 

53 FilterLabel, 

54 ImageFitsReader, 

55 MaskedImageFitsReader, 

56 MaskFitsReader, 

57) 

58 

59# Needed for ApCorrMap to resolve properly 

60from lsst.afw.math import BoundedField # noqa: F401 

61from lsst.daf.base import PropertyList 

62from lsst.daf.butler import DatasetProvenance, FormatterV2 

63from lsst.resources import ResourcePath 

64from lsst.utils.classes import cached_getter 

65from lsst.utils.introspection import find_outside_stacklevel 

66 

67from ..utils import add_provenance_to_fits_header 

68 

69if TYPE_CHECKING: 

70 import lsst.afw.cameraGeom 

71 

72 

73_LOG = logging.getLogger(__name__) 

74 

75_ALWAYS_USE_ASTROPY_FOR_COMPONENT_READ = False 

76"""If True, the astropy code will always be used to read component and cutouts 

77even if the file is local, the cutout is too large, or the dataset type is 

78wrong. This should mostly be used for testing. 

79""" 

80 

81 

82class _ReaderClassLike(Protocol): 

83 def __init__(self, path: str) -> None: ... 

84 def readBBox(self) -> lsst.geom.Box2I: ... 

85 def read(self, bbox: lsst.geom.Box2I = lsst.geom.Box2I(), dtype: Any = None) -> Any: ... 

86 def readImage(self, bbox: lsst.geom.Box2I = lsst.geom.Box2I(), dtype: Any = None) -> Any: ... 

87 def readMask(self, bbox: lsst.geom.Box2I = lsst.geom.Box2I(), dtype: Any = None) -> Any: ... 

88 def readVariance(self, bbox: lsst.geom.Box2I = lsst.geom.Box2I(), dtype: Any = None) -> Any: ... 

89 def readDetector(self) -> lsst.afw.cameraGeom.Detector: ... 

90 def readComponent(self, component: str) -> Any: ... 

91 def readMetadata(self) -> PropertyList: ... 

92 def readSerializationVersion(self) -> int: ... 

93 

94 

95class FitsImageFormatterBase(FormatterV2): 

96 """Base class formatter for image-like storage classes stored via FITS. 

97 

98 Notes 

99 ----- 

100 This class makes no assumptions about how many HDUs are used to represent 

101 the image on disk, and includes no support for writing. It's really just a 

102 collection of miscellaneous boilerplate common to all FITS image 

103 formatters. 

104 

105 Concrete subclasses must implement `readComponent`, `readFull`, and 

106 `write_local_file` (even if just to disable them by raising an exception). 

107 """ 

108 

109 can_read_from_local_file = True 

110 default_extension = ".fits" 

111 supported_extensions: ClassVar[Set[str]] = frozenset({".fits", ".fits.gz", ".fits.fz", ".fz", ".fit"}) 

112 

113 unsupported_parameters: ClassVar[Set[str]] = frozenset() 

114 """Support all parameters.""" 

115 

116 _reader = None 

117 _reader_path: str | None = None 

118 

119 ReaderClass: type[_ReaderClassLike] # must be set by concrete subclasses 

120 """Class to use for reading FITS files in the expected way. 

121 (e.g., `type` [`lsst.afw.image.ImageFitsReader]) 

122 """ 

123 

124 @property 

125 def reader(self) -> _ReaderClassLike: 

126 """The reader object that backs this formatter's read operations. 

127 

128 This is computed on first use and then cached. It should never be 

129 accessed when writing. Currently assumes a local file. 

130 """ 

131 if self._reader is None: 

132 if self._reader_path is None: 

133 raise RuntimeError("Internal error in formatter; failing to set path.") 

134 self._reader = self.ReaderClass(self._reader_path) 

135 return self._reader 

136 

137 @property 

138 @cached_getter 

139 def checked_parameters(self) -> dict[str, Any]: 

140 """The parameters passed by the butler user, after checking them 

141 against the storage class and transforming `None` into an empty `dict` 

142 (`dict`). 

143 

144 This is computed on first use and then cached. It should never be 

145 accessed when writing. Subclasses that need additional checking should 

146 delegate to `super` and then check the result before returning it. 

147 """ 

148 parameters = self.file_descriptor.parameters 

149 if parameters is None: 

150 parameters = {} 

151 self.file_descriptor.storageClass.validateParameters(parameters) 

152 return parameters 

153 

154 @property 

155 def storageClass_dtype(self) -> np.dtype | None: 

156 """The numpy data type associated with the storage class.""" 

157 dtype: np.dtype | None = None 

158 try: 

159 # lsst.afw.image.Exposure is generic base class and does not have 

160 # the dtype attribute. 

161 dtype = np.dtype(self.file_descriptor.storageClass.pytype.dtype) # type: ignore[attr-defined] 

162 except AttributeError: 

163 pass 

164 return dtype 

165 

166 def read_from_local_file(self, path: str, component: str | None = None, expected_size: int = -1) -> Any: 

167 # Docstring inherited. 

168 if _is_future_visit_image(self.file_descriptor.readStorageClass.name, component): 

169 from lsst.images import VisitImage 

170 

171 return VisitImage.read_legacy( 

172 path, 

173 component=component, 

174 preserve_quantization=self.checked_parameters.get("preserve_quantization", False), 

175 ) 

176 elif self.checked_parameters.get("preserve_quantization", False): 

177 raise NotImplementedError( 

178 "preserve_quantization=True only works when converting to VisitImage on read." 

179 ) 

180 

181 # The methods doing the reading all currently assume local file 

182 # and assume that the file descriptor refers to a local file. 

183 # With FormatterV2 that file descriptor does not refer to a local 

184 # file. 

185 self._reader_path = path 

186 self._reader = None # Ensure the reader class is reset. 

187 try: 

188 if component is not None: 

189 in_memory_dataset = self.readComponent(component) 

190 else: 

191 in_memory_dataset = self.readFull() 

192 finally: 

193 self._reader = None # Release the file handle. 

194 return in_memory_dataset 

195 

196 @abstractmethod 

197 def readComponent(self, component: str) -> Any: 

198 """Read a component dataset. 

199 

200 Parameters 

201 ---------- 

202 component : `str`, optional 

203 Component to read from the file. 

204 

205 Returns 

206 ------- 

207 obj : `typing.Any` 

208 In-memory component object. 

209 

210 Raises 

211 ------ 

212 KeyError 

213 Raised if the requested component cannot be handled. 

214 """ 

215 raise NotImplementedError() 

216 

217 @abstractmethod 

218 def readFull(self) -> Any: 

219 """Read the full dataset (while still accounting for parameters). 

220 

221 Returns 

222 ------- 

223 obj : `typing.Any` 

224 In-memory component object. 

225 

226 """ 

227 raise NotImplementedError() 

228 

229 

230class StandardFitsImageFormatterBase(FitsImageFormatterBase): 

231 """Base class interface for image-like storage stored via FITS, 

232 written using LSST code. 

233 

234 Notes 

235 ----- 

236 Concrete subclasses must provide at least the ``ReaderClass`` attribute. 

237 

238 The provided implementation of `readComponent` handles only the 'bbox', 

239 'dimensions', and 'xy0' components common to all image-like storage 

240 classes. Subclasses with additional components should handle them first, 

241 then delegate to ``super()`` for these (or, if necessary, delegate first 

242 and catch `KeyError`). 

243 

244 The provided implementation of `readFull` handles only parameters that 

245 can be forwarded directly to the reader class (usually ``bbox`` and 

246 ``origin``). Concrete subclasses that need to handle additional parameters 

247 should generally reimplement without delegating (the implementation is 

248 trivial). 

249 

250 This Formatter supports write recipes, and assumes its in-memory type has 

251 ``writeFits`` and (for write recipes) ``writeFitsWithOptions`` methods. 

252 

253 Each ``StandardFitsImageFormatterBase`` recipe for FITS compression should 

254 define ``image``, ``mask`` and ``variance`` entries, each of which may 

255 contain entries supported by 

256 `lsst.afw.fits.CompressionOptions.from_mapping` (``null`` disables 

257 compression). 

258 

259 A very simple example YAML recipe (for the `lsst.afw.image.Exposure` 

260 specialization): 

261 

262 .. code-block:: yaml 

263 

264 lsst.obs.base.fitsExposureFormatter.FitsExposureFormatter: 

265 default: 

266 image: &default 

267 algorithm: GZIP_2 

268 mask: *default 

269 variance: *default 

270 

271 """ 

272 

273 supported_write_parameters = frozenset({"recipe"}) 

274 

275 def readComponent(self, component: str) -> Any: 

276 # Docstring inherited. 

277 if component in ("bbox", "dimensions", "xy0"): 

278 bbox = self.reader.readBBox() 

279 if component == "dimensions": 

280 return bbox.getDimensions() 

281 elif component == "xy0": 

282 return bbox.getMin() 

283 else: 

284 return bbox 

285 else: 

286 raise KeyError(f"Unknown component requested: {component}") 

287 

288 def readFull(self) -> Any: 

289 # Docstring inherited. 

290 return self.reader.read(**self.checked_parameters, dtype=self.storageClass_dtype) 

291 

292 def write_local_file(self, in_memory_dataset: Any, uri: ResourcePath) -> None: 

293 """Serialize the image to FITS. 

294 

295 Parameters 

296 ---------- 

297 in_memory_dataset : `object` 

298 Image to write. Must support a ``writeFits`` or 

299 ``writeFitsWithOptions`` interface. 

300 uri : `lsst.resources.ResourcePath` 

301 Location to write the local file. 

302 """ 

303 # check to see if we have a recipe requested 

304 recipeName = self.write_parameters.get("recipe") 

305 recipe = self.get_image_compression_settings(recipeName) 

306 if recipe: 

307 in_memory_dataset.writeFitsWithOptions(uri.ospath, options=recipe) 

308 else: 

309 in_memory_dataset.writeFits(uri.ospath) 

310 

311 def get_image_compression_settings(self, recipeName: str | None) -> dict: 

312 """Retrieve the relevant compression settings for this recipe. 

313 

314 Parameters 

315 ---------- 

316 recipeName : `str` or `None` 

317 Label associated with the collection of compression parameters 

318 to select. 

319 

320 Returns 

321 ------- 

322 settings : `dict` 

323 The selected settings. 

324 """ 

325 # if no recipe has been provided and there is no default 

326 # return immediately 

327 if not recipeName: 

328 if "default" not in self.write_recipes: 

329 return {} 

330 recipeName = "default" 

331 

332 if recipeName not in self.write_recipes: 

333 raise RuntimeError(f"Unrecognized recipe option given for compression: {recipeName}") 

334 

335 recipe = self.write_recipes[recipeName] 

336 if recipe is None: 

337 return {} 

338 seed: int | None = None 

339 for plane in ("image", "mask", "variance"): 

340 if plane in recipe and (quantization := recipe[plane].get("quantization")) is not None: 

341 if quantization.get("seed", 0) == 0: 

342 if seed is None: 

343 # Set the seed based on data ID. We can't just use 

344 # 'hash', since like 'set' that's not deterministic. 

345 # And we can't rely on a DimensionPacker because those 

346 # are only defined for certain combinations of 

347 # dimensions. Doing an MD5 of the JSON feels like 

348 # overkill but I don't really see anything much 

349 # simpler. 

350 hash_bytes = hashlib.md5( 

351 json.dumps(list(self.data_id.required_values)).encode(), 

352 usedforsecurity=False, 

353 ).digest() 

354 # And it *really* feels like overkill when we squash 

355 # that into the [1, 10000] range allowed by FITS. 

356 seed = 1 + int.from_bytes(hash_bytes) % 9999 

357 _LOG.debug( 

358 "Setting compression quantization seed for %s %s %s to %s.", 

359 self.data_id, 

360 self.dataset_ref.datasetType.name, 

361 plane, 

362 seed, 

363 ) 

364 quantization["seed"] = seed 

365 else: 

366 _LOG.warning( 

367 "Compression quantization seed for %s %s %s was set explicitly to %s.", 

368 self.dataset_ref.datasetType.name, 

369 self.data_id, 

370 plane, 

371 quantization["seed"], 

372 ) 

373 else: 

374 _LOG.debug( 

375 "No quantization found for %s %s %s.", 

376 self.dataset_ref.datasetType.name, 

377 self.data_id, 

378 plane, 

379 ) 

380 return recipe 

381 

382 @classmethod 

383 def validate_write_recipes(cls, recipes: Mapping[str, Any] | None) -> Mapping[str, Any] | None: 

384 """Validate supplied recipes for this formatter. 

385 

386 The recipes are supplemented with default values where appropriate. 

387 

388 Parameters 

389 ---------- 

390 recipes : `dict` or `None` 

391 Recipes to validate. Can be empty dict or `None`. 

392 

393 Returns 

394 ------- 

395 validated : `dict` 

396 Validated recipes. Returns what was given if there are no 

397 recipes listed. 

398 

399 Raises 

400 ------ 

401 RuntimeError 

402 Raised if validation fails. 

403 """ 

404 if not recipes: 

405 # We can not insist on recipes being specified. 

406 return recipes 

407 

408 validated: dict[str, Any] = {} 

409 for name, recipe in recipes.items(): 

410 if recipe is not None: 

411 validated[name] = {} 

412 for plane in ["image", "mask", "variance"]: 

413 try: 

414 options = CompressionOptions.from_mapping(recipe[plane]) 

415 except Exception as err: 

416 err.add_note(f"Validating write recipe {name!r} ({plane!r} section).") 

417 raise 

418 validated[name][plane] = options.to_dict() 

419 else: 

420 validated[name] = None 

421 return validated 

422 

423 

424class FitsImageFormatter(StandardFitsImageFormatterBase): 

425 """Concrete formatter for reading/writing `~lsst.afw.image.Image` 

426 from/to FITS. 

427 """ 

428 

429 ReaderClass = ImageFitsReader 

430 

431 

432class FitsMaskFormatter(StandardFitsImageFormatterBase): 

433 """Concrete formatter for reading/writing `~lsst.afw.image.Mask` 

434 from/to FITS. 

435 """ 

436 

437 ReaderClass = MaskFitsReader 

438 

439 

440class FitsMaskedImageFormatter(StandardFitsImageFormatterBase): 

441 """Concrete formatter for reading/writing `~lsst.afw.image.MaskedImage` 

442 from/to FITS. 

443 """ 

444 

445 ReaderClass = MaskedImageFitsReader 

446 

447 def readComponent(self, component: str) -> Any: 

448 # Docstring inherited. 

449 if component == "image": 

450 return self.reader.readImage(**self.checked_parameters, dtype=self.storageClass_dtype) 

451 elif component == "mask": 

452 return self.reader.readMask(**self.checked_parameters) 

453 elif component == "variance": 

454 return self.reader.readVariance(**self.checked_parameters, dtype=self.storageClass_dtype) 

455 else: 

456 # Delegate to base for bbox, dimensions, xy0. 

457 return super().readComponent(component) 

458 

459 

460def standardizeAmplifierParameters( 

461 parameters: dict[str, Any], on_disk_detector: lsst.afw.cameraGeom.Detector | None 

462) -> tuple[lsst.afw.cameraGeom.Amplifier, lsst.afw.cameraGeom.Detector, bool]: 

463 """Preprocess the Exposure storage class's "amp" and "detector" parameters. 

464 

465 This checks the given objects for consistency with the on-disk geometry and 

466 converts amplifier IDs/names to Amplifier instances. 

467 

468 Parameters 

469 ---------- 

470 parameters : `dict` 

471 Dictionary of parameters passed to formatter. See the Exposure storage 

472 class definition in daf_butler for allowed keys and values. 

473 on_disk_detector : `lsst.afw.cameraGeom.Detector` or `None` 

474 Detector that represents the on-disk image being loaded, or `None` if 

475 this is unknown (and hence the user must provide one in 

476 ``parameters`` if "amp" is in ``parameters``). 

477 

478 Returns 

479 ------- 

480 amplifier : `lsst.afw.cameraGeom.Amplifier` or `None` 

481 An amplifier object that defines a subimage to load, or `None` if there 

482 was no "amp" parameter. 

483 detector : `lsst.afw.cameraGeom.Detector` or `None` 

484 A detector object whose amplifiers are in the same s/orientation 

485 state as the on-disk image. If there is no "amp" parameter, 

486 ``on_disk_detector`` is simply passed through. 

487 regions_differ : `bool` 

488 `True` if the on-disk detector and the detector given in the parameters 

489 had different bounding boxes for one or more regions. This can happen 

490 if the true overscan region sizes can only be determined when the image 

491 is actually read, but otherwise it should be considered user error. 

492 """ 

493 if (amplifier := parameters.get("amp")) is None: 

494 return None, on_disk_detector, False 

495 if "bbox" in parameters or "origin" in parameters: 

496 raise ValueError("Cannot pass 'amp' with 'bbox' or 'origin'.") 

497 if isinstance(amplifier, int | str): 

498 amp_key = amplifier 

499 target_amplifier = None 

500 else: 

501 amp_key = amplifier.getName() 

502 target_amplifier = amplifier 

503 if (detector := parameters.get("detector")) is not None: 

504 if on_disk_detector is not None: 

505 # User passed a detector and we also found one on disk. Check them 

506 # for consistency. Note that we are checking the amps we'd get 

507 # from the two detectors against each other, not the amplifier we 

508 # got directly from the user, as the latter is allowed to differ in 

509 # assembly/orientation state. 

510 comparison = on_disk_detector[amp_key].compareGeometry(detector[amp_key]) 

511 if comparison & comparison.ASSEMBLY_DIFFERS: 

512 raise ValueError( 

513 "The given 'detector' has a different assembly state and/or orientation from " 

514 f"the on-disk one for amp {amp_key}." 

515 ) 

516 else: 

517 if on_disk_detector is None: 

518 raise ValueError( 

519 f"No on-disk detector and no detector given; cannot load amplifier from key {amp_key}. " 

520 "Please provide either a 'detector' parameter or an Amplifier instance in the " 

521 "'amp' parameter." 

522 ) 

523 comparison = AmplifierGeometryComparison.EQUAL 

524 detector = on_disk_detector 

525 if target_amplifier is None: 

526 target_amplifier = detector[amp_key] 

527 return target_amplifier, detector, comparison & comparison.REGIONS_DIFFER 

528 

529 

530class _ComponentCache(NamedTuple): 

531 id_: uuid.UUID | None = None 

532 reader: ExposureFitsReader | None = None 

533 bbox: lsst.geom.Box2I | None = None 

534 mem: MemFileManager | None = None 

535 

536 

537class FitsExposureFormatter(FitsMaskedImageFormatter): 

538 """Concrete formatter for reading/writing `~lsst.afw.image.Exposure` 

539 from/to FITS. 

540 

541 Notes 

542 ----- 

543 This class inherits from `FitsMaskedImageFormatter` even though 

544 `lsst.afw.image.Exposure` doesn't inherit from 

545 `lsst.afw.image.MaskedImage`; this is just an easy way to be able to 

546 delegate to `FitsMaskedImageFormatter` for component-handling, and 

547 should be replaced with e.g. both calling a free function if that slight 

548 type covariance violation ever becomes a practical problem. 

549 """ 

550 

551 can_read_from_uri = True 

552 ReaderClass = ExposureFitsReader 

553 # TODO: Remove MemFileManager from cache when DM-49640 is fixed. 

554 _lock = threading.Lock() 

555 _cached_fits: _ComponentCache = _ComponentCache() 

556 

557 def read_from_uri(self, uri: ResourcePath, component: str | None = None, expected_size: int = -1) -> Any: 

558 # For now only support small non-pixel components. In future 

559 # could work with cutouts. 

560 self._reader = None # Guarantee things are reset. 

561 

562 parameters = self.checked_parameters.copy() 

563 preserve_quantization = parameters.pop("preserve_quantization", False) 

564 # Full read, always use local file read. 

565 if preserve_quantization or (not component and not parameters): 

566 return NotImplemented 

567 

568 if not _ALWAYS_USE_ASTROPY_FOR_COMPONENT_READ and uri.isLocal: 

569 # For a local URI allow afw to read it directly. 

570 return NotImplemented 

571 pixel_components = ("mask", "image", "variance") 

572 

573 if component in pixel_components: 

574 # For pixel access currently this can not be cached in memory 

575 # and the performance gains are unclear. Assume local file 

576 # read with file caching for now. 

577 return NotImplemented 

578 

579 # With current file layouts the non-pixel extensions account for 1/3 

580 # of the file size and it is more efficient to download the entire 

581 # file. 

582 if not ( 

583 _ALWAYS_USE_ASTROPY_FOR_COMPONENT_READ 

584 or self._dataset_ref.dataId.mapping.keys().isdisjoint({"tract", "patch"}) 

585 ): 

586 return NotImplemented 

587 

588 # Cutouts can be optimized. For now only use this optimization 

589 # if bbox is the only parameter and the number of pixels in the 

590 # bounding box is reasonable. 

591 bbox = None 

592 origin = lsst.afw.image.PARENT 

593 if not component: 

594 # Try to support PARENT and LOCAL origin but if there are any 

595 # other parameters do not attempt a cutout. 

596 if parameters.keys() - {"bbox", "origin"}: 

597 return NotImplemented 

598 bbox = parameters["bbox"] 

599 origin = parameters.get("origin", lsst.afw.image.PARENT) 

600 # For larger cutouts use the full file. 

601 max_cutout_size = 500 * 500 

602 if not _ALWAYS_USE_ASTROPY_FOR_COMPONENT_READ and bbox.width * bbox.height > max_cutout_size: 

603 return NotImplemented 

604 

605 # We only cache component reads since those are small. 

606 if component: 

607 with self._lock: 

608 cache = type(self)._cached_fits 

609 if self.dataset_ref.id == cache.id_: 

610 if component in {"xy0", "dimensions", "bbox"} and cache.bbox is not None: 

611 match component: 

612 case "xy0": 

613 return cache.bbox.getMin() 

614 case "dimensions": 

615 return cache.bbox.getDimensions() 

616 case "bbox": 

617 return cache.bbox 

618 else: 

619 self._reader = cache.reader 

620 return self.readComponent(component) 

621 

622 try: 

623 fs, fspath = uri.to_fsspec() 

624 except Exception: 

625 # fsspec cannot be initialized, fall back to downloading the file. 

626 return NotImplemented 

627 

628 bbox_component = None 

629 try: 

630 hdul = [] 

631 with fs.open(fspath) as f, astropy.io.fits.open(f) as fits_obj: 

632 # Read all non-pixel components and cache. 

633 for hdu in fits_obj: 

634 hdr = hdu.header 

635 extname = hdr.get("EXTNAME") 

636 # Older files have IMAGE in EXTNAME in PRIMARY so check 

637 # for EXTEND=T. 

638 extend = hdr.get("EXTEND") 

639 if not extend and extname and extname.lower() in pixel_components: 

640 # Calculate the dimensional components for later 

641 # caching. Do not derive from cached FITS reader 

642 # because they depend on the dimensionality of the 

643 # pixel data and we do not want to cache the pixel 

644 # data. 

645 if bbox_component is None: 

646 shape = hdu.shape 

647 dimensions = lsst.geom.Extent2I(shape[1], shape[0]) 

648 

649 # XY0 is defined in the A WCS. 

650 pl = PropertyList() 

651 pl.update(hdr) 

652 xy0 = getImageXY0FromMetadata(pl, "A", strip=False) 

653 

654 # This is the PARENT bbox. 

655 bbox_component = lsst.geom.Box2I(xy0, dimensions) 

656 

657 # Handle cutout request. 

658 if bbox: 

659 if origin == lsst.afw.image.PARENT: 

660 full_bbox = bbox_component 

661 else: 

662 full_bbox = lsst.geom.Box2I( 

663 lsst.geom.Point2I(0, 0), bbox_component.getDimensions 

664 ) 

665 minX = bbox.getBeginX() - full_bbox.getBeginX() 

666 maxX = bbox.getEndX() - full_bbox.getBeginX() 

667 minY = bbox.getBeginY() - full_bbox.getBeginY() 

668 maxY = bbox.getEndY() - full_bbox.getBeginY() 

669 data = hdu.section[minY:maxY, minX:maxX] 

670 

671 # Must correct the header WCS to take into 

672 # account the offset. 

673 if (k := "CRPIX1") in hdr: 

674 hdr[k] -= minX 

675 if (k := "CRPIX2") in hdr: 

676 hdr[k] -= minY 

677 if (k := "LTV1") in hdr: 

678 hdr[k] = -bbox.getBeginX() 

679 if (k := "LTV2") in hdr: 

680 hdr[k] = -bbox.getBeginY() 

681 if (k := "CRVAL1A") in hdr: 

682 hdr[k] = bbox.getBeginX() 

683 if (k := "CRVAL2A") in hdr: 

684 hdr[k] = bbox.getBeginY() 

685 else: 

686 data = np.zeros([1, 1], dtype=np.int32) 

687 

688 # Construct a new HDU and copy the header. 

689 stripped_hdu = astropy.io.fits.ImageHDU(data=data, header=hdr) 

690 hdul.append(stripped_hdu) 

691 else: 

692 hdul.append(hdu) 

693 stripped_fits = astropy.io.fits.HDUList(hdus=hdul) 

694 # Write the FITS file to in-memory FITS. 

695 buffer = BytesIO() 

696 stripped_fits.writeto(buffer) 

697 except Exception as e: 

698 # For some reason we can't open the remote file so fall back. 

699 _LOG.debug( 

700 "Attempted remote read of components but encountered an error. " 

701 "Falling back to file download. Error: %s", 

702 str(e), 

703 ) 

704 return NotImplemented 

705 

706 # Pass the new FITS buffer to the reader class without going through 

707 # a temporary file. We can assume this is relatively small for 

708 # components. 

709 fits_data = buffer.getvalue() 

710 mem = MemFileManager(len(fits_data)) 

711 mem.setData(fits_data, len(fits_data)) 

712 self._reader = self.ReaderClass(mem) 

713 

714 if component: 

715 with self._lock: 

716 type(self)._cached_fits = _ComponentCache( 

717 id_=self.dataset_ref.id, 

718 reader=self._reader, 

719 mem=mem, 

720 bbox=bbox_component, 

721 ) 

722 match component: 

723 case "xy0": 

724 if bbox_component is None: # For mypy. 

725 return None 

726 return bbox_component.getMin() 

727 case "dimensions": 

728 if bbox_component is None: 

729 return None 

730 return bbox_component.getDimensions() 

731 case "bbox": 

732 return bbox_component 

733 case _: 

734 return self.readComponent(component) 

735 else: 

736 # Must be a cutout. We have applied the bbox parameter so no 

737 # parameters should be passed here. 

738 cutout = self.reader.read(dtype=self.storageClass_dtype) 

739 cutout.getInfo().setFilter(self._fixFilterLabels(cutout.getInfo().getFilter())) 

740 return cutout 

741 

742 def add_provenance( 

743 self, in_memory_dataset: Any, /, *, provenance: DatasetProvenance | None = None 

744 ) -> Any: 

745 # Add provenance via FITS headers. 

746 add_provenance_to_fits_header(in_memory_dataset.metadata, self.dataset_ref, provenance) 

747 return in_memory_dataset 

748 

749 def readComponent(self, component: str) -> Any: 

750 # Docstring inherited. 

751 # Generic components can be read via a string name; DM-27754 will make 

752 # this mapping larger at the expense of the following one. 

753 genericComponents = { 

754 "summaryStats": ExposureInfo.KEY_SUMMARY_STATS, 

755 } 

756 if (genericComponentName := genericComponents.get(component)) is not None: 

757 return self.reader.readComponent(genericComponentName) 

758 # Other components have hard-coded method names, but don't take 

759 # parameters. 

760 standardComponents = { 

761 "id": "readExposureId", 

762 "metadata": "readMetadata", 

763 "wcs": "readWcs", 

764 "coaddInputs": "readCoaddInputs", 

765 "psf": "readPsf", 

766 "photoCalib": "readPhotoCalib", 

767 "filter": "readFilter", 

768 "validPolygon": "readValidPolygon", 

769 "apCorrMap": "readApCorrMap", 

770 "visitInfo": "readVisitInfo", 

771 "transmissionCurve": "readTransmissionCurve", 

772 "detector": "readDetector", 

773 "exposureInfo": "readExposureInfo", 

774 } 

775 if (methodName := standardComponents.get(component)) is not None: 

776 result = getattr(self.reader, methodName)() 

777 if component == "filter": 

778 return self._fixFilterLabels(result) 

779 return result 

780 # Delegate to MaskedImage and ImageBase implementations for the rest. 

781 return super().readComponent(component) 

782 

783 def readFull(self) -> Any: 

784 # Docstring inherited. 

785 amplifier, detector, _ = standardizeAmplifierParameters( 

786 self.checked_parameters, 

787 self.reader.readDetector(), 

788 ) 

789 if amplifier is not None: 

790 amplifier_isolator = AmplifierIsolator( 

791 amplifier, 

792 self.reader.readBBox(), 

793 detector, 

794 ) 

795 result = amplifier_isolator.transform_subimage( 

796 self.reader.read(bbox=amplifier_isolator.subimage_bbox, dtype=self.storageClass_dtype) 

797 ) 

798 result.setDetector(amplifier_isolator.make_detector()) 

799 else: 

800 result = self.reader.read(**self.checked_parameters, dtype=self.storageClass_dtype) 

801 result.getInfo().setFilter(self._fixFilterLabels(result.getInfo().getFilter())) 

802 return result 

803 

804 def _fixFilterLabels( 

805 self, file_filter_label: lsst.afw.image.FilterLabel, should_be_standardized: bool | None = None 

806 ) -> lsst.afw.image.FilterLabel: 

807 """Compare the filter label read from the file with the one in the 

808 data ID. 

809 

810 Parameters 

811 ---------- 

812 file_filter_label : `lsst.afw.image.FilterLabel` or `None` 

813 Filter label read from the file, if there was one. 

814 should_be_standardized : `bool`, optional 

815 If `True`, expect ``file_filter_label`` to be consistent with the 

816 data ID and warn only if it is not. If `False`, expect it to be 

817 inconsistent and warn only if the data ID is incomplete and hence 

818 the `FilterLabel` cannot be fixed. If `None` (default) guess 

819 whether the file should be standardized by looking at the 

820 serialization version number in file, which requires this method to 

821 have been run after `readFull` or `readComponent`. 

822 

823 Returns 

824 ------- 

825 filter_label : `lsst.afw.image.FilterLabel` or `None` 

826 The preferred filter label; may be the given one or one built from 

827 the data ID. `None` is returned if there should never be any 

828 filters associated with this dataset type. 

829 

830 Notes 

831 ----- 

832 Most test coverage for this method is in ci_hsc_gen3, where we have 

833 much easier access to test data that exhibits the problems it attempts 

834 to solve. 

835 """ 

836 # Remember filter data ID keys that weren't in this particular data ID, 

837 # so we can warn about them later. 

838 missing = [] 

839 band = None 

840 physical_filter = None 

841 if "band" in self.data_id.dimensions.names: 

842 band = self.data_id.get("band") 

843 # band isn't in the data ID; is that just because this data ID 

844 # hasn't been filled in with everything the Registry knows, or 

845 # because this dataset is never associated with a band? 

846 if band is None and not self.data_id.hasFull() and "band" in self.data_id.dimensions.implied: 

847 missing.append("band") 

848 if "physical_filter" in self.data_id.dimensions.names: 

849 physical_filter = self.data_id.get("physical_filter") 

850 # Same check as above for band, but for physical_filter. 

851 if ( 

852 physical_filter is None 

853 and not self.data_id.hasFull() 

854 and "physical_filter" in self.data_id.dimensions.implied 

855 ): 

856 missing.append("physical_filter") 

857 if should_be_standardized is None: 

858 version = self.reader.readSerializationVersion() 

859 should_be_standardized = version >= 2 

860 if missing: 

861 # Data ID identifies a filter but the actual filter label values 

862 # haven't been fetched from the database; we have no choice but 

863 # to use the one in the file. 

864 # Warn if that's more likely than not to be bad, because the file 

865 # predates filter standardization. 

866 if not should_be_standardized: 

867 warnings.warn( 

868 f"Data ID {self.data_id} is missing (implied) value(s) for {missing}; " 

869 "the correctness of this Exposure's FilterLabel cannot be guaranteed. " 

870 "Call Registry.expandDataId before Butler.get to avoid this.", 

871 # Report the warning from outside of middleware or the 

872 # relevant runQuantum method. 

873 stacklevel=find_outside_stacklevel( 

874 "lsst.obs.base", "lsst.pipe.base", "lsst.daf.butler", allow_methods={"runQuantum"} 

875 ), 

876 ) 

877 return file_filter_label 

878 if band is None and physical_filter is None: 

879 data_id_filter_label = None 

880 else: 

881 data_id_filter_label = FilterLabel(band=band, physical=physical_filter) 

882 if data_id_filter_label != file_filter_label and should_be_standardized: 

883 # File was written after FilterLabel and standardization, but its 

884 # FilterLabel doesn't agree with the data ID: this indicates a bug 

885 # in whatever code produced the Exposure (though it may be one that 

886 # has been fixed since the file was written). 

887 warnings.warn( 

888 f"Reading {self.file_descriptor.location} with data ID {self.data_id}: " 

889 f"filter label mismatch (file is {file_filter_label}, data ID is " 

890 f"{data_id_filter_label}). This is probably a bug in the code that produced it.", 

891 stacklevel=find_outside_stacklevel( 

892 "lsst.obs.base", "lsst.pipe.base", "lsst.daf.butler", allow_methods={"runQuantum"} 

893 ), 

894 ) 

895 return data_id_filter_label 

896 

897 

898def _is_future_visit_image(storage_class_name: str, component: str | None) -> bool: 

899 match storage_class_name, component: 

900 case ("VisitImage", None): 

901 return True 

902 case ("ImageV2", "image" | "variance"): 

903 return True 

904 case ("MaskV2", "mask"): 

905 return True 

906 case ("BoxV2", "bbox"): 

907 return True 

908 case ("PointSpreadFunction", "psf"): 

909 return True 

910 case ("DetectorV2", "detector"): 

911 return True 

912 # The components below can't be used unless we fix daf_butler 

913 # restrictions on component names (they're checked against the original 

914 # storage class component names). 

915 case ("ObservationSummaryStats", "summary_stats"): 

916 return True 

917 case ("ObservationInfo", "obs_info"): 

918 return True 

919 case ("StructuredDataDict", "aperture_corrections"): 

920 return True 

921 case ("ImageField", "photometric_scaling"): 

922 return True 

923 return False