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

338 statements  

« prev     ^ index     » next       coverage.py v7.13.2, created at 2026-02-13 01:17 -0800

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 # The methods doing the reading all currently assume local file 

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

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

171 # file. 

172 self._reader_path = path 

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

174 try: 

175 if component is not None: 

176 in_memory_dataset = self.readComponent(component) 

177 else: 

178 in_memory_dataset = self.readFull() 

179 finally: 

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

181 return in_memory_dataset 

182 

183 @abstractmethod 

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

185 """Read a component dataset. 

186 

187 Parameters 

188 ---------- 

189 component : `str`, optional 

190 Component to read from the file. 

191 

192 Returns 

193 ------- 

194 obj : `typing.Any` 

195 In-memory component object. 

196 

197 Raises 

198 ------ 

199 KeyError 

200 Raised if the requested component cannot be handled. 

201 """ 

202 raise NotImplementedError() 

203 

204 @abstractmethod 

205 def readFull(self) -> Any: 

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

207 

208 Returns 

209 ------- 

210 obj : `typing.Any` 

211 In-memory component object. 

212 

213 """ 

214 raise NotImplementedError() 

215 

216 

217class StandardFitsImageFormatterBase(FitsImageFormatterBase): 

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

219 written using LSST code. 

220 

221 Notes 

222 ----- 

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

224 

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

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

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

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

229 and catch `KeyError`). 

230 

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

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

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

234 should generally reimplement without delegating (the implementation is 

235 trivial). 

236 

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

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

239 

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

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

242 contain entries supported by 

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

244 compression). 

245 

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

247 specialization): 

248 

249 .. code-block:: yaml 

250 

251 lsst.obs.base.fitsExposureFormatter.FitsExposureFormatter: 

252 default: 

253 image: &default 

254 algorithm: GZIP_2 

255 mask: *default 

256 variance: *default 

257 

258 """ 

259 

260 supported_write_parameters = frozenset({"recipe"}) 

261 

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

263 # Docstring inherited. 

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

265 bbox = self.reader.readBBox() 

266 if component == "dimensions": 

267 return bbox.getDimensions() 

268 elif component == "xy0": 

269 return bbox.getMin() 

270 else: 

271 return bbox 

272 else: 

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

274 

275 def readFull(self) -> Any: 

276 # Docstring inherited. 

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

278 

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

280 """Serialize the image to FITS. 

281 

282 Parameters 

283 ---------- 

284 in_memory_dataset : `object` 

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

286 ``writeFitsWithOptions`` interface. 

287 uri : `lsst.resources.ResourcePath` 

288 Location to write the local file. 

289 """ 

290 # check to see if we have a recipe requested 

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

292 recipe = self.get_image_compression_settings(recipeName) 

293 if recipe: 

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

295 else: 

296 in_memory_dataset.writeFits(uri.ospath) 

297 

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

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

300 

301 Parameters 

302 ---------- 

303 recipeName : `str` or `None` 

304 Label associated with the collection of compression parameters 

305 to select. 

306 

307 Returns 

308 ------- 

309 settings : `dict` 

310 The selected settings. 

311 """ 

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

313 # return immediately 

314 if not recipeName: 

315 if "default" not in self.write_recipes: 

316 return {} 

317 recipeName = "default" 

318 

319 if recipeName not in self.write_recipes: 

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

321 

322 recipe = self.write_recipes[recipeName] 

323 if recipe is None: 

324 return {} 

325 seed: int | None = None 

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

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

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

329 if seed is None: 

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

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

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

333 # are only defined for certain combinations of 

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

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

336 # simpler. 

337 hash_bytes = hashlib.md5( 

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

339 usedforsecurity=False, 

340 ).digest() 

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

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

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

344 _LOG.debug( 

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

346 self.data_id, 

347 self.dataset_ref.datasetType.name, 

348 plane, 

349 seed, 

350 ) 

351 quantization["seed"] = seed 

352 else: 

353 _LOG.warning( 

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

355 self.dataset_ref.datasetType.name, 

356 self.data_id, 

357 plane, 

358 quantization["seed"], 

359 ) 

360 else: 

361 _LOG.debug( 

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

363 self.dataset_ref.datasetType.name, 

364 self.data_id, 

365 plane, 

366 ) 

367 return recipe 

368 

369 @classmethod 

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

371 """Validate supplied recipes for this formatter. 

372 

373 The recipes are supplemented with default values where appropriate. 

374 

375 Parameters 

376 ---------- 

377 recipes : `dict` or `None` 

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

379 

380 Returns 

381 ------- 

382 validated : `dict` 

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

384 recipes listed. 

385 

386 Raises 

387 ------ 

388 RuntimeError 

389 Raised if validation fails. 

390 """ 

391 if not recipes: 

392 # We can not insist on recipes being specified. 

393 return recipes 

394 

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

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

397 if recipe is not None: 

398 validated[name] = {} 

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

400 try: 

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

402 except Exception as err: 

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

404 raise 

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

406 else: 

407 validated[name] = None 

408 return validated 

409 

410 

411class FitsImageFormatter(StandardFitsImageFormatterBase): 

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

413 from/to FITS. 

414 """ 

415 

416 ReaderClass = ImageFitsReader 

417 

418 

419class FitsMaskFormatter(StandardFitsImageFormatterBase): 

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

421 from/to FITS. 

422 """ 

423 

424 ReaderClass = MaskFitsReader 

425 

426 

427class FitsMaskedImageFormatter(StandardFitsImageFormatterBase): 

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

429 from/to FITS. 

430 """ 

431 

432 ReaderClass = MaskedImageFitsReader 

433 

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

435 # Docstring inherited. 

436 if component == "image": 

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

438 elif component == "mask": 

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

440 elif component == "variance": 

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

442 else: 

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

444 return super().readComponent(component) 

445 

446 

447def standardizeAmplifierParameters( 

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

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

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

451 

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

453 converts amplifier IDs/names to Amplifier instances. 

454 

455 Parameters 

456 ---------- 

457 parameters : `dict` 

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

459 class definition in daf_butler for allowed keys and values. 

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

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

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

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

464 

465 Returns 

466 ------- 

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

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

469 was no "amp" parameter. 

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

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

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

473 ``on_disk_detector`` is simply passed through. 

474 regions_differ : `bool` 

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

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

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

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

479 """ 

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

481 return None, on_disk_detector, False 

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

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

484 if isinstance(amplifier, int | str): 

485 amp_key = amplifier 

486 target_amplifier = None 

487 else: 

488 amp_key = amplifier.getName() 

489 target_amplifier = amplifier 

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

491 if on_disk_detector is not None: 

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

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

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

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

496 # assembly/orientation state. 

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

498 if comparison & comparison.ASSEMBLY_DIFFERS: 

499 raise ValueError( 

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

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

502 ) 

503 else: 

504 if on_disk_detector is None: 

505 raise ValueError( 

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

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

508 "'amp' parameter." 

509 ) 

510 comparison = AmplifierGeometryComparison.EQUAL 

511 detector = on_disk_detector 

512 if target_amplifier is None: 

513 target_amplifier = detector[amp_key] 

514 return target_amplifier, detector, comparison & comparison.REGIONS_DIFFER 

515 

516 

517class _ComponentCache(NamedTuple): 

518 id_: uuid.UUID | None = None 

519 reader: ExposureFitsReader | None = None 

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

521 mem: MemFileManager | None = None 

522 

523 

524class FitsExposureFormatter(FitsMaskedImageFormatter): 

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

526 from/to FITS. 

527 

528 Notes 

529 ----- 

530 This class inherits from `FitsMaskedImageFormatter` even though 

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

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

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

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

535 type covariance violation ever becomes a practical problem. 

536 """ 

537 

538 can_read_from_uri = True 

539 ReaderClass = ExposureFitsReader 

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

541 _lock = threading.Lock() 

542 _cached_fits: _ComponentCache = _ComponentCache() 

543 

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

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

546 # could work with cutouts. 

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

548 

549 # Full read, always use local file read. 

550 if not component and not self.checked_parameters: 

551 return NotImplemented 

552 

553 if not _ALWAYS_USE_ASTROPY_FOR_COMPONENT_READ and uri.isLocal: 

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

555 return NotImplemented 

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

557 

558 if component in pixel_components: 

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

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

561 # read with file caching for now. 

562 return NotImplemented 

563 

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

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

566 # file. 

567 if not ( 

568 _ALWAYS_USE_ASTROPY_FOR_COMPONENT_READ 

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

570 ): 

571 return NotImplemented 

572 

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

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

575 # bounding box is reasonable. 

576 bbox = None 

577 origin = lsst.afw.image.PARENT 

578 if not component: 

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

580 # other parameters do not attempt a cutout. 

581 if self.checked_parameters.keys() - {"bbox", "origin"}: 

582 return NotImplemented 

583 bbox = self.checked_parameters["bbox"] 

584 origin = self.checked_parameters.get("origin", lsst.afw.image.PARENT) 

585 # For larger cutouts use the full file. 

586 max_cutout_size = 500 * 500 

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

588 return NotImplemented 

589 

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

591 if component: 

592 with self._lock: 

593 cache = type(self)._cached_fits 

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

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

596 match component: 

597 case "xy0": 

598 return cache.bbox.getMin() 

599 case "dimensions": 

600 return cache.bbox.getDimensions() 

601 case "bbox": 

602 return cache.bbox 

603 else: 

604 self._reader = cache.reader 

605 return self.readComponent(component) 

606 

607 try: 

608 fs, fspath = uri.to_fsspec() 

609 except Exception: 

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

611 return NotImplemented 

612 

613 bbox_component = None 

614 try: 

615 hdul = [] 

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

617 # Read all non-pixel components and cache. 

618 for hdu in fits_obj: 

619 hdr = hdu.header 

620 extname = hdr.get("EXTNAME") 

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

622 # for EXTEND=T. 

623 extend = hdr.get("EXTEND") 

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

625 # Calculate the dimensional components for later 

626 # caching. Do not derive from cached FITS reader 

627 # because they depend on the dimensionality of the 

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

629 # data. 

630 if bbox_component is None: 

631 shape = hdu.shape 

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

633 

634 # XY0 is defined in the A WCS. 

635 pl = PropertyList() 

636 pl.update(hdr) 

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

638 

639 # This is the PARENT bbox. 

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

641 

642 # Handle cutout request. 

643 if bbox: 

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

645 full_bbox = bbox_component 

646 else: 

647 full_bbox = lsst.geom.Box2I( 

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

649 ) 

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

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

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

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

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

655 

656 # Must correct the header WCS to take into 

657 # account the offset. 

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

659 hdr[k] -= minX 

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

661 hdr[k] -= minY 

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

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

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

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

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

667 hdr[k] = bbox.getBeginX() 

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

669 hdr[k] = bbox.getBeginY() 

670 else: 

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

672 

673 # Construct a new HDU and copy the header. 

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

675 hdul.append(stripped_hdu) 

676 else: 

677 hdul.append(hdu) 

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

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

680 buffer = BytesIO() 

681 stripped_fits.writeto(buffer) 

682 except Exception as e: 

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

684 _LOG.debug( 

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

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

687 str(e), 

688 ) 

689 return NotImplemented 

690 

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

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

693 # components. 

694 fits_data = buffer.getvalue() 

695 mem = MemFileManager(len(fits_data)) 

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

697 self._reader = self.ReaderClass(mem) 

698 

699 if component: 

700 with self._lock: 

701 type(self)._cached_fits = _ComponentCache( 

702 id_=self.dataset_ref.id, 

703 reader=self._reader, 

704 mem=mem, 

705 bbox=bbox_component, 

706 ) 

707 match component: 

708 case "xy0": 

709 if bbox_component is None: # For mypy. 

710 return None 

711 return bbox_component.getMin() 

712 case "dimensions": 

713 if bbox_component is None: 

714 return None 

715 return bbox_component.getDimensions() 

716 case "bbox": 

717 return bbox_component 

718 case _: 

719 return self.readComponent(component) 

720 else: 

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

722 # parameters should be passed here. 

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

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

725 return cutout 

726 

727 def add_provenance( 

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

729 ) -> Any: 

730 # Add provenance via FITS headers. 

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

732 return in_memory_dataset 

733 

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

735 # Docstring inherited. 

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

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

738 genericComponents = { 

739 "summaryStats": ExposureInfo.KEY_SUMMARY_STATS, 

740 } 

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

742 return self.reader.readComponent(genericComponentName) 

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

744 # parameters. 

745 standardComponents = { 

746 "id": "readExposureId", 

747 "metadata": "readMetadata", 

748 "wcs": "readWcs", 

749 "coaddInputs": "readCoaddInputs", 

750 "psf": "readPsf", 

751 "photoCalib": "readPhotoCalib", 

752 "filter": "readFilter", 

753 "validPolygon": "readValidPolygon", 

754 "apCorrMap": "readApCorrMap", 

755 "visitInfo": "readVisitInfo", 

756 "transmissionCurve": "readTransmissionCurve", 

757 "detector": "readDetector", 

758 "exposureInfo": "readExposureInfo", 

759 } 

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

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

762 if component == "filter": 

763 return self._fixFilterLabels(result) 

764 return result 

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

766 return super().readComponent(component) 

767 

768 def readFull(self) -> Any: 

769 # Docstring inherited. 

770 amplifier, detector, _ = standardizeAmplifierParameters( 

771 self.checked_parameters, 

772 self.reader.readDetector(), 

773 ) 

774 if amplifier is not None: 

775 amplifier_isolator = AmplifierIsolator( 

776 amplifier, 

777 self.reader.readBBox(), 

778 detector, 

779 ) 

780 result = amplifier_isolator.transform_subimage( 

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

782 ) 

783 result.setDetector(amplifier_isolator.make_detector()) 

784 else: 

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

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

787 return result 

788 

789 def _fixFilterLabels( 

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

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

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

793 data ID. 

794 

795 Parameters 

796 ---------- 

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

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

799 should_be_standardized : `bool`, optional 

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

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

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

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

804 whether the file should be standardized by looking at the 

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

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

807 

808 Returns 

809 ------- 

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

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

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

813 filters associated with this dataset type. 

814 

815 Notes 

816 ----- 

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

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

819 to solve. 

820 """ 

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

822 # so we can warn about them later. 

823 missing = [] 

824 band = None 

825 physical_filter = None 

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

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

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

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

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

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

832 missing.append("band") 

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

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

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

836 if ( 

837 physical_filter is None 

838 and not self.data_id.hasFull() 

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

840 ): 

841 missing.append("physical_filter") 

842 if should_be_standardized is None: 

843 version = self.reader.readSerializationVersion() 

844 should_be_standardized = version >= 2 

845 if missing: 

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

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

848 # to use the one in the file. 

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

850 # predates filter standardization. 

851 if not should_be_standardized: 

852 warnings.warn( 

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

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

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

856 # Report the warning from outside of middleware or the 

857 # relevant runQuantum method. 

858 stacklevel=find_outside_stacklevel( 

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

860 ), 

861 ) 

862 return file_filter_label 

863 if band is None and physical_filter is None: 

864 data_id_filter_label = None 

865 else: 

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

867 if data_id_filter_label != file_filter_label and should_be_standardized: 

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

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

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

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

872 warnings.warn( 

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

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

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

876 stacklevel=find_outside_stacklevel( 

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

878 ), 

879 ) 

880 return data_id_filter_label