Coverage for python/astro_metadata_translator/translators/decam.py: 32%
162 statements
« prev ^ index » next coverage.py v7.10.6, created at 2025-10-09 03:26 -0700
1# This file is part of astro_metadata_translator.
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 LICENSE 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.
12"""Metadata translation code for DECam FITS headers."""
14from __future__ import annotations
16__all__ = ("DecamTranslator",)
18import logging
19import posixpath
20import re
21from collections.abc import Iterator, Mapping, MutableMapping
22from typing import TYPE_CHECKING, Any
24import astropy.time
25import astropy.units as u
26from astropy.coordinates import Angle, EarthLocation, UnknownSiteException
27from astropy.io import fits
29from ..translator import CORRECTIONS_RESOURCE_ROOT, cache_translation
30from .fits import FitsTranslator
31from .helpers import altaz_from_degree_headers, is_non_science, tracking_from_degree_headers
33if TYPE_CHECKING:
34 import astropy.coordinates
36log = logging.getLogger(__name__)
38_CTIO_FALLBACK_LOCATION = EarthLocation.from_geocentric(
39 1814303.74553723, -5214365.7436216, -3187340.56598756, unit=u.m
40)
43class DecamTranslator(FitsTranslator):
44 """Metadata translator for DECam standard headers."""
46 name = "DECam"
47 """Name of this translation class"""
49 supported_instrument = "DECam"
50 """Supports the DECam instrument."""
52 default_resource_root = posixpath.join(CORRECTIONS_RESOURCE_ROOT, "DECam")
53 """Default resource path root to use to locate header correction files."""
55 # DECam has no rotator, and the instrument angle on sky is set to +Y=East,
56 # +X=South which we define as a 90 degree rotation and an X-flip.
57 _const_map = {
58 "boresight_rotation_angle": Angle(90 * u.deg),
59 "boresight_rotation_coord": "sky",
60 }
62 _trivial_map: dict[str, str | list[str] | tuple[Any, ...]] = {
63 "exposure_time": ("EXPTIME", {"unit": u.s}),
64 "exposure_time_requested": ("EXPREQ", {"unit": u.s}),
65 "dark_time": ("DARKTIME", {"unit": u.s}),
66 "boresight_airmass": ("AIRMASS", {"checker": is_non_science}),
67 "observation_id": "OBSID",
68 "object": "OBJECT",
69 "science_program": "PROPID",
70 "detector_num": "CCDNUM",
71 "detector_serial": "DETECTOR",
72 "detector_unique_name": "DETPOS",
73 "telescope": ("TELESCOP", {"default": "CTIO 4.0-m telescope"}),
74 "instrument": ("INSTRUME", {"default": "DECam"}),
75 # Ensure that reasonable values are always available
76 "relative_humidity": ("HUMIDITY", {"default": 40.0, "minimum": 0, "maximum": 100.0}),
77 "temperature": ("OUTTEMP", {"unit": u.deg_C, "default": 10.0, "minimum": -10.0, "maximum": 40.0}),
78 # Header says torr but seems to be mbar. Use hPa unit
79 # which is the SI equivalent of mbar.
80 "pressure": ("PRESSURE", {"unit": u.hPa, "default": 771.611, "minimum": 700.0, "maximum": 850.0}),
81 }
83 # DECam has no formal concept of an observing day. To ensure that
84 # observations from a single night all have the same observing_day, adopt
85 # the same offset used by the Vera Rubin Observatory of 12 hours.
86 _observing_day_offset = astropy.time.TimeDelta(12 * 3600, format="sec", scale="tai")
88 # List from Frank Valdes (2024-03-21).
89 _sky_observation_types: tuple[str, ...] = ("science", "object", "standard", "sky flat")
90 _non_sky_observation_types: tuple[str, ...] = ("zero", "dark", "dome flat")
92 # Unique detector names are currently not used but are read directly from
93 # header.
94 # The detector_group could be N or S with detector_name corresponding
95 # to the number in that group.
96 detector_names = {
97 1: "S29",
98 2: "S30",
99 3: "S31",
100 4: "S25",
101 5: "S26",
102 6: "S27",
103 7: "S28",
104 8: "S20",
105 9: "S21",
106 10: "S22",
107 11: "S23",
108 12: "S24",
109 13: "S14",
110 14: "S15",
111 15: "S16",
112 16: "S17",
113 17: "S18",
114 18: "S19",
115 19: "S8",
116 20: "S9",
117 21: "S10",
118 22: "S11",
119 23: "S12",
120 24: "S13",
121 25: "S1",
122 26: "S2",
123 27: "S3",
124 28: "S4",
125 29: "S5",
126 30: "S6",
127 31: "S7",
128 32: "N1",
129 33: "N2",
130 34: "N3",
131 35: "N4",
132 36: "N5",
133 37: "N6",
134 38: "N7",
135 39: "N8",
136 40: "N9",
137 41: "N10",
138 42: "N11",
139 43: "N12",
140 44: "N13",
141 45: "N14",
142 46: "N15",
143 47: "N16",
144 48: "N17",
145 49: "N18",
146 50: "N19",
147 51: "N20",
148 52: "N21",
149 53: "N22",
150 54: "N23",
151 55: "N24",
152 56: "N25",
153 57: "N26",
154 58: "N27",
155 59: "N28",
156 60: "N29",
157 62: "N31",
158 }
160 @classmethod
161 def can_translate(cls, header: Mapping[str, Any], filename: str | None = None) -> bool:
162 """Indicate whether this translation class can translate the
163 supplied header.
165 Checks the INSTRUME and FILTER headers.
167 Parameters
168 ----------
169 header : `dict`-like
170 Header to convert to standardized form.
171 filename : `str`, optional
172 Name of file being translated.
174 Returns
175 -------
176 can : `bool`
177 `True` if the header is recognized by this class. `False`
178 otherwise.
179 """
180 # Use INSTRUME. Because of defaulting behavior only do this
181 # if we really have an INSTRUME header
182 if "INSTRUME" in header:
183 via_instrume = super().can_translate(header, filename=filename)
184 if via_instrume:
185 return via_instrume
186 if cls.is_keyword_defined(header, "FILTER") and "DECam" in header["FILTER"]:
187 return True
188 return False
190 @cache_translation
191 def to_exposure_id(self) -> int:
192 """Calculate exposure ID.
194 Returns
195 -------
196 id : `int`
197 ID of exposure.
198 """
199 value = self._header["EXPNUM"]
200 self._used_these_cards("EXPNUM")
201 return value
203 @cache_translation
204 def to_observation_counter(self) -> int:
205 """Return the lifetime exposure number.
207 Returns
208 -------
209 sequence : `int`
210 The observation counter.
211 """
212 return self.to_exposure_id()
214 @cache_translation
215 def to_visit_id(self) -> int:
216 # Docstring will be inherited. Property defined in properties.py
217 return self.to_exposure_id()
219 @cache_translation
220 def to_datetime_end(self) -> astropy.time.Time:
221 # Docstring will be inherited. Property defined in properties.py
222 # Instcals have no DATE-END or DTUTC
223 datetime_end = self._from_fits_date("DTUTC", scale="utc")
224 if datetime_end is None:
225 datetime_end = self.to_datetime_begin() + self.to_exposure_time()
226 return datetime_end
228 def _translate_from_calib_id(self, field: str) -> str:
229 """Fetch the ID from the CALIB_ID header.
231 Calibration products made with constructCalibs have some metadata
232 saved in its FITS header CALIB_ID.
234 Parameters
235 ----------
236 field : `str`
237 Field to extract from the ``CALIB_ID`` header.
239 Returns
240 -------
241 value : `str`
242 The value extracted from the calibration header for that field.
243 """
244 data = self._header["CALIB_ID"]
245 match = re.search(rf".*{field}=(\S+)", data)
246 if not match:
247 raise RuntimeError(f"Header CALIB_ID with value '{data}' has not field '{field}'")
248 self._used_these_cards("CALIB_ID")
249 return match.groups()[0]
251 @cache_translation
252 def to_physical_filter(self) -> str | None:
253 """Calculate physical filter.
255 Return `None` if the keyword FILTER does not exist in the header,
256 which can happen for some valid Community Pipeline products.
258 Returns
259 -------
260 filter : `str`
261 The full filter name.
262 """
263 if self.is_key_ok("FILTER"):
264 value = self._header["FILTER"].strip()
265 self._used_these_cards("FILTER")
266 return value
267 elif self.is_key_ok("CALIB_ID"):
268 return self._translate_from_calib_id("filter")
269 else:
270 return None
272 @cache_translation
273 def to_location(self) -> astropy.coordinates.EarthLocation:
274 """Calculate the observatory location.
276 Returns
277 -------
278 location : `astropy.coordinates.EarthLocation`
279 An object representing the location of the telescope.
280 """
281 if self.is_key_ok("OBS-LONG"):
282 # OBS-LONG has west-positive sign so must be flipped
283 lon = self._header["OBS-LONG"] * -1.0
284 value = EarthLocation.from_geodetic(lon, self._header["OBS-LAT"], self._header["OBS-ELEV"])
285 self._used_these_cards("OBS-LONG", "OBS-LAT", "OBS-ELEV")
286 else:
287 # Look up the value since some files do not have location
288 try:
289 value = EarthLocation.of_site("ctio")
290 except UnknownSiteException:
291 value = _CTIO_FALLBACK_LOCATION
293 return value
295 @cache_translation
296 def to_observation_type(self) -> str:
297 """Calculate the observation type.
299 Returns
300 -------
301 typ : `str`
302 Observation type. Normalized to standard set.
303 """
304 if not self.is_key_ok("OBSTYPE"):
305 return "none"
306 obstype = self._header["OBSTYPE"].strip().lower()
307 self._used_these_cards("OBSTYPE")
308 if obstype == "object":
309 return "science"
310 return obstype
312 @cache_translation
313 def to_tracking_radec(self) -> astropy.coordinates.SkyCoord:
314 # Docstring will be inherited. Property defined in properties.py
315 radecsys = ("RADESYS",)
316 radecpairs = (("TELRA", "TELDEC"),)
317 return tracking_from_degree_headers(self, radecsys, radecpairs, unit=(u.hourangle, u.deg))
319 @cache_translation
320 def to_altaz_begin(self) -> astropy.coordinates.AltAz:
321 # Docstring will be inherited. Property defined in properties.py
322 return altaz_from_degree_headers(self, (("ZD", "AZ"),), self.to_datetime_begin(), is_zd={"ZD"})
324 @cache_translation
325 def to_detector_exposure_id(self) -> int | None:
326 # Docstring will be inherited. Property defined in properties.py
327 exposure_id = self.to_exposure_id()
328 if exposure_id is None:
329 return None
330 return int(f"{exposure_id:07d}{self.to_detector_num():02d}")
332 @cache_translation
333 def to_detector_group(self) -> str:
334 # Docstring will be inherited. Property defined in properties.py
335 name = self.to_detector_unique_name()
336 return name[0]
338 @cache_translation
339 def to_detector_name(self) -> str:
340 # Docstring will be inherited. Property defined in properties.py
341 name = self.to_detector_unique_name()
342 return name[1:]
344 @cache_translation
345 def to_focus_z(self) -> u.Quantity:
346 # Docstring will be inherited. Property defined in properties.py
347 # ``TELFOCUS`` is a comma-separated string with six focus offsets
348 # (fx, fy, fz, tx, ty, tz) recorded in units of microns.
349 tel_focus_list = self._header["TELFOCUS"].split(",")
350 return float(tel_focus_list[2]) * u.um
352 @classmethod
353 def fix_header(
354 cls, header: MutableMapping[str, Any], instrument: str, obsid: str, filename: str | None = None
355 ) -> bool:
356 """Fix DECam headers.
358 Parameters
359 ----------
360 header : `dict`
361 The header to update. Updates are in place.
362 instrument : `str`
363 The name of the instrument.
364 obsid : `str`
365 Unique observation identifier associated with this header.
366 Will always be provided.
367 filename : `str`, optional
368 Filename associated with this header. May not be set since headers
369 can be fixed independently of any filename being known.
371 Returns
372 -------
373 modified = `bool`
374 Returns `True` if the header was updated.
376 Notes
377 -----
378 Fixes the following issues:
380 * If OBSTYPE contains "zero" or "bias",
381 update the FILTER keyword to "solid plate 0.0 0.0".
383 Corrections are reported as debug level log messages.
384 """
385 modified = False
387 # Calculate the standard label to use for log messages
388 log_label = cls._construct_log_prefix(obsid, filename)
390 obstype = header.get("OBSTYPE", "unknown")
392 if "bias" in obstype.lower() or "zero" in obstype.lower():
393 header["FILTER"] = "solid plate 0.0 0.0"
394 modified = True
395 log.debug("%s: Set FILTER to %s because OBSTYPE is %s", log_label, header["FILTER"], obstype)
397 # The following provides corrections to the boresight header info for
398 # the raw data files obtained from http://astroarchive.noirlab.edu/,
399 # program number 2013A-0351, of the Trifid and Lagoon Nebulae region.
400 if "Trifid" in header.get("OBJECT", "unknown") and obstype == "object":
401 # Create a translator since we need the date
402 translator = cls(header)
403 date_obs = translator.to_datetime_begin()
404 date_min = astropy.time.Time("2013-02-10", format="isot", scale="utc")
405 date_max = astropy.time.Time("2013-02-14", format="isot", scale="utc")
406 if date_obs > date_min and date_obs < date_max:
407 # Median differences (in deg) between the raw header
408 # RA/DEC values & the boresight values obtianed from
409 # https://astroarchive.noirlab.edu for these data.
410 delta_ra = -0.027417
411 delta_dec = -0.002833
412 ra_hour = Angle(header["RA"], unit="hour")
413 dec_deg = Angle(header["DEC"], unit="degree")
414 header["TELRA"] = (ra_hour + delta_ra * u.degree).hour
415 header["TELDEC"] = (dec_deg + delta_dec * u.degree).degree
416 modified = True
417 log.debug(
418 "Found Trifid in OBJECT header. Changing TELRA and TELDEC headers to "
419 "RA + %.5f deg and DEC + %.5f deg, respectively.",
420 delta_ra,
421 delta_dec,
422 )
424 return modified
426 @classmethod
427 def determine_translatable_headers(
428 cls, filename: str, primary: MutableMapping[str, Any] | None = None
429 ) -> Iterator[MutableMapping[str, Any]]:
430 """Given a file return all the headers usable for metadata translation.
432 DECam files are multi-extension FITS with a primary header and
433 each detector stored in a subsequent extension. DECam uses
434 ``INHERIT=T`` and each detector header will be merged with the
435 primary header.
437 Guide headers are not returned.
439 Parameters
440 ----------
441 filename : `str`
442 Path to a file in a format understood by this translator.
443 primary : `dict`-like, optional
444 The primary header obtained by the caller. This is sometimes
445 already known, for example if a system is trying to bootstrap
446 without already knowing what data is in the file. Will be
447 merged with detector headers if supplied, else will be read
448 from the file.
450 Yields
451 ------
452 headers : iterator of `dict`-like
453 Each detector header in turn. The supplied header will be merged
454 with the contents of each detector header.
456 Notes
457 -----
458 This translator class is specifically tailored to raw DECam data and
459 is not designed to work with general FITS files. The normal paradigm
460 is for the caller to have read the first header and then called
461 `determine_translator()` on the result to work out which translator
462 class to then call to obtain the real headers to be used for
463 translation.
464 """
465 # Circular dependency so must defer import.
466 from ..headers import merge_headers
468 # This is convoluted because we need to turn an Optional variable
469 # to a Dict so that mypy is happy.
470 primary_hdr = primary if primary else {}
472 # Since we want to scan many HDUs we use astropy directly to keep
473 # the file open rather than continually opening and closing it
474 # as we go to each HDU.
475 with fits.open(filename) as fits_file:
476 # Astropy does not automatically handle the INHERIT=T in
477 # DECam headers so the primary header must be merged.
478 first_pass = True
480 for hdu in fits_file:
481 if first_pass:
482 if not primary_hdr:
483 primary_hdr = hdu.header
484 first_pass = False
485 continue
487 header = hdu.header
488 if "CCDNUM" not in header: # Primary does not have CCDNUM
489 continue
490 if header["CCDNUM"] > 62: # ignore guide CCDs
491 continue
492 yield merge_headers([primary_hdr, header], mode="overwrite")
494 @classmethod
495 def observing_date_to_offset(cls, observing_date: astropy.time.Time) -> astropy.time.TimeDelta | None:
496 """Return the offset to use when calculating the observing day.
498 Parameters
499 ----------
500 observing_date : `astropy.time.Time`
501 The date of the observation. Unused.
503 Returns
504 -------
505 offset : `astropy.time.TimeDelta`
506 The offset to apply. The offset is always 12 hours. DECam has
507 no defined observing day concept in its headers. To ensure that
508 observations from a single night all have the same observing_day,
509 adopt the same offset used by the Vera Rubin Observatory of
510 12 hours.
511 """
512 return cls._observing_day_offset