Coverage for python/lsst/images/fits/_common.py: 54%

263 statements  

« prev     ^ index     » next       coverage.py v7.16.2, created at 2026-10-09 10:54 +0000

1# This file is part of lsst-images. 

2# 

3# Developed for the LSST Data Management System. 

4# This product includes software developed by the LSST Project 

5# (https://www.lsst.org). 

6# See the COPYRIGHT file at the top-level directory of this distribution 

7# for details of code ownership. 

8# 

9# Use of this source code is governed by a 3-clause BSD-style 

10# license that can be found in the LICENSE file. 

11 

12 

13from __future__ import annotations 

14 

15__all__ = ( 

16 "FITS_SOURCE_REGEX", 

17 "JSON_COLUMN", 

18 "JSON_EXTNAME", 

19 "ExtensionHDU", 

20 "ExtensionKey", 

21 "FitsCompressionAlgorithm", 

22 "FitsCompressionOptions", 

23 "FitsDitherAlgorithm", 

24 "FitsExternalMetadata", 

25 "FitsOpaqueMetadata", 

26 "FitsQuantizationOptions", 

27 "InvalidFitsArchiveError", 

28 "PointerModel", 

29 "PrecompressedImage", 

30 "add_offset_wcs", 

31 "parse_legacy_bunit", 

32 "read_offset_wcs", 

33 "read_yx0", 

34 "strip_butler_cards", 

35 "strip_legacy_exposure_cards", 

36 "strip_wcs_cards", 

37 "suppress_fits_card_warnings", 

38) 

39 

40import contextlib 

41import dataclasses 

42import enum 

43import itertools 

44import re 

45import string 

46import warnings 

47from collections.abc import Iterator 

48from typing import Any, ClassVar, Self, final 

49 

50import astropy.io.fits 

51import astropy.io.fits.verify 

52import numpy as np 

53import pydantic 

54 

55from .._geom import YX, Box 

56from ..serialization import ( 

57 ArchiveReadError, 

58 ExternalMetadata, 

59 ExternalMetadataValue, 

60 OpaqueArchiveMetadata, 

61 TableColumnModel, 

62) 

63 

64type ExtensionHDU = astropy.io.fits.ImageHDU | astropy.io.fits.CompImageHDU | astropy.io.fits.BinTableHDU 

65 

66FITS_SOURCE_REGEX = re.compile(r"fits:(?P<extname>[\w/\-]+)(,(?P<extver>\d+))?(\[\d+\])?") 

67 

68JSON_EXTNAME: str = "JSON" 

69JSON_COLUMN: str = "JSON" 

70 

71 

72@contextlib.contextmanager 

73def suppress_fits_card_warnings() -> Iterator[None]: 

74 """Silence Astropy warnings for FITS cards it fixes up automatically. 

75 

76 When a header is written, Astropy converts keywords longer than eight 

77 characters into ``HIERARCH`` cards and truncates comments that overflow a 

78 card. Both emit a `~astropy.io.fits.verify.VerifyWarning` that carries no 

79 actionable information, so they are silenced around the code that writes 

80 these headers. 

81 """ 

82 with warnings.catch_warnings(): 

83 warnings.filterwarnings( 

84 "ignore", 

85 message="Keyword name .* is greater than 8 characters", 

86 category=astropy.io.fits.verify.VerifyWarning, 

87 ) 

88 warnings.filterwarnings( 

89 "ignore", 

90 message=".*comment will be truncated", 

91 category=astropy.io.fits.verify.VerifyWarning, 

92 ) 

93 yield 

94 

95 

96@dataclasses.dataclass(frozen=True) 

97class ExtensionKey: 

98 """The identifiers for a single FITS extension HDU within a file.""" 

99 

100 name: str = "" 

101 """Value of the EXTNAME keyword, or an empty string for the primary HDU.""" 

102 

103 ver: int = 1 

104 """Value of the EXTVER keyword (which may be absent if it is ``1``).""" 

105 

106 @classmethod 

107 def from_index_row(cls, row: np.ndarray) -> ExtensionKey: 

108 """Construct from a row of the index binary table appended to the 

109 FITS files written by the package. 

110 

111 Parameters 

112 ---------- 

113 row 

114 Structured-array row with ``EXTNAME`` and ``EXTVER`` fields. 

115 """ 

116 return cls(str(row["EXTNAME"]).upper(), int(row["EXTVER"])) 

117 

118 @classmethod 

119 def from_str(cls, source: str) -> ExtensionKey: 

120 """Construct from the `str` coercion of this type, which is used 

121 as the 'source' field in various Pydantic models that serve as 

122 references to other HDUs. 

123 

124 Parameters 

125 ---------- 

126 source 

127 ``fits:<extname>[,<extver>]`` reference string, as produced by 

128 ``str(self)``. 

129 """ 

130 if (m := FITS_SOURCE_REGEX.fullmatch(source)) is None: 130 ↛ 131line 130 didn't jump to line 131 because the condition on line 130 was never true

131 raise ArchiveReadError(f"Bad 'source' string for FITS: {source!r}.") 

132 ver = 1 

133 if m.group("extver") is not None: 

134 ver = int(m.group("extver")) 

135 return cls(m.group("extname"), ver) 

136 

137 def check(self) -> None: 

138 if not FITS_SOURCE_REGEX.match(str(self)): 138 ↛ 139line 138 didn't jump to line 139 because the condition on line 138 was never true

139 raise ValueError( 

140 f"Invalid source key: '{str(self)}'; name characters must be alphanumeric, '-', '_', or '/'." 

141 ) 

142 

143 def __str__(self) -> str: 

144 if self.ver > 1: 

145 return f"fits:{self.name},{self.ver}" 

146 else: 

147 return f"fits:{self.name}" 

148 

149 

150class PointerModel(pydantic.BaseModel): 

151 column: TableColumnModel = pydantic.Field(description="Table column for this cell this pointer targets.") 

152 row: int = pydantic.Field(description="Zero-indexed row for cell this pointer targets.") 

153 

154 

155class InvalidFitsArchiveError(RuntimeError): 

156 """The error type raised when the content of a FITS file presumed to have 

157 been written by FitsOutputArchive is not self-consistent. 

158 """ 

159 

160 

161class FitsCompressionAlgorithm(enum.StrEnum): 

162 """FITS compression algorithms supported by this package. 

163 

164 See the FITS standard for definitions. 

165 """ 

166 

167 GZIP_1 = "GZIP_1" 

168 GZIP_2 = "GZIP_2" 

169 RICE_1 = "RICE_1" 

170 

171 

172class FitsDitherAlgorithm(enum.StrEnum): 

173 """FITS quantization dither algorithms supported by this package. 

174 

175 See the FITS standard for definitions. 

176 """ 

177 

178 NO_DITHER = "NO_DITHER" 

179 SUBTRACTIVE_DITHER_1 = "SUBTRACTIVE_DITHER_1" 

180 SUBTRACTIVE_DITHER_2 = "SUBTRACTIVE_DITHER_2" 

181 

182 def to_astropy_quantize_method(self) -> int: 

183 """Convert to the integer code used by Astropy.""" 

184 match self: 

185 case self.NO_DITHER: 185 ↛ 186line 185 didn't jump to line 186 because the pattern on line 185 never matched

186 return -1 

187 case self.SUBTRACTIVE_DITHER_1: 187 ↛ 188line 187 didn't jump to line 188 because the pattern on line 187 never matched

188 return 1 

189 case self.SUBTRACTIVE_DITHER_2: 

190 return 2 

191 raise AssertionError("Invalid enum value.") 

192 

193 

194class FitsQuantizationOptions(pydantic.BaseModel, frozen=True): 

195 """Quantization options for FITS compression.""" 

196 

197 dither: FitsDitherAlgorithm 

198 """How to add random noise during quantization to reduce biases.""" 

199 

200 level: float 

201 """Quantization level. 

202 

203 When positive, this is the fraction of the measured standard deviation that 

204 corresponds to an integer step. When negative, it is ``-ZSCALE``, the 

205 scaling to apply directly to the original pixels before quantization. 

206 """ 

207 

208 seed: int | None = None 

209 """Random number seed to use for dithering. 

210 

211 Values between 1 and 10000 (inclusive) are used directly. ``0`` will 

212 generate a value from the current time, and ``-1`` will generate a value 

213 from the checksum of the image. 

214 

215 If `None`, the ``compression_seed`` parameter must be passed to 

216 `FitsOutputArchive.open` if any quantized compression is configured. 

217 """ 

218 

219 

220class FitsCompressionOptions(pydantic.BaseModel, frozen=True): 

221 """Configuration options for FITS compression.""" 

222 

223 algorithm: FitsCompressionAlgorithm = FitsCompressionAlgorithm.GZIP_2 

224 """Compression algorithm to use.""" 

225 

226 tile_shape: tuple[int, ...] | None = None 

227 """Shape ``(..., y, x)`` of independently compressed tiles. 

228 

229 The default of `None` leaves the tile shape up to the ``tile_shape`` 

230 argument to `~lsst.images.serialization.OutputArchive.add_array`. 

231 """ 

232 

233 quantization: FitsQuantizationOptions | None = None 

234 """Quantization to apply before compression, if any.""" 

235 

236 DEFAULT: ClassVar[FitsCompressionOptions | None] 

237 """Default compression options (lossless ``GZIP_2``).""" 

238 

239 LOSSY: ClassVar[FitsCompressionOptions] 

240 """Default lossy compression options.""" 

241 

242 def make_hdu(self, data: np.ndarray, name: str) -> astropy.io.fits.CompImageHDU: 

243 """Make an `astropy.io.fits.CompImageHDU` object from these options. 

244 

245 Parameters 

246 ---------- 

247 data 

248 Pixel data to store in the HDU. 

249 name 

250 ``EXTNAME`` of the HDU. 

251 """ 

252 if self.quantization is not None: 

253 return astropy.io.fits.CompImageHDU( 

254 data, 

255 name=name, 

256 compression_type=self.algorithm.value, 

257 tile_shape=self.tile_shape, 

258 quantize_method=self.quantization.dither.to_astropy_quantize_method(), 

259 quantize_level=self.quantization.level, 

260 dither_seed=self.quantization.seed, 

261 ) 

262 else: 

263 return astropy.io.fits.CompImageHDU( 

264 data, 

265 name=name, 

266 compression_type=self.algorithm.value, 

267 tile_shape=self.tile_shape, 

268 quantize_level=0.0, 

269 ) 

270 

271 

272FitsCompressionOptions.DEFAULT = FitsCompressionOptions() 

273FitsCompressionOptions.LOSSY = FitsCompressionOptions( 

274 algorithm=FitsCompressionAlgorithm.RICE_1, 

275 quantization=FitsQuantizationOptions(dither=FitsDitherAlgorithm.SUBTRACTIVE_DITHER_2, level=16.0), 

276) 

277 

278 

279_COMPRESSION_KEYS = frozenset( 

280 ( 

281 "ZIMAGE", 

282 "ZCMPTYPE", 

283 "ZBITPIX", 

284 "ZNAXIS", 

285 "ZMASKCMP", 

286 "ZQUANTIZ", 

287 "ZDITHER0", 

288 "ZCHECKSUM", 

289 "ZDATASUM", 

290 ) 

291) 

292_COMPRESSION_PREFIX_KEYS = ("ZNAXIS", "ZTILE", "ZNAME", "ZVAL") 

293 

294 

295@dataclasses.dataclass 

296class PrecompressedImage: 

297 """Already-compressed FITS HDUs that are attached to high-level objects 

298 via `FitsOpaqueMetadata`, allowing lossy-compressed pixel values to be 

299 round-tripped exactly. 

300 """ 

301 

302 header: astropy.io.fits.Header 

303 """Header for the HDU. 

304 

305 This contains only FITS tile-compression keywords. 

306 """ 

307 

308 data: astropy.io.fits.FITS_rec 

309 """FITS binary table data that serves as the low-level representation of a 

310 tile-compressed image HDU. 

311 """ 

312 

313 @classmethod 

314 def from_bintable(cls, hdu: astropy.io.fits.BinTableHDU) -> Self: 

315 """Construct from a binary table HDU. 

316 

317 Parameters 

318 ---------- 

319 hdu 

320 Binary table HDU, typically read from a FITS file opened with 

321 ``disable_image_compression=True``. 

322 

323 Returns 

324 ------- 

325 PrecompressedImage 

326 A `PrecompressedImage` instance. 

327 """ 

328 header = astropy.io.fits.Header( 

329 [ 

330 card 

331 for card in hdu.header.cards 

332 if card.keyword in _COMPRESSION_KEYS 

333 or any(card.keyword.startswith(k) for k in _COMPRESSION_PREFIX_KEYS) 

334 ] 

335 ) 

336 # This is an opportunity to fix CFITSIO's non-standard writing of the 

337 # old RICE_ONE value instead of RICE_1. 

338 if header["ZCMPTYPE"] == "RICE_ONE": 

339 header["ZCMPTYPE"] = "RICE_1" 

340 return cls(header=header, data=hdu.data) 

341 

342 

343@final 

344class FitsExternalMetadata(ExternalMetadata): 

345 """Read-only, by-keyword access to the values in a FITS header. 

346 

347 Parameters 

348 ---------- 

349 header 

350 Header to read. It is held by reference, so later changes to it are 

351 seen. `None` is equivalent to an empty header. 

352 

353 Notes 

354 ----- 

355 ``COMMENT``, ``HISTORY``, and blank-keyword cards are not visible. 

356 Keywords are reported in upper case, with HIERARCH keywords lacking the 

357 ``HIERARCH`` prefix (e.g. ``LSST ISR UNITS``). Lookups ignore case and 

358 accept an optional ``HIERARCH`` prefix. Cards with no value are 

359 reported as `None`. 

360 """ 

361 

362 _HIDDEN_KEYWORDS: ClassVar[frozenset[str]] = frozenset({"", "COMMENT", "HISTORY"}) 

363 """Keywords of FITS commentary cards, which are not visible.""" 

364 

365 def __init__(self, header: astropy.io.fits.Header | None) -> None: 

366 self._header = header if header is not None else astropy.io.fits.Header() 

367 

368 @staticmethod 

369 def _from_card_value(value: Any) -> ExternalMetadataValue: 

370 """Convert a FITS card value to an `ExternalMetadataValue`, mapping 

371 the value of a card with no value to `None`. 

372 """ 

373 if isinstance(value, astropy.io.fits.card.Undefined): 373 ↛ 374line 373 didn't jump to line 374 because the condition on line 373 was never true

374 return None 

375 return value 

376 

377 def _is_visible(self, key: object) -> bool: 

378 """Return whether the header has a visible card for ``key``.""" 

379 return ( 

380 isinstance(key, str) 

381 and key.upper().removeprefix("HIERARCH ") not in self._HIDDEN_KEYWORDS 

382 and key in self._header 

383 ) 

384 

385 def __getitem__(self, key: str) -> ExternalMetadataValue: 

386 if not self._is_visible(key): 

387 raise KeyError(key) 

388 return self._from_card_value(self._header[key]) 

389 

390 def __contains__(self, key: object) -> bool: 

391 return self._is_visible(key) 

392 

393 def __iter__(self) -> Iterator[str]: 

394 return iter( 

395 dict.fromkeys( 

396 keyword 

397 for keyword in (k.upper() for k in self._header.keys()) 

398 if keyword not in self._HIDDEN_KEYWORDS 

399 ) 

400 ) 

401 

402 def __len__(self) -> int: 

403 return sum(1 for _ in self) 

404 

405 def get_all(self, key: str) -> tuple[ExternalMetadataValue, ...]: 

406 # Docstring inherited. 

407 if not self._is_visible(key): 

408 raise KeyError(key) 

409 return tuple(self._from_card_value(self._header[key, i]) for i in range(self._header.count(key))) 

410 

411 

412@final 

413@dataclasses.dataclass 

414class FitsOpaqueMetadata(OpaqueArchiveMetadata): 

415 """Opaque metadata that may be carried around by a serializable type to 

416 propagate serialization options and opaque information without that type 

417 knowing how it was serialized. 

418 """ 

419 

420 headers: dict[ExtensionKey, astropy.io.fits.Header] = dataclasses.field(default_factory=dict) 

421 """FITS headers found (but not interpreted and stripped) when reading, to 

422 be propagated on write. 

423 

424 Keys are EXTNAME/EXTVER combinations, or ("", 1) for the primary header. 

425 Header information in opaque metadata is considered immutable, allowing it 

426 to be transferred by reference to copies and subsets of the object it is 

427 attached to. 

428 """ 

429 

430 precompressed: dict[str, PrecompressedImage] = dataclasses.field(default_factory=dict) 

431 """FITS tile-compressed HDUs that should be written out directly instead 

432 of the in-memory data provided. 

433 

434 Keys are EXTNAME values, which must be unique (no EXTVER disambiguation). 

435 Precompressed pixel values are never copied or transferred to subsets. 

436 """ 

437 

438 def add_header( 

439 self, 

440 header: astropy.io.fits.Header, 

441 name: str | None = None, 

442 ver: int | None = None, 

443 key: ExtensionKey | None = None, 

444 ) -> None: 

445 """Add a header to the opaque metadata if it is not already present, 

446 and strip EXTNAME, EXTVER, and DATE if present. 

447 

448 Parameters 

449 ---------- 

450 header 

451 Header to add. May be modified in place. 

452 name 

453 EXTNAME (all caps). If not provided, the EXTNAME card must be 

454 present in the header. Use ``""`` for the primary header. 

455 ver 

456 EXTVER. If not provided and the EXTVER card is not present, 

457 defaults to 1. 

458 key 

459 Combination of EXTNAME and EXTVER; used instead of both ``name`` 

460 and ``ver`` if provided. 

461 """ 

462 if key is None: 

463 if name is None: 

464 name = header["EXTNAME"] 

465 if ver is None: 

466 ver = header.get("EXTVER", 1) 

467 key = ExtensionKey(name, ver) 

468 strip_butler_cards(header) 

469 if key not in self.headers: 469 ↛ exitline 469 didn't return from function 'add_header' because the condition on line 469 was always true

470 header.remove("EXTNAME", ignore_missing=True) 

471 header.remove("EXTVER", ignore_missing=True) 

472 # DATE is written fresh on every HDU, so it is regenerated on write 

473 # rather than propagated as opaque metadata. 

474 header.remove("DATE", ignore_missing=True) 

475 self.headers[key] = header 

476 

477 def maybe_use_precompressed(self, name: str) -> astropy.io.fits.BinTableHDU | None: 

478 """Look up the given EXTNAME to see if there is a tile compressed image 

479 HDU that should be used directly, instead of requantizing. 

480 

481 Parameters 

482 ---------- 

483 name 

484 EXTNAME (all caps). 

485 

486 Returns 

487 ------- 

488 `astropy.io.fits.BinTableHDU` | `None` 

489 An already-compressed HDU, in binary table form, or `None` if there 

490 is no precompressed HDU for this EXTNAME. 

491 """ 

492 if (precompressed := self.precompressed.get(name)) is None: 492 ↛ 494line 492 didn't jump to line 494 because the condition on line 492 was always true

493 return None 

494 return astropy.io.fits.BinTableHDU(precompressed.data, header=precompressed.header.copy(), name=name) 

495 

496 def extract_legacy_primary_header(self, header: astropy.io.fits.Header) -> dict[str, Any]: 

497 """Update the opaque metadata with the header of the primary HDU 

498 of a legacy (`lsst.afw.image`) FITS file, stripping cards we know we 

499 don't need and extracting any ``LSST IMAGES ...`` cards into a 

500 dictionary we return. 

501 

502 Parameters 

503 ---------- 

504 header 

505 Primary HDU header of the legacy FITS file. 

506 """ 

507 metadata: dict[str, Any] = {} 

508 for n in itertools.count(): 508 ↛ 513line 508 didn't jump to line 513 because the loop on line 508 didn't complete

509 if (key := header.pop(f"LSST IMAGES KEY {n + 1}", ...)) is ...: 

510 break 

511 value = header.pop(f"LSST IMAGES VALUE {n + 1}") 

512 metadata[key] = value 

513 primary_header = header.copy(strip=True) 

514 # No idea what these spare TAN-SIP headers are doing in the afw 

515 # FITS files, but we'll strip them here: 

516 primary_header.remove("A_ORDER", ignore_missing=True) 

517 primary_header.remove("B_ORDER", ignore_missing=True) 

518 primary_header.remove("DATE", ignore_missing=True) 

519 strip_legacy_exposure_cards(primary_header) 

520 strip_butler_cards(primary_header) 

521 self.headers[ExtensionKey()] = primary_header 

522 return metadata 

523 

524 def add_cutdown_primary_header(self, header: astropy.io.fits.Header) -> None: 

525 """Add the primary header of a cut-down ``lsst.images`` FITS file as 

526 opaque metadata. 

527 

528 A cut-down file is one whose JSON-tree, index, and nested-archive HDUs 

529 have been dropped, leaving only the primary, image, mask, and variance 

530 HDUs (e.g. as written by ``dax_images_cutout``). Only the 

531 container-layout cards that would be misleading after re-serialization 

532 are stripped; all other cards are retained so they survive the round 

533 trip. 

534 

535 Parameters 

536 ---------- 

537 header 

538 Primary header to read. Not modified. 

539 """ 

540 primary_header = header.copy(strip=True) 

541 for keyword in ("FMTVER", "INDXADDR", "INDXSIZE", "JSONADDR", "JSONSIZE", "DATAMODL"): 

542 primary_header.remove(keyword, ignore_missing=True) 

543 self.add_header(primary_header, name="", ver=1) 

544 

545 def copy(self) -> FitsOpaqueMetadata: 

546 # Docstring inherited. 

547 return FitsOpaqueMetadata(headers=self.headers) 

548 

549 def subset(self, bbox: Box) -> FitsOpaqueMetadata: 

550 # Docstring inherited. 

551 return FitsOpaqueMetadata(headers=self.headers) 

552 

553 def external_metadata(self) -> FitsExternalMetadata: 

554 # Docstring inherited. 

555 return FitsExternalMetadata(self.headers.get(ExtensionKey())) 

556 

557 def get_instrumental_unit(self) -> astropy.units.UnitBase | None: 

558 """Extract the ``LSST ISR UNIT`` key from the primary header (if it 

559 exists) and wrap it with Astropy. 

560 """ 

561 if (primary_header := self.headers.get(ExtensionKey())) is not None: 561 ↛ 564line 561 didn't jump to line 564 because the condition on line 561 was always true

562 if (instrumental_unit_str := primary_header.get("LSST ISR UNITS")) is not None: 562 ↛ 563line 562 didn't jump to line 563 because the condition on line 562 was never true

563 return astropy.units.Unit(instrumental_unit_str) 

564 return None 

565 

566 

567def add_offset_wcs(header: astropy.io.fits.Header, *, x: int | float, y: int | float, key: str = "A") -> None: 

568 """Add a trivial FITS WCS to a header that applies the appropriate offset 

569 to map FITS array coordinates to a logical pixel grid. 

570 

571 Parameters 

572 ---------- 

573 header 

574 Header to update in-place. 

575 x 

576 Logical coordinate of the first column. 

577 y 

578 Logical coordinate of the first row. 

579 key 

580 Single-character suffix for this WCS. 

581 """ 

582 header.set(f"CTYPE1{key}", "LINEAR") 

583 header.set(f"CTYPE2{key}", "LINEAR") 

584 header.set(f"CRPIX1{key}", 1.0) 

585 header.set(f"CRPIX2{key}", 1.0) 

586 header.set(f"CRVAL1{key}", float(x)) 

587 header.set(f"CRVAL2{key}", float(y)) 

588 header.set(f"CUNIT1{key}", "PIXEL") 

589 header.set(f"CUNIT2{key}", "PIXEL") 

590 

591 

592def read_offset_wcs(header: astropy.io.fits.Header, *, key: str = "A") -> tuple[int, int] | None: 

593 """Recover the logical pixel origin from a linear offset WCS written by 

594 `add_offset_wcs`. 

595 

596 Parameters 

597 ---------- 

598 header 

599 Header to read from. 

600 key 

601 Single-character suffix for the WCS to read. 

602 

603 Returns 

604 ------- 

605 `tuple` [`int`, `int`] | `None` 

606 The ``(x, y)`` logical coordinates of the first column and row, or 

607 `None` if the header has no ``LINEAR`` WCS with this suffix. 

608 """ 

609 if header.get(f"CTYPE1{key}") != "LINEAR": 

610 return None 

611 x = header[f"CRVAL1{key}"] - (header.get(f"CRPIX1{key}", 1.0) - 1.0) 

612 y = header[f"CRVAL2{key}"] - (header.get(f"CRPIX2{key}", 1.0) - 1.0) 

613 return (round(x), round(y)) 

614 

615 

616def read_yx0(header: astropy.io.fits.Header) -> YX[int]: 

617 """Recover the logical origin of an array from its FITS header. 

618 

619 Parameters 

620 ---------- 

621 header 

622 Header to read from. Not modified. 

623 

624 Returns 

625 ------- 

626 `~lsst.images.YX` [`int`] 

627 Logical coordinate of the first pixel, ordered ``(y, x)``. 

628 

629 Raises 

630 ------ 

631 ValueError 

632 Raised if the header records no origin via either the ``LTV1``/``LTV2`` 

633 cards (written by `lsst.afw.image`) or a linear offset WCS (written by 

634 `add_offset_wcs`). 

635 

636 Notes 

637 ----- 

638 The ``LTV1``/``LTV2`` convention is preferred when present, falling back to 

639 the ``A`` offset WCS otherwise. 

640 """ 

641 if "LTV1" in header: 

642 return YX(y=-round(header["LTV2"]), x=-round(header["LTV1"])) 

643 if (xy0 := read_offset_wcs(header)) is not None: 

644 return YX(y=xy0[1], x=xy0[0]) 

645 raise ValueError("Header records no LTV1/LTV2 cards or linear offset WCS.") 

646 

647 

648_WCS_VECTOR_KEYS = ("CUNIT", "CRPIX", "CRPIX", "CRVAL", "CRDELT", "CROTA", "CRDER", "CSYER", "CDELT") 

649_WCS_MATRIX_KEYS = ("CD{0}_{1}", "PC{0}_{1}") 

650 

651 

652def strip_wcs_cards(header: astropy.io.fits.Header) -> None: 

653 """Strip WCS cards from a FITS header. 

654 

655 This does *not* attempt to cover all possible FITS WCS forms; it focuses on 

656 the ones we actually plan to write (simple undistorted ones + TAN-SIP). 

657 

658 Parameters 

659 ---------- 

660 header 

661 FITS header to strip WCS cards from in place. 

662 """ 

663 wcsaxes = header.pop("WCSAXES", 2) 

664 for wcsname in [""] + list(string.ascii_uppercase): 

665 header.remove("RADESYS" + wcsname, ignore_missing=True) 

666 if "CTYPE1" + wcsname in header: 

667 ctype: str = "" # just for linters that can't figure out that the loop always executes 

668 for n in range(wcsaxes): 

669 suffix = f"{n + 1}{wcsname}" 

670 ctype = header.pop("CTYPE" + suffix) 

671 for key in _WCS_VECTOR_KEYS: 

672 header.remove(key + suffix, ignore_missing=True) 

673 for m in range(wcsaxes): 

674 for tmpl in _WCS_MATRIX_KEYS: 

675 header.remove(tmpl.format(m + 1, suffix), ignore_missing=True) 

676 if ctype.endswith("-SIP"): 676 ↛ 677line 676 didn't jump to line 677 because the condition on line 676 was never true

677 _strip_sip_poly(header, wcsname, "A") 

678 _strip_sip_poly(header, wcsname, "B") 

679 _strip_sip_poly(header, wcsname, "AP") 

680 _strip_sip_poly(header, wcsname, "BP") 

681 header.remove("LONPOLE", ignore_missing=True) 

682 header.remove("LATPOLE", ignore_missing=True) 

683 header.remove("MJDREF", ignore_missing=True) 

684 header.remove("PV1_3", ignore_missing=True) 

685 header.remove("PV1_4", ignore_missing=True) 

686 

687 

688def _strip_sip_poly(header: astropy.io.fits.Header, wcsname: str, which: str) -> None: 

689 order: int | None = header.pop(f"{which}_ORDER{wcsname}", None) 

690 if order is not None: 

691 for i, j in itertools.product(range(order + 1), range(order + 1)): 

692 header.remove(f"{which}_{i}_{j}{wcsname}", ignore_missing=True) 

693 

694 

695def strip_legacy_exposure_cards(header: astropy.io.fits.Header) -> None: 

696 """Strip header keywords added by lsst.afw.image.Exposure. 

697 

698 Parameters 

699 ---------- 

700 header 

701 FITS header to strip exposure cards from in place. 

702 """ 

703 header.remove("AR_HDU", ignore_missing=True) 

704 for name in ( 

705 "FILTER", 

706 "DETECTOR", 

707 "VALID_POLYGON", 

708 "SKYWCS", 

709 "PSF", 

710 "SUMMARYSTATS", 

711 "AP_CORR_MAP", 

712 "PHOTOCALIB", 

713 ): 

714 header.remove(f"{name}_ID", ignore_missing=True) 

715 header.remove(f"ARCHIVE_ID_{name}", ignore_missing=True) 

716 

717 

718def strip_butler_cards(header: astropy.io.fits.Header) -> None: 

719 """Strip header keywords added by butler provenance that would be 

720 incorrect if propagated to a downstream file. 

721 

722 Parameters 

723 ---------- 

724 header 

725 FITS header to strip butler provenance cards from in place. 

726 """ 

727 for key in list(header): 

728 if key.startswith("LSST BUTLER"): 

729 del header[key] 

730 

731 

732def parse_legacy_bunit( 

733 bunit: str, instrumental_unit: astropy.units.UnitBase | None = None 

734) -> astropy.units.UnitBase: 

735 """Handle a BUNIT value written by lsst.afw or LSST-supported instruments. 

736 

737 Parameters 

738 ---------- 

739 bunit 

740 BUNIT FITS header card value. 

741 instrumental_unit 

742 The correct instrumental unit for this instrument's post-ISR images. 

743 Should only be passed if this *is* a post-ISR image. 

744 

745 Notes 

746 ----- 

747 This works around two problems with the FITS BUNIT headers: 

748 

749 - At least according to Astropy, BUNIT='ADU' is not valid but 'adu' is. 

750 - Some files are written with BUNIT='adu' when the LSSTCam ISR's output 

751 units are 'electron'. 

752 """ 

753 bunit = bunit.replace("ADU", "adu") 

754 try: 

755 unit = astropy.units.Unit(bunit, format="fits") 

756 except ValueError: 

757 # Accept non-FITS units by assuming Astropy can still figure 

758 # them out if we don't specify the format. 

759 unit = astropy.units.Unit(bunit) 

760 # Fix incorrect BUNIT='adu' in LSST preliminary_visit_image, accounting 

761 # for variance planes with squared units as well. 

762 if unit == astropy.units.adu and instrumental_unit == astropy.units.electron: 

763 unit = astropy.units.electron 

764 if unit == astropy.units.adu**2 and instrumental_unit == astropy.units.electron: 764 ↛ 765line 764 didn't jump to line 765 because the condition on line 764 was never true

765 unit = astropy.units.electron**2 

766 return unit