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

341 statements  

« prev     ^ index     » next       coverage.py v7.13.5, created at 2026-05-07 08:21 +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 self.dataset_ref.datasetType.storageClass_name == "VisitImage": 

169 from lsst.images import VisitImage 

170 

171 return VisitImage.read_legacy(path, component=component, preserve_quantization=True) 

172 

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

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

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

176 # file. 

177 self._reader_path = path 

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

179 try: 

180 if component is not None: 

181 in_memory_dataset = self.readComponent(component) 

182 else: 

183 in_memory_dataset = self.readFull() 

184 finally: 

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

186 return in_memory_dataset 

187 

188 @abstractmethod 

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

190 """Read a component dataset. 

191 

192 Parameters 

193 ---------- 

194 component : `str`, optional 

195 Component to read from the file. 

196 

197 Returns 

198 ------- 

199 obj : `typing.Any` 

200 In-memory component object. 

201 

202 Raises 

203 ------ 

204 KeyError 

205 Raised if the requested component cannot be handled. 

206 """ 

207 raise NotImplementedError() 

208 

209 @abstractmethod 

210 def readFull(self) -> Any: 

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

212 

213 Returns 

214 ------- 

215 obj : `typing.Any` 

216 In-memory component object. 

217 

218 """ 

219 raise NotImplementedError() 

220 

221 

222class StandardFitsImageFormatterBase(FitsImageFormatterBase): 

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

224 written using LSST code. 

225 

226 Notes 

227 ----- 

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

229 

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

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

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

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

234 and catch `KeyError`). 

235 

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

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

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

239 should generally reimplement without delegating (the implementation is 

240 trivial). 

241 

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

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

244 

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

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

247 contain entries supported by 

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

249 compression). 

250 

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

252 specialization): 

253 

254 .. code-block:: yaml 

255 

256 lsst.obs.base.fitsExposureFormatter.FitsExposureFormatter: 

257 default: 

258 image: &default 

259 algorithm: GZIP_2 

260 mask: *default 

261 variance: *default 

262 

263 """ 

264 

265 supported_write_parameters = frozenset({"recipe"}) 

266 

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

268 # Docstring inherited. 

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

270 bbox = self.reader.readBBox() 

271 if component == "dimensions": 

272 return bbox.getDimensions() 

273 elif component == "xy0": 

274 return bbox.getMin() 

275 else: 

276 return bbox 

277 else: 

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

279 

280 def readFull(self) -> Any: 

281 # Docstring inherited. 

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

283 

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

285 """Serialize the image to FITS. 

286 

287 Parameters 

288 ---------- 

289 in_memory_dataset : `object` 

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

291 ``writeFitsWithOptions`` interface. 

292 uri : `lsst.resources.ResourcePath` 

293 Location to write the local file. 

294 """ 

295 # check to see if we have a recipe requested 

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

297 recipe = self.get_image_compression_settings(recipeName) 

298 if recipe: 

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

300 else: 

301 in_memory_dataset.writeFits(uri.ospath) 

302 

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

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

305 

306 Parameters 

307 ---------- 

308 recipeName : `str` or `None` 

309 Label associated with the collection of compression parameters 

310 to select. 

311 

312 Returns 

313 ------- 

314 settings : `dict` 

315 The selected settings. 

316 """ 

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

318 # return immediately 

319 if not recipeName: 

320 if "default" not in self.write_recipes: 

321 return {} 

322 recipeName = "default" 

323 

324 if recipeName not in self.write_recipes: 

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

326 

327 recipe = self.write_recipes[recipeName] 

328 if recipe is None: 

329 return {} 

330 seed: int | None = None 

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

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

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

334 if seed is None: 

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

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

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

338 # are only defined for certain combinations of 

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

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

341 # simpler. 

342 hash_bytes = hashlib.md5( 

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

344 usedforsecurity=False, 

345 ).digest() 

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

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

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

349 _LOG.debug( 

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

351 self.data_id, 

352 self.dataset_ref.datasetType.name, 

353 plane, 

354 seed, 

355 ) 

356 quantization["seed"] = seed 

357 else: 

358 _LOG.warning( 

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

360 self.dataset_ref.datasetType.name, 

361 self.data_id, 

362 plane, 

363 quantization["seed"], 

364 ) 

365 else: 

366 _LOG.debug( 

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

368 self.dataset_ref.datasetType.name, 

369 self.data_id, 

370 plane, 

371 ) 

372 return recipe 

373 

374 @classmethod 

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

376 """Validate supplied recipes for this formatter. 

377 

378 The recipes are supplemented with default values where appropriate. 

379 

380 Parameters 

381 ---------- 

382 recipes : `dict` or `None` 

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

384 

385 Returns 

386 ------- 

387 validated : `dict` 

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

389 recipes listed. 

390 

391 Raises 

392 ------ 

393 RuntimeError 

394 Raised if validation fails. 

395 """ 

396 if not recipes: 

397 # We can not insist on recipes being specified. 

398 return recipes 

399 

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

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

402 if recipe is not None: 

403 validated[name] = {} 

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

405 try: 

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

407 except Exception as err: 

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

409 raise 

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

411 else: 

412 validated[name] = None 

413 return validated 

414 

415 

416class FitsImageFormatter(StandardFitsImageFormatterBase): 

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

418 from/to FITS. 

419 """ 

420 

421 ReaderClass = ImageFitsReader 

422 

423 

424class FitsMaskFormatter(StandardFitsImageFormatterBase): 

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

426 from/to FITS. 

427 """ 

428 

429 ReaderClass = MaskFitsReader 

430 

431 

432class FitsMaskedImageFormatter(StandardFitsImageFormatterBase): 

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

434 from/to FITS. 

435 """ 

436 

437 ReaderClass = MaskedImageFitsReader 

438 

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

440 # Docstring inherited. 

441 if component == "image": 

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

443 elif component == "mask": 

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

445 elif component == "variance": 

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

447 else: 

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

449 return super().readComponent(component) 

450 

451 

452def standardizeAmplifierParameters( 

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

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

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

456 

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

458 converts amplifier IDs/names to Amplifier instances. 

459 

460 Parameters 

461 ---------- 

462 parameters : `dict` 

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

464 class definition in daf_butler for allowed keys and values. 

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

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

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

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

469 

470 Returns 

471 ------- 

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

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

474 was no "amp" parameter. 

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

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

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

478 ``on_disk_detector`` is simply passed through. 

479 regions_differ : `bool` 

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

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

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

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

484 """ 

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

486 return None, on_disk_detector, False 

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

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

489 if isinstance(amplifier, int | str): 

490 amp_key = amplifier 

491 target_amplifier = None 

492 else: 

493 amp_key = amplifier.getName() 

494 target_amplifier = amplifier 

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

496 if on_disk_detector is not None: 

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

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

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

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

501 # assembly/orientation state. 

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

503 if comparison & comparison.ASSEMBLY_DIFFERS: 

504 raise ValueError( 

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

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

507 ) 

508 else: 

509 if on_disk_detector is None: 

510 raise ValueError( 

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

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

513 "'amp' parameter." 

514 ) 

515 comparison = AmplifierGeometryComparison.EQUAL 

516 detector = on_disk_detector 

517 if target_amplifier is None: 

518 target_amplifier = detector[amp_key] 

519 return target_amplifier, detector, comparison & comparison.REGIONS_DIFFER 

520 

521 

522class _ComponentCache(NamedTuple): 

523 id_: uuid.UUID | None = None 

524 reader: ExposureFitsReader | None = None 

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

526 mem: MemFileManager | None = None 

527 

528 

529class FitsExposureFormatter(FitsMaskedImageFormatter): 

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

531 from/to FITS. 

532 

533 Notes 

534 ----- 

535 This class inherits from `FitsMaskedImageFormatter` even though 

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

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

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

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

540 type covariance violation ever becomes a practical problem. 

541 """ 

542 

543 can_read_from_uri = True 

544 ReaderClass = ExposureFitsReader 

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

546 _lock = threading.Lock() 

547 _cached_fits: _ComponentCache = _ComponentCache() 

548 

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

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

551 # could work with cutouts. 

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

553 

554 # Full read, always use local file read. 

555 if not component and not self.checked_parameters: 

556 return NotImplemented 

557 

558 if not _ALWAYS_USE_ASTROPY_FOR_COMPONENT_READ and uri.isLocal: 

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

560 return NotImplemented 

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

562 

563 if component in pixel_components: 

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

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

566 # read with file caching for now. 

567 return NotImplemented 

568 

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

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

571 # file. 

572 if not ( 

573 _ALWAYS_USE_ASTROPY_FOR_COMPONENT_READ 

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

575 ): 

576 return NotImplemented 

577 

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

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

580 # bounding box is reasonable. 

581 bbox = None 

582 origin = lsst.afw.image.PARENT 

583 if not component: 

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

585 # other parameters do not attempt a cutout. 

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

587 return NotImplemented 

588 bbox = self.checked_parameters["bbox"] 

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

590 # For larger cutouts use the full file. 

591 max_cutout_size = 500 * 500 

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

593 return NotImplemented 

594 

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

596 if component: 

597 with self._lock: 

598 cache = type(self)._cached_fits 

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

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

601 match component: 

602 case "xy0": 

603 return cache.bbox.getMin() 

604 case "dimensions": 

605 return cache.bbox.getDimensions() 

606 case "bbox": 

607 return cache.bbox 

608 else: 

609 self._reader = cache.reader 

610 return self.readComponent(component) 

611 

612 try: 

613 fs, fspath = uri.to_fsspec() 

614 except Exception: 

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

616 return NotImplemented 

617 

618 bbox_component = None 

619 try: 

620 hdul = [] 

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

622 # Read all non-pixel components and cache. 

623 for hdu in fits_obj: 

624 hdr = hdu.header 

625 extname = hdr.get("EXTNAME") 

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

627 # for EXTEND=T. 

628 extend = hdr.get("EXTEND") 

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

630 # Calculate the dimensional components for later 

631 # caching. Do not derive from cached FITS reader 

632 # because they depend on the dimensionality of the 

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

634 # data. 

635 if bbox_component is None: 

636 shape = hdu.shape 

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

638 

639 # XY0 is defined in the A WCS. 

640 pl = PropertyList() 

641 pl.update(hdr) 

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

643 

644 # This is the PARENT bbox. 

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

646 

647 # Handle cutout request. 

648 if bbox: 

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

650 full_bbox = bbox_component 

651 else: 

652 full_bbox = lsst.geom.Box2I( 

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

654 ) 

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

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

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

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

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

660 

661 # Must correct the header WCS to take into 

662 # account the offset. 

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

664 hdr[k] -= minX 

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

666 hdr[k] -= minY 

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

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

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

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

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

672 hdr[k] = bbox.getBeginX() 

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

674 hdr[k] = bbox.getBeginY() 

675 else: 

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

677 

678 # Construct a new HDU and copy the header. 

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

680 hdul.append(stripped_hdu) 

681 else: 

682 hdul.append(hdu) 

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

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

685 buffer = BytesIO() 

686 stripped_fits.writeto(buffer) 

687 except Exception as e: 

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

689 _LOG.debug( 

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

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

692 str(e), 

693 ) 

694 return NotImplemented 

695 

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

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

698 # components. 

699 fits_data = buffer.getvalue() 

700 mem = MemFileManager(len(fits_data)) 

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

702 self._reader = self.ReaderClass(mem) 

703 

704 if component: 

705 with self._lock: 

706 type(self)._cached_fits = _ComponentCache( 

707 id_=self.dataset_ref.id, 

708 reader=self._reader, 

709 mem=mem, 

710 bbox=bbox_component, 

711 ) 

712 match component: 

713 case "xy0": 

714 if bbox_component is None: # For mypy. 

715 return None 

716 return bbox_component.getMin() 

717 case "dimensions": 

718 if bbox_component is None: 

719 return None 

720 return bbox_component.getDimensions() 

721 case "bbox": 

722 return bbox_component 

723 case _: 

724 return self.readComponent(component) 

725 else: 

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

727 # parameters should be passed here. 

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

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

730 return cutout 

731 

732 def add_provenance( 

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

734 ) -> Any: 

735 # Add provenance via FITS headers. 

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

737 return in_memory_dataset 

738 

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

740 # Docstring inherited. 

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

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

743 genericComponents = { 

744 "summaryStats": ExposureInfo.KEY_SUMMARY_STATS, 

745 } 

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

747 return self.reader.readComponent(genericComponentName) 

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

749 # parameters. 

750 standardComponents = { 

751 "id": "readExposureId", 

752 "metadata": "readMetadata", 

753 "wcs": "readWcs", 

754 "coaddInputs": "readCoaddInputs", 

755 "psf": "readPsf", 

756 "photoCalib": "readPhotoCalib", 

757 "filter": "readFilter", 

758 "validPolygon": "readValidPolygon", 

759 "apCorrMap": "readApCorrMap", 

760 "visitInfo": "readVisitInfo", 

761 "transmissionCurve": "readTransmissionCurve", 

762 "detector": "readDetector", 

763 "exposureInfo": "readExposureInfo", 

764 } 

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

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

767 if component == "filter": 

768 return self._fixFilterLabels(result) 

769 return result 

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

771 return super().readComponent(component) 

772 

773 def readFull(self) -> Any: 

774 # Docstring inherited. 

775 amplifier, detector, _ = standardizeAmplifierParameters( 

776 self.checked_parameters, 

777 self.reader.readDetector(), 

778 ) 

779 if amplifier is not None: 

780 amplifier_isolator = AmplifierIsolator( 

781 amplifier, 

782 self.reader.readBBox(), 

783 detector, 

784 ) 

785 result = amplifier_isolator.transform_subimage( 

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

787 ) 

788 result.setDetector(amplifier_isolator.make_detector()) 

789 else: 

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

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

792 return result 

793 

794 def _fixFilterLabels( 

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

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

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

798 data ID. 

799 

800 Parameters 

801 ---------- 

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

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

804 should_be_standardized : `bool`, optional 

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

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

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

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

809 whether the file should be standardized by looking at the 

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

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

812 

813 Returns 

814 ------- 

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

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

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

818 filters associated with this dataset type. 

819 

820 Notes 

821 ----- 

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

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

824 to solve. 

825 """ 

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

827 # so we can warn about them later. 

828 missing = [] 

829 band = None 

830 physical_filter = None 

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

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

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

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

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

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

837 missing.append("band") 

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

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

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

841 if ( 

842 physical_filter is None 

843 and not self.data_id.hasFull() 

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

845 ): 

846 missing.append("physical_filter") 

847 if should_be_standardized is None: 

848 version = self.reader.readSerializationVersion() 

849 should_be_standardized = version >= 2 

850 if missing: 

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

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

853 # to use the one in the file. 

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

855 # predates filter standardization. 

856 if not should_be_standardized: 

857 warnings.warn( 

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

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

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

861 # Report the warning from outside of middleware or the 

862 # relevant runQuantum method. 

863 stacklevel=find_outside_stacklevel( 

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

865 ), 

866 ) 

867 return file_filter_label 

868 if band is None and physical_filter is None: 

869 data_id_filter_label = None 

870 else: 

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

872 if data_id_filter_label != file_filter_label and should_be_standardized: 

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

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

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

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

877 warnings.warn( 

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

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

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

881 stacklevel=find_outside_stacklevel( 

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

883 ), 

884 ) 

885 return data_id_filter_label