Skip to content

Currently available interfaces to spectrographs

Data loading interfaces for data from different spectrographs. Each spectrograph is represented as a different class, with all providing the same user-interface (after they load all the necessary information from disk).

.. note:: During the years, the supported spectrographs have been (or will be) subjected to technical interventions, leading to changes in the instrumental profile. If SBART creates a stellar template that mixes information from before and after (such interventions), we will be adding noise into the model (as the spectral shape will be different)

To avoid such issues, we introduce the concept of "subInstruments", which define individual data-sets within a
spectrograph's lifespan. Thus, for each sub-Instrument, we create individual templates and extract independent
radial velocities. For more details on the dates used to delimit the data-sets, we refer to the docs of the individual
instruments.

instrument_dict = {'ESPRESSO': ESPRESSO, 'HARPS': HARPS, 'HARPSN': HARPSN, 'MAROONX': MAROONX, 'SimulatedSpirou': SimulatedSpirou, 'CARMENES': CARMENES, 'POET': POET} module-attribute

CARMENES

CARMENES

Bases: Frame

Interface to handle KOBE-CARMENES data (optical arm only).

Using this class implies passing an input file (shaq_output_folder), from where we will load the following information:

  • Column 0 -> BJD
  • Column 5 -> CCF RV in km/s
  • Column 3 -> CCF RV Error in km/s
  • Column 10 -> BERV in km/s (used if we configure the override_BERV option)
  • Column 7 -> Instrumental drift in m/s
  • Column 8 -> Error in Instrumental drift in m/s
  • Column 9 -> QC flag on instrumental drift

Any observation that can't be cross-matched to the BJD will be marked as invalid

Source code in src/ASTRA/Instruments/CARMENES.py 
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
class CARMENES(Frame):
    """Interface to handle KOBE-CARMENES data (optical arm only).

    Using this class implies passing an input file (shaq_output_folder), from where we will load
    the following information:

    - Column 0 -> BJD
    - Column 5 -> CCF RV in km/s
    - Column 3 -> CCF RV Error in km/s
    - Column 10 -> BERV in km/s (used if we configure the override_BERV option)
    - Column 7 -> Instrumental drift in m/s
    - Column 8 -> Error in Instrumental drift in m/s
    - Column 9 -> QC flag on instrumental drift

    Any observation that can't be cross-matched to the BJD will be marked as invalid
    """

    _default_params = Frame._default_params + DefaultValues(
        shaq_output_folder=UserParam(None, constraint=PathValue, mandatory=True),
        override_BERV=UserParam(True, constraint=BooleanValue, mandatory=False),
        is_KOBE_data=UserParam(True, constraint=BooleanValue, mandatory=False),
        max_hours_to_calibration=UserParam(
            100,
            constraint=NumericValue,
            mandatory=False,
            description="Maximum number of hours between observation and calibration. If exceeded, frame is invalid",
        ),
        sigma_clip_flux=UserParam(
            -1,
            constraint=NumericValue,  # -1 means that it is disabled,
            mandatory=False,
            description="If positive, then sigma clip the flux",
        ),
    )
    _default_params.update(
        "IS_SA_CORRECTED", UserParam(True, constraint=BooleanValue, mandatory=False)
    )

    sub_instruments = {
        "CARMENES": datetime.datetime.max,
    }
    _name = "CARMENES"

    KW_map = {
        "OBJECT": "OBJECT",
        "BJD": "HIERARCH CARACAL BJD",
        "MJD": "MJD-OBS",
        "ISO-DATE": "DATE-OBS",  # TODO: to check this KW name
        "DRS-VERSION": "HIERARCH CARACAL FOX VERSION",
        "RA": "RA",
        "DEC": "DEC",
    }

    def __init__(
        self,
        file_path,
        user_configs: Optional[Dict[str, Any]] = None,
        reject_subInstruments=None,
        frameID=None,
        quiet_user_params: bool = True,
    ):
        """

        Parameters
        ----------
        file_path
            Path to the S2D (or S1D) file.
        user_configs
            Dictionary whose keys are the configurable options of ESPRESSO (check above)
        reject_subInstruments
            Iterable of subInstruments to fully reject
        frameID
            ID for this observation. Only used for organization purposes by :class:`~SBART.data_objects.DataClass`
        """

        self._blaze_corrected = True

        super().__init__(
            inst_name=self._name,
            array_size={"S2D": [61, 4096]},
            file_path=file_path,
            frameID=frameID,
            KW_map=self.KW_map,
            available_indicators=("FWHM", "BIS SPAN"),
            user_configs=user_configs,
            reject_subInstruments=reject_subInstruments,
            need_external_data_load=True,
            quiet_user_params=quiet_user_params,
        )
        coverage = [500, 1750]
        self.instrument_properties["wavelength_coverage"] = coverage
        self.instrument_properties["is_drift_corrected"] = False

        # TODO: ensure that there are no problem when using the NIR arm of CARMENES!!!!!
        self.instrument_properties["resolution"] = 94600

        # lat/lon from: https://geohack.toolforge.org/geohack.php?pagename=Calar_Alto_Observatory&params=37_13_25_N_2_32_46_W_type:landmark_region:ES
        # height from: https://en.wikipedia.org/wiki/Calar_Alto_Observatory
        lat, lon = 37.223611, -2.546111
        self.instrument_properties["EarthLocation"] = EarthLocation.from_geodetic(
            lat=lat, lon=lon, height=2168
        )

        # from https://www.mide.com/air-pressure-at-altitude-calculator
        # and convert to Pa to mbar
        self.instrument_properties["site_pressure"] = 778.5095

        self.is_BERV_corrected = False

    def get_spectral_type(self):
        name_lowercase = self.file_path.stem
        if "vis_A" in name_lowercase:
            return "S2D"
        else:
            raise custom_exceptions.InternalError(
                f"{self.name} can't recognize the file that it received ( - {self.file_path.stem})!"
            )

    def load_instrument_specific_KWs(self, header):
        self.observation_info["airmass"] = header[f"AIRMASS"]

        # Load BERV info + previous RV
        self.observation_info["MAX_BERV"] = 30 * kilometer_second
        self.observation_info["BERV"] = (
            header["HIERARCH CARACAL BERV "] * kilometer_second
        )

        # TODO: check ambient temperature on CARMENES data TO SEE IF IT IS THE "REAL ONE"
        # Environmental KWs for telfit (also needs airmassm previously loaded)
        ambi_KWs = {
            "relative_humidity": "AMBI RHUM",
            "ambient_temperature": "AMBI TEMPERATURE",
        }

        for name, endKW in ambi_KWs.items():
            self.observation_info[name] = header[f"HIERARCH CAHA GEN {endKW}"]
            if "temperature" in name:  # store temperature in KELVIN for TELFIT
                self.observation_info[name] = convert_temperature(
                    self.observation_info[name], old_scale="Celsius", new_scale="Kelvin"
                )
        for order in range(self.N_orders):
            self.observation_info["orderwise_SNRs"].append(
                header[f"HIERARCH CARACAL FOX SNR {order}"]
            )

        try:
            self.observation_info["MOON PHASE"] = header[
                "HIERARCH CAHA INS SCHEDULER MOON PHASE"
            ]
            self.observation_info["MOON DISTANCE"] = header[
                "HIERARCH CAHA INS SCHEDULER MOON DISTANCE"
            ]
        except KeyError:
            self.observation_info["MOON PHASE"] = 0
            self.observation_info["MOON DISTANCE"] = 0

    def check_header_QC(self, header: fits.Header):
        """Header QC checks for CARMENES.

        1) Drift calibration was done (CARACAL DRIFT FP REF exists)
        2) Time between observation and calibration is smaller than "max_hours_to_calibration"

        Can add the following status:

        - KW_WARNING("Drift flag of KOBE is greater than 1")
            If th drift value is greater than one. This is actually set in the DataClass.load_CARMENES_extra_information()

        Args:
            header (fits.Header): _description_
        """
        kill_messages = []

        try:
            drift_file = header["CARACAL DRIFT FP REF"]

            night_drift = drift_file.split("s-")[0].split("car-")[-1]
            night_drift = datetime.datetime.strptime(night_drift, r"%Y%m%dT%Hh%Mm%S")
            night_data = datetime.datetime.strptime(
                header["CARACAL DATE-OBS"], r"%Y-%m-%dT%H:%M:%S"
            )
            msg = f"\tDistance between calibration (FP) and OB has surpassed the limit of {self._internal_configs['max_hours_to_calibration']}"
            if (
                abs(night_drift - night_data).total_seconds()
                > datetime.timedelta(
                    hours=self._internal_configs["max_hours_to_calibration"]
                ).total_seconds()
            ):
                kill_messages.append(msg)

        except KeyError:
            kill_messages.append("Observation missing the CARACAL DRIFT FP REF keyword")

        if self._internal_configs["is_KOBE_data"]:
            # we can have skysubs outside of KOBE
            for msg in kill_messages:
                logger.critical(msg)
                self.add_to_status(FATAL_KW(msg.replace("\t", "")))

    def load_S2D_data(self):
        if self.is_open:
            return
        super().load_S2D_data()

        with fits.open(self.file_path) as hdulist:
            s2d_data = hdulist["SPEC"].data * 100000  # spetra from all olders
            err_data = hdulist["SIG"].data * 100000
            wavelengths = hdulist["WAVE"].data  # vacuum wavelengths; no BERV correction

        self.wavelengths = wavelengths
        self.spectra = s2d_data
        self.uncertainties = err_data
        self.build_mask(bypass_QualCheck=True)
        return 1

    def finalize_data_load(self, bad_flag: Optional[Flag] = None):
        bad_flag = MISSING_SHAQ_RVS
        super().finalize_data_load(bad_flag)

        moon_sep = self.observation_info["MOON DISTANCE"]
        moon_illum = self.observation_info["MOON PHASE"]
        rv = self.observation_info["DRS_RV"]
        berv = self.observation_info["BERV"]
        fwhm = self.observation_info["FWHM"]

        curr_rv_diff = (rv - berv).to(kilometer_second).value
        conditions = (
            moon_sep > 80,
            (30 < moon_sep < 80) and moon_illum < 0.6,
            abs(curr_rv_diff) > 5 * fwhm,
        )
        if not any(conditions):
            logger.warning("Target is close to the moon, setting warning flag")
            self.add_to_status(KW_WARNING("Target is close to moon"))

    def load_S1D_data(self) -> Mask:
        raise NotImplementedError

    def build_mask(self, bypass_QualCheck: bool = False):
        # We evaluate the bad orders all at once
        super().build_mask(bypass_QualCheck, assess_bad_orders=False)

        bpmap0 = np.zeros((61, 4096), dtype=np.uint64)
        bpmap0[
            14:38, [2453 - 3, 2453 - 2, 2453 - 1, 2453, 2453 + 1, 2453 + 2, 2453 + 3]
        ] |= 1
        bpmap0[14:38, 1643] |= 1  # ghost of hotspot tail
        bpmap0[14:38, 2459] |= (
            1  # spikes of hotspot satellite (bug not correct due to bug in v2.00)
        )
        bpmap0[15:41, 3374] |= 1  # displaced column; ignore by marking as nan
        bpmap0[28, 3395:3400] |= 1  # car-20160701T00h49m36s-sci-gtoc-vis.fits
        bpmap0[34, 838:850] |= 1  # car-20160803T22h46m41s-sci-gtoc-vis.fits
        bpmap0[34, 2035:2044] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
        bpmap0[34, 3150:3161] |= 1  # car-20160803T22h46m41s-sci-gtoc-vis.fits
        bpmap0[35, 403:410] |= 1  # car-20160803T22h46m41s-sci-gtoc-vis
        bpmap0[35, 754:759] |= 1  # car-20170419T03h27m48s-sci-gtoc-vis
        bpmap0[35, 1083:1093] |= 1  # car-20160803T22h46m41s-sci-gtoc-vis
        bpmap0[35, 1944:1956] |= 1  # car-20160803T22h46m41s-sci-gtoc-vis
        bpmap0[35, 2710:2715] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
        bpmap0[35, 3050:3070] |= 1  # car-20160803T22h46m41s-sci-gtoc-vis
        bpmap0[35, 3706:3717] |= 1  # car-20160803T22h46m41s-sci-gtoc-vis
        bpmap0[35, 3706:3717] |= 1  # car-20160803T22h46m41s-sci-gtoc-vis
        bpmap0[36, 303:308] |= 1  # car-20170419T03h27m48s-sci-gtoc-vis
        bpmap0[36, 312:317] |= 1  # car-20170419T03h27m48s-sci-gtoc-vis
        bpmap0[36, 1311:1315] |= 1  # car-20170419T03h27m48s-sci-gtoc-vis
        bpmap0[36, 1325:1329] |= 1  # car-20170419T03h27m48s-sci-gtoc-vis
        bpmap0[37, 1326:1343] |= 1  # car-20170419T03h27m48s-sci-gtoc-vis
        bpmap0[39, 1076:1082] |= 1  # car-20170626T02h00m17s-sci-gtoc-vis
        bpmap0[39, 1204:1212] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
        bpmap0[39, 1236:1243] |= 1  # car-20170419T03h27m48s-sci-gtoc-vis
        bpmap0[39, 1463:1468] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
        bpmap0[39, 2196:2203] |= 1  # car-20160520T03h10m13s-sci-gtoc-vis.fits
        bpmap0[39, 2493:2504] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
        bpmap0[39, 3705:3717] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
        bpmap0[40, 2765:2773] |= 1  # car-20170419T03h27m48s-sci-gtoc-vis
        bpmap0[40, 3146:3153] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
        bpmap0[40, 3556:3564] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
        bpmap0[41, 486:491] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
        bpmap0[41, 495:501] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
        bpmap0[41, 1305:1315] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
        bpmap0[42, 480:490] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
        bpmap0[42, 1316:1330] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
        bpmap0[42, 2363:2368] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
        bpmap0[42, 2375:2382] |= 1  # car-20170509T03h05m21s-sci-gtoc-vis
        bpmap0[44, 3355:3361] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
        bpmap0[46, 311:321] |= 1  # car-20160701T00h49m36s-sci-gtoc-vis.fits
        bpmap0[46, 835:845] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
        bpmap0[46, 1156:1171] |= 1  # car-20160701T00h49m36s-sci-gtoc-vis.fits
        bpmap0[46, 1895:1905] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
        bpmap0[46, 2212:2232] |= 1  # car-20160701T00h49m36s-sci-gtoc-vis.fits
        bpmap0[47, 2127:2133] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
        bpmap0[47, 2218:2223] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
        bpmap0[47, 2260:2266] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
        bpmap0[47, 2313:2319] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
        bpmap0[47, 3111:3116] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
        bpmap0[47, 3267:3272] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
        bpmap0[47, 3316:3321] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
        bpmap0[47, 3432:3438] |= 1  # car-20170509T03h05m21s-sci-gtoc-vis
        bpmap0[47, 3480:3488] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
        bpmap0[47, 3658:3665] |= 1  # car-20170509T03h05m21s-sci-gtoc-vis
        bpmap0[49, 1008:1017] |= 1  # car-20160701T00h49m36s-sci-gtoc-vis.fits
        bpmap0[49, 2532:2544] |= 1  # car-20160701T00h49m36s-sci-gtoc-vis.fits
        bpmap0[49, 3046:3056] |= 1  # car-20160701T00h49m36s-sci-gtoc-vis.fits
        bpmap0[49, 3574:3588] |= 1  # car-20160701T00h49m36s-sci-gtoc-vis.fits

        # Constructing the bad pixel map, as defined by SERVAL in here:
        # https://github.com/mzechmeister/serval/blob/c2f47b26f1102333dfe76f93c2a686807cda02ce/src/inst_CARM_VIS.py#L95

        bpmap0[25:41, 3606:3588] |= 1  # Our rejection

        if self._internal_configs["sigma_clip_flux"] > 0:
            sigma = self._internal_configs["sigma_clip_flux"]
            for order_number in range(self.N_orders):
                cont = median_filter(self.spectra[order_number], size=500)
                inds = np.where(
                    self.spectra[order_number]
                    >= cont + sigma * self.uncertainties[order_number]
                )
                bpmap0[order_number, inds] |= 1
        self.spectral_mask.add_indexes_to_mask(np.where(bpmap0 != 0), QUAL_DATA)

        # remove extremely negative points!
        self.spectral_mask.add_indexes_to_mask(
            np.where(self.spectra < -3 * self.uncertainties), MISSING_DATA
        )
        self.spectral_mask.add_indexes_to_mask(
            np.where(self.uncertainties == 0), MISSING_DATA
        )

        self.assess_bad_orders()

    def close_arrays(self):
        super().close_arrays()
        self.is_BERV_corrected = False

sub_instruments = {'CARMENES': datetime.datetime.max} class-attribute instance-attribute

KW_map = {'OBJECT': 'OBJECT', 'BJD': 'HIERARCH CARACAL BJD', 'MJD': 'MJD-OBS', 'ISO-DATE': 'DATE-OBS', 'DRS-VERSION': 'HIERARCH CARACAL FOX VERSION', 'RA': 'RA', 'DEC': 'DEC'} class-attribute instance-attribute

is_BERV_corrected = False instance-attribute

get_spectral_type()

Source code in src/ASTRA/Instruments/CARMENES.py 
142
143
144
145
146
147
148
149
def get_spectral_type(self):
    name_lowercase = self.file_path.stem
    if "vis_A" in name_lowercase:
        return "S2D"
    else:
        raise custom_exceptions.InternalError(
            f"{self.name} can't recognize the file that it received ( - {self.file_path.stem})!"
        )

load_instrument_specific_KWs(header)

Source code in src/ASTRA/Instruments/CARMENES.py 
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
def load_instrument_specific_KWs(self, header):
    self.observation_info["airmass"] = header[f"AIRMASS"]

    # Load BERV info + previous RV
    self.observation_info["MAX_BERV"] = 30 * kilometer_second
    self.observation_info["BERV"] = (
        header["HIERARCH CARACAL BERV "] * kilometer_second
    )

    # TODO: check ambient temperature on CARMENES data TO SEE IF IT IS THE "REAL ONE"
    # Environmental KWs for telfit (also needs airmassm previously loaded)
    ambi_KWs = {
        "relative_humidity": "AMBI RHUM",
        "ambient_temperature": "AMBI TEMPERATURE",
    }

    for name, endKW in ambi_KWs.items():
        self.observation_info[name] = header[f"HIERARCH CAHA GEN {endKW}"]
        if "temperature" in name:  # store temperature in KELVIN for TELFIT
            self.observation_info[name] = convert_temperature(
                self.observation_info[name], old_scale="Celsius", new_scale="Kelvin"
            )
    for order in range(self.N_orders):
        self.observation_info["orderwise_SNRs"].append(
            header[f"HIERARCH CARACAL FOX SNR {order}"]
        )

    try:
        self.observation_info["MOON PHASE"] = header[
            "HIERARCH CAHA INS SCHEDULER MOON PHASE"
        ]
        self.observation_info["MOON DISTANCE"] = header[
            "HIERARCH CAHA INS SCHEDULER MOON DISTANCE"
        ]
    except KeyError:
        self.observation_info["MOON PHASE"] = 0
        self.observation_info["MOON DISTANCE"] = 0

check_header_QC(header: fits.Header)

Header QC checks for CARMENES.

1) Drift calibration was done (CARACAL DRIFT FP REF exists) 2) Time between observation and calibration is smaller than "max_hours_to_calibration"

Can add the following status:

  • KW_WARNING("Drift flag of KOBE is greater than 1") If th drift value is greater than one. This is actually set in the DataClass.load_CARMENES_extra_information()

Args: header (fits.Header): description

Source code in src/ASTRA/Instruments/CARMENES.py 
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
def check_header_QC(self, header: fits.Header):
    """Header QC checks for CARMENES.

    1) Drift calibration was done (CARACAL DRIFT FP REF exists)
    2) Time between observation and calibration is smaller than "max_hours_to_calibration"

    Can add the following status:

    - KW_WARNING("Drift flag of KOBE is greater than 1")
        If th drift value is greater than one. This is actually set in the DataClass.load_CARMENES_extra_information()

    Args:
        header (fits.Header): _description_
    """
    kill_messages = []

    try:
        drift_file = header["CARACAL DRIFT FP REF"]

        night_drift = drift_file.split("s-")[0].split("car-")[-1]
        night_drift = datetime.datetime.strptime(night_drift, r"%Y%m%dT%Hh%Mm%S")
        night_data = datetime.datetime.strptime(
            header["CARACAL DATE-OBS"], r"%Y-%m-%dT%H:%M:%S"
        )
        msg = f"\tDistance between calibration (FP) and OB has surpassed the limit of {self._internal_configs['max_hours_to_calibration']}"
        if (
            abs(night_drift - night_data).total_seconds()
            > datetime.timedelta(
                hours=self._internal_configs["max_hours_to_calibration"]
            ).total_seconds()
        ):
            kill_messages.append(msg)

    except KeyError:
        kill_messages.append("Observation missing the CARACAL DRIFT FP REF keyword")

    if self._internal_configs["is_KOBE_data"]:
        # we can have skysubs outside of KOBE
        for msg in kill_messages:
            logger.critical(msg)
            self.add_to_status(FATAL_KW(msg.replace("\t", "")))

load_S2D_data()

Source code in src/ASTRA/Instruments/CARMENES.py 
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
def load_S2D_data(self):
    if self.is_open:
        return
    super().load_S2D_data()

    with fits.open(self.file_path) as hdulist:
        s2d_data = hdulist["SPEC"].data * 100000  # spetra from all olders
        err_data = hdulist["SIG"].data * 100000
        wavelengths = hdulist["WAVE"].data  # vacuum wavelengths; no BERV correction

    self.wavelengths = wavelengths
    self.spectra = s2d_data
    self.uncertainties = err_data
    self.build_mask(bypass_QualCheck=True)
    return 1

finalize_data_load(bad_flag: Optional[Flag] = None)

Source code in src/ASTRA/Instruments/CARMENES.py 
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
def finalize_data_load(self, bad_flag: Optional[Flag] = None):
    bad_flag = MISSING_SHAQ_RVS
    super().finalize_data_load(bad_flag)

    moon_sep = self.observation_info["MOON DISTANCE"]
    moon_illum = self.observation_info["MOON PHASE"]
    rv = self.observation_info["DRS_RV"]
    berv = self.observation_info["BERV"]
    fwhm = self.observation_info["FWHM"]

    curr_rv_diff = (rv - berv).to(kilometer_second).value
    conditions = (
        moon_sep > 80,
        (30 < moon_sep < 80) and moon_illum < 0.6,
        abs(curr_rv_diff) > 5 * fwhm,
    )
    if not any(conditions):
        logger.warning("Target is close to the moon, setting warning flag")
        self.add_to_status(KW_WARNING("Target is close to moon"))

load_S1D_data() -> Mask

Source code in src/ASTRA/Instruments/CARMENES.py 
267
268
def load_S1D_data(self) -> Mask:
    raise NotImplementedError

build_mask(bypass_QualCheck: bool = False)

Source code in src/ASTRA/Instruments/CARMENES.py 
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
def build_mask(self, bypass_QualCheck: bool = False):
    # We evaluate the bad orders all at once
    super().build_mask(bypass_QualCheck, assess_bad_orders=False)

    bpmap0 = np.zeros((61, 4096), dtype=np.uint64)
    bpmap0[
        14:38, [2453 - 3, 2453 - 2, 2453 - 1, 2453, 2453 + 1, 2453 + 2, 2453 + 3]
    ] |= 1
    bpmap0[14:38, 1643] |= 1  # ghost of hotspot tail
    bpmap0[14:38, 2459] |= (
        1  # spikes of hotspot satellite (bug not correct due to bug in v2.00)
    )
    bpmap0[15:41, 3374] |= 1  # displaced column; ignore by marking as nan
    bpmap0[28, 3395:3400] |= 1  # car-20160701T00h49m36s-sci-gtoc-vis.fits
    bpmap0[34, 838:850] |= 1  # car-20160803T22h46m41s-sci-gtoc-vis.fits
    bpmap0[34, 2035:2044] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
    bpmap0[34, 3150:3161] |= 1  # car-20160803T22h46m41s-sci-gtoc-vis.fits
    bpmap0[35, 403:410] |= 1  # car-20160803T22h46m41s-sci-gtoc-vis
    bpmap0[35, 754:759] |= 1  # car-20170419T03h27m48s-sci-gtoc-vis
    bpmap0[35, 1083:1093] |= 1  # car-20160803T22h46m41s-sci-gtoc-vis
    bpmap0[35, 1944:1956] |= 1  # car-20160803T22h46m41s-sci-gtoc-vis
    bpmap0[35, 2710:2715] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
    bpmap0[35, 3050:3070] |= 1  # car-20160803T22h46m41s-sci-gtoc-vis
    bpmap0[35, 3706:3717] |= 1  # car-20160803T22h46m41s-sci-gtoc-vis
    bpmap0[35, 3706:3717] |= 1  # car-20160803T22h46m41s-sci-gtoc-vis
    bpmap0[36, 303:308] |= 1  # car-20170419T03h27m48s-sci-gtoc-vis
    bpmap0[36, 312:317] |= 1  # car-20170419T03h27m48s-sci-gtoc-vis
    bpmap0[36, 1311:1315] |= 1  # car-20170419T03h27m48s-sci-gtoc-vis
    bpmap0[36, 1325:1329] |= 1  # car-20170419T03h27m48s-sci-gtoc-vis
    bpmap0[37, 1326:1343] |= 1  # car-20170419T03h27m48s-sci-gtoc-vis
    bpmap0[39, 1076:1082] |= 1  # car-20170626T02h00m17s-sci-gtoc-vis
    bpmap0[39, 1204:1212] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
    bpmap0[39, 1236:1243] |= 1  # car-20170419T03h27m48s-sci-gtoc-vis
    bpmap0[39, 1463:1468] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
    bpmap0[39, 2196:2203] |= 1  # car-20160520T03h10m13s-sci-gtoc-vis.fits
    bpmap0[39, 2493:2504] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
    bpmap0[39, 3705:3717] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
    bpmap0[40, 2765:2773] |= 1  # car-20170419T03h27m48s-sci-gtoc-vis
    bpmap0[40, 3146:3153] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
    bpmap0[40, 3556:3564] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
    bpmap0[41, 486:491] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
    bpmap0[41, 495:501] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
    bpmap0[41, 1305:1315] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
    bpmap0[42, 480:490] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
    bpmap0[42, 1316:1330] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
    bpmap0[42, 2363:2368] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
    bpmap0[42, 2375:2382] |= 1  # car-20170509T03h05m21s-sci-gtoc-vis
    bpmap0[44, 3355:3361] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
    bpmap0[46, 311:321] |= 1  # car-20160701T00h49m36s-sci-gtoc-vis.fits
    bpmap0[46, 835:845] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
    bpmap0[46, 1156:1171] |= 1  # car-20160701T00h49m36s-sci-gtoc-vis.fits
    bpmap0[46, 1895:1905] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
    bpmap0[46, 2212:2232] |= 1  # car-20160701T00h49m36s-sci-gtoc-vis.fits
    bpmap0[47, 2127:2133] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
    bpmap0[47, 2218:2223] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
    bpmap0[47, 2260:2266] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
    bpmap0[47, 2313:2319] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
    bpmap0[47, 3111:3116] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
    bpmap0[47, 3267:3272] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
    bpmap0[47, 3316:3321] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
    bpmap0[47, 3432:3438] |= 1  # car-20170509T03h05m21s-sci-gtoc-vis
    bpmap0[47, 3480:3488] |= 1  # car-20160714T00h18m29s-sci-gtoc-vis
    bpmap0[47, 3658:3665] |= 1  # car-20170509T03h05m21s-sci-gtoc-vis
    bpmap0[49, 1008:1017] |= 1  # car-20160701T00h49m36s-sci-gtoc-vis.fits
    bpmap0[49, 2532:2544] |= 1  # car-20160701T00h49m36s-sci-gtoc-vis.fits
    bpmap0[49, 3046:3056] |= 1  # car-20160701T00h49m36s-sci-gtoc-vis.fits
    bpmap0[49, 3574:3588] |= 1  # car-20160701T00h49m36s-sci-gtoc-vis.fits

    # Constructing the bad pixel map, as defined by SERVAL in here:
    # https://github.com/mzechmeister/serval/blob/c2f47b26f1102333dfe76f93c2a686807cda02ce/src/inst_CARM_VIS.py#L95

    bpmap0[25:41, 3606:3588] |= 1  # Our rejection

    if self._internal_configs["sigma_clip_flux"] > 0:
        sigma = self._internal_configs["sigma_clip_flux"]
        for order_number in range(self.N_orders):
            cont = median_filter(self.spectra[order_number], size=500)
            inds = np.where(
                self.spectra[order_number]
                >= cont + sigma * self.uncertainties[order_number]
            )
            bpmap0[order_number, inds] |= 1
    self.spectral_mask.add_indexes_to_mask(np.where(bpmap0 != 0), QUAL_DATA)

    # remove extremely negative points!
    self.spectral_mask.add_indexes_to_mask(
        np.where(self.spectra < -3 * self.uncertainties), MISSING_DATA
    )
    self.spectral_mask.add_indexes_to_mask(
        np.where(self.uncertainties == 0), MISSING_DATA
    )

    self.assess_bad_orders()

close_arrays()

Source code in src/ASTRA/Instruments/CARMENES.py 
364
365
366
def close_arrays(self):
    super().close_arrays()
    self.is_BERV_corrected = False

load_CARMENES_extra_information(self: DataClass) -> None

CARMENES pipeline does not give RVs, we have to do an external load of the information

Parameters:

Name Type Description Default
shaq_folder str

Path to the main folder of shaq-outputs. where all the KOBE-*** targets live

 required
Source code in src/ASTRA/Instruments/CARMENES.py 
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
def load_CARMENES_extra_information(self: DataClass) -> None:
    """CARMENES pipeline does not give RVs, we have to do an external load of the information

    Parameters
    ----------
    shaq_folder : str
        Path to the main folder of shaq-outputs. where all the KOBE-*** targets live
    """

    name_to_search = self.Target.true_name

    if self.observations[0]._internal_configs["is_KOBE_data"]:
        if "KOBE-" not in name_to_search:
            name_to_search = (
                "KOBE-" + name_to_search
            )  # temporary fix for naming problem!
    else:
        logger.info(
            f"Not loading KOBE data, searching for {name_to_search} dat file with Rvs"
        )

    shaq_folder = Path(self.observations[0]._internal_configs["shaq_output_folder"])
    override_BERV = self.observations[0]._internal_configs["override_BERV"]

    if shaq_folder.name.endswith("dat"):
        logger.info("Received the previous RV file, not searching for outputs")
        shaqfile = shaq_folder
    else:
        logger.info("Searching for outputs of previous RV extraction")
        shaqfile = shaq_folder / name_to_search / f"{name_to_search}_RVs.dat"

    if shaqfile.exists():
        logger.info("Loading extra CARMENES data from {}", shaqfile)
    else:
        logger.critical(f"RV file does not exist on {shaqfile}")
        raise custom_exceptions.InvalidConfiguration("Missing RV file for data")

    number_loads = 0
    locs = []
    loaded_BJDs = [frame.get_KW_value("BJD") for frame in self.observations]
    with open(shaqfile) as file:
        for line in file:
            if "#" in line:  # header or other "BAD" files
                continue
            # TODO: implement a more thorough check in here, to mark the "bad" frames as invalid!
            ll = line.strip().split()
            if len(ll) == 0:
                logger.warning(f"shaq RV from {name_to_search} has empty line")
                continue
            bjd = round(float(ll[1]) - 2400000.0, 7)  # we have the full bjd date

            try:
                index = loaded_BJDs.index(
                    bjd
                )  # to make sure that everything is loaded in the same order
                locs.append(index)
            except ValueError:
                logger.warning("RV shaq has entry that does not exist in the S2D files")
                continue

            self.observations[index].import_KW_from_outside(
                "DRS_RV", float(ll[5]) * kilometer_second, optional=False
            )
            self.observations[index].import_KW_from_outside(
                "DRS_RV_ERR", float(ll[3]) * kilometer_second, optional=False
            )
            if override_BERV:
                self.observations[index].import_KW_from_outside(
                    "BERV", float(ll[10]) * kilometer_second, optional=False
                )
                berv_factor = 1 + float(ll[10]) / SPEED_OF_LIGHT
                self.observations[index].import_KW_from_outside(
                    "BERV_FACTOR", berv_factor, optional=False
                )

            self.observations[index].import_KW_from_outside(
                "FWHM", float(ll[11]), optional=True
            )
            self.observations[index].import_KW_from_outside(
                "BIS SPAN", float(ll[13]), optional=True
            )

            drift_val = np.nan_to_num(float(ll[7])) * meter_second
            drift_err = np.nan_to_num(float(ll[8])) * meter_second
            drift_flag = float(ll[9])
            if drift_flag > 1:
                self.observations[index].add_to_status(
                    KW_WARNING("Drift flag of KOBE is greater than 1")
                )

            self.observations[index].import_KW_from_outside(
                "drift", drift_val, optional=False
            )
            self.observations[index].import_KW_from_outside(
                "drift_ERR", drift_err, optional=False
            )

            number_loads += 1
            self.observations[index].finalized_external_data_load()

    if number_loads < len(self.observations):
        msg = f"RV shaq outputs does not have value for all S2D files of {name_to_search} ({number_loads}/{len(self.observations)})"

        logger.critical(msg)

CustomFits

Interface to MAROON-X data.

SimulatedSpirou

Bases: Frame

Interface to handle simulated data.

Source code in src/ASTRA/Instruments/CustomFits.py 
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
class SimulatedSpirou(Frame):
    """Interface to handle simulated data."""

    _default_params = Frame._default_params

    # Adding one more day than the official time, so that we ensure to include any observations
    # from the last day

    sub_instruments = {
        "SimulatedSpirou": datetime.datetime.max,
    }
    _name = "SimulatedSpirou"

    order_intervals: dict[DETECTOR_DEFINITION, slice] = {
        DETECTOR_DEFINITION.WHITE_LIGHT: list(range(48)),
    }
    KW_map = {}

    def __init__(
        self,
        file_path: Path,
        user_configs: Optional[Dict[str, Any]] = None,
        reject_subInstruments=None,
        frameID=None,
        quiet_user_params: bool = True,
    ):
        """Construct MAROON-X object.

        Parameters
        ----------
        file_path
            Path to the S2D (or S1D) file.
        user_configs
            Dictionary whose keys are the configurable options of ESPRESSO (check above)
        reject_subInstruments
            Iterable of subInstruments to fully reject
        frameID
            ID for this observation. Only used for organization purposes by :class:`~SBART.data_objects.DataClass`

        """
        self._blaze_corrected = True

        super().__init__(
            inst_name=self._name,
            array_size={"S2D": [48, 4088]},
            file_path=file_path,
            frameID=frameID,
            KW_map=self.KW_map,
            available_indicators=[],
            user_configs=user_configs,
            reject_subInstruments=reject_subInstruments,
            need_external_data_load=False,
            quiet_user_params=quiet_user_params,
        )
        coverage = [955, 2438]
        self.instrument_properties["wavelength_coverage"] = coverage
        self.instrument_properties["is_drift_corrected"] = True

        self.instrument_properties["resolution"] = 100_000

        # lat/lon from: https://geohack.toolforge.org/geohack.php?params=19_49_25_N_155_28_9_W
        lat, lon = 19.820667, -155.468056
        self.instrument_properties["EarthLocation"] = EarthLocation.from_geodetic(lat=lat, lon=lon, height=4214)

        # from https://www.mide.com/air-pressure-at-altitude-calculator
        # and convert from Pa to mbar
        self.instrument_properties["site_pressure"] = 599.4049

        self.is_BERV_corrected = True

    def get_spectral_type(self):
        return "S2D"

    def load_instrument_specific_KWs(self, header):
        # Load BERV info + previous RV
        self.observation_info["MAX_BERV"] = 30 * kilometer_second
        self.observation_info["BERV"] = header["BERV"] * meter_second

        self.observation_info["DRS_RV"] = header["DRS_RV"] * meter_second
        self.observation_info["DRS_RV_ERR"] = header["DRS_RV_ERR"] * meter_second
        self.observation_info["BJD"] = header["BJD"]

        t = Time(header["BJD"], format="jd", scale="tdb")
        iso_string = t.iso
        self.observation_info["ISO-DATE"] = "T".join(iso_string.split(" "))

        for order in range(self.N_orders):
            self.observation_info["orderwise_SNRs"].append(header[f"SNR {order}"])

    def check_header_QC(self, header: fits.Header):
        """Header QC checks for CARMENES.

        1) Drift calibration was done (CARACAL DRIFT FP REF exists)
        2) Time between observation and calibration is smaller than "max_hours_to_calibration"

        Can add the following status:

        - KW_WARNING("Drift flag of KOBE is greater than 1")
            If th drift value is greater than one. This is actually set in the DataClass.load_CARMENES_extra_information()

        Args:
            header (fits.Header): _description_
        """
        ...

    def load_S2D_data(self):
        if self.is_open:
            return
        super().load_S2D_data()

        with fits.open(self.file_path) as hdulist:
            s2d_data = hdulist["SPEC"].data * 100000  # spetra from all olders
            err_data = hdulist["SIG"].data * 100000
            wavelengths = hdulist["WAVE"].data

        self.wavelengths = wavelengths
        self.spectra = s2d_data
        self.uncertainties = err_data
        self.build_mask(bypass_QualCheck=True)
        return 1

    def load_S1D_data(self) -> Mask:
        raise NotImplementedError

    def build_mask(self, bypass_QualCheck: bool = False):
        # We evaluate the bad orders all at once
        super().build_mask(bypass_QualCheck, assess_bad_orders=False)

        # remove extremely negative points!
        self.spectral_mask.add_indexes_to_mask(np.where(self.spectra < -3 * self.uncertainties), MISSING_DATA)
        self.spectral_mask.add_indexes_to_mask(np.where(self.uncertainties == 0), MISSING_DATA)

        self.assess_bad_orders()

    def close_arrays(self):
        super().close_arrays()

sub_instruments = {'SimulatedSpirou': datetime.datetime.max} class-attribute instance-attribute

order_intervals: dict[DETECTOR_DEFINITION, slice] = {DETECTOR_DEFINITION.WHITE_LIGHT: list(range(48))} class-attribute instance-attribute

KW_map = {} class-attribute instance-attribute

is_BERV_corrected = True instance-attribute

get_spectral_type()

Source code in src/ASTRA/Instruments/CustomFits.py 
92
93
def get_spectral_type(self):
    return "S2D"

load_instrument_specific_KWs(header)

Source code in src/ASTRA/Instruments/CustomFits.py 
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
def load_instrument_specific_KWs(self, header):
    # Load BERV info + previous RV
    self.observation_info["MAX_BERV"] = 30 * kilometer_second
    self.observation_info["BERV"] = header["BERV"] * meter_second

    self.observation_info["DRS_RV"] = header["DRS_RV"] * meter_second
    self.observation_info["DRS_RV_ERR"] = header["DRS_RV_ERR"] * meter_second
    self.observation_info["BJD"] = header["BJD"]

    t = Time(header["BJD"], format="jd", scale="tdb")
    iso_string = t.iso
    self.observation_info["ISO-DATE"] = "T".join(iso_string.split(" "))

    for order in range(self.N_orders):
        self.observation_info["orderwise_SNRs"].append(header[f"SNR {order}"])

check_header_QC(header: fits.Header)

Header QC checks for CARMENES.

1) Drift calibration was done (CARACAL DRIFT FP REF exists) 2) Time between observation and calibration is smaller than "max_hours_to_calibration"

Can add the following status:

  • KW_WARNING("Drift flag of KOBE is greater than 1") If th drift value is greater than one. This is actually set in the DataClass.load_CARMENES_extra_information()

Args: header (fits.Header): description

Source code in src/ASTRA/Instruments/CustomFits.py 
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
def check_header_QC(self, header: fits.Header):
    """Header QC checks for CARMENES.

    1) Drift calibration was done (CARACAL DRIFT FP REF exists)
    2) Time between observation and calibration is smaller than "max_hours_to_calibration"

    Can add the following status:

    - KW_WARNING("Drift flag of KOBE is greater than 1")
        If th drift value is greater than one. This is actually set in the DataClass.load_CARMENES_extra_information()

    Args:
        header (fits.Header): _description_
    """
    ...

load_S2D_data()

Source code in src/ASTRA/Instruments/CustomFits.py 
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
def load_S2D_data(self):
    if self.is_open:
        return
    super().load_S2D_data()

    with fits.open(self.file_path) as hdulist:
        s2d_data = hdulist["SPEC"].data * 100000  # spetra from all olders
        err_data = hdulist["SIG"].data * 100000
        wavelengths = hdulist["WAVE"].data

    self.wavelengths = wavelengths
    self.spectra = s2d_data
    self.uncertainties = err_data
    self.build_mask(bypass_QualCheck=True)
    return 1

load_S1D_data() -> Mask

Source code in src/ASTRA/Instruments/CustomFits.py 
143
144
def load_S1D_data(self) -> Mask:
    raise NotImplementedError

build_mask(bypass_QualCheck: bool = False)

Source code in src/ASTRA/Instruments/CustomFits.py 
146
147
148
149
150
151
152
153
154
def build_mask(self, bypass_QualCheck: bool = False):
    # We evaluate the bad orders all at once
    super().build_mask(bypass_QualCheck, assess_bad_orders=False)

    # remove extremely negative points!
    self.spectral_mask.add_indexes_to_mask(np.where(self.spectra < -3 * self.uncertainties), MISSING_DATA)
    self.spectral_mask.add_indexes_to_mask(np.where(self.uncertainties == 0), MISSING_DATA)

    self.assess_bad_orders()

close_arrays()

Source code in src/ASTRA/Instruments/CustomFits.py 
156
157
def close_arrays(self):
    super().close_arrays()

ESO_PIPELINE

ESO_PIPELINE

Bases: Frame

Interface to handle data products (S2D and S1D) of the ESO pipeline.

User parameters:

================================ ================ ================ ================ ================ Parameter name Mandatory Default Value Valid Values Comment ================================ ================ ================ ================ ================ apply_FluxCorr False False boolean Apply polynomial flux correction to flux to preserve red-blue balance Telluric_Corrected False False boolean Use S2D data in the format given by R. Allart's correction tool ================================ ================ ================ ================ ================

.. note:: Also check the User parameters of the parent classes for further customization options of SBART

Subclassed by:
Source code in src/ASTRA/Instruments/ESO_PIPELINE.py 
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
class ESO_PIPELINE(Frame):
    """Interface to handle data products (S2D and S1D) of the ESO pipeline.

    **User parameters:**

    ================================ ================ ================ ================ ================
    Parameter name                      Mandatory      Default Value    Valid Values    Comment
    ================================ ================ ================ ================ ================
    apply_FluxCorr                    False              False         boolean          Apply polynomial flux correction to flux to preserve red-blue balance
    Telluric_Corrected                False              False         boolean          Use S2D data in the format given by R. Allart's correction tool
    ================================ ================ ================ ================ ================

    .. note::
        Also check the **User parameters** of the parent classes for further customization options of SBART

    """

    _default_params = Frame._default_params + DefaultValues(
        Telluric_Corrected=UserParam(
            False,
            constraint=BooleanValue,
            description="The Frame was already corrected from telluric features",
        ),
        UseMolecfit=UserParam(False, constraint=BooleanValue),
        use_old_pipeline=UserParam(
            default_value=False,
            constraint=BooleanValue,
            description="Use data from the old pipeline. Only available to selected instruments",
        ),
        SCIRED_CHECK_IS_FATAL=UserParam(
            default_value=True,
            constraint=BooleanValue,
            description="Automatically reject frames with QC SCIRED CHECK = 0 ",
        ),
    )

    _default_params.update(
        "apply_FluxCorr",
        UserParam(
            False,
            constraint=BooleanValue,
            description="Apply the blue-red flux correction due to the wavelength dependence of the atmospheric extinction. Only available on data from ESO pipeline (ESPRESSO)",
        ),
    )

    _default_params.update(
        "apply_FluxBalance_Norm",
        UserParam(False, constraint=BooleanValue),
    )

    def __init__(
        self,
        inst_name: str,
        array_size: dict,
        file_path,
        KW_identifier: str,
        user_configs: Optional[Dict[str, Any]] = None,
        reject_subInstruments: Optional[Iterable[str]] = None,
        frameID: Optional[int] = None,
        quiet_user_params: bool = True,
        override_KW_map=None,
        override_indicators=None,
        need_external_data_load=False,
        is_solar: bool = False,
    ):
        """Prepare data load.

        Parameters
        ----------
        file_path
            Path to the S2D (or S1D) file.
        user_configs
            Dictionary whose keys are the configurable options of ESPRESSO (check above)
        reject_subInstruments
            Iterable of subInstruments to fully reject
        frameID
            ID for this observation. Only used for organization purposes by :class:`~SBART.data_objects.DataClass`

        """
        # Wavelength coverage

        self.UT_number = None
        self.KW_identifier = KW_identifier

        KW_map = {
            "OBJECT": "OBJECT",
            "BJD": f"HIERARCH {KW_identifier} QC BJD",
            "MJD": "MJD-OBS",
            "ISO-DATE": "DATE-OBS",
            "DRS-VERSION": "HIERARCH ESO PRO REC1 PIPE ID",
            "MD5-CHECK": "DATAMD5",
            "RA": "RA",
            "DEC": "DEC",
            "SPEC_TYPE": f"HIERARCH {KW_identifier} QC CCF MASK",
            "EXPTIME": "EXPTIME",
            "DRS_CCF_MASK": f"HIERARCH {KW_identifier} QC CCF MASK",
            "DRS_FLUX_CORRECTION_TEMPLATE": "HIERARCH ESO PRO REC1 CAL8 NAME",
            "INS NAME": "INSTRUME",
            "INS MODE": f"HIERARCH {KW_identifier} INS MODE",
            "PROG ID": f"HIERARCH {KW_identifier} OBS PROG ID",
            "OBS NAME": f"HIERARCH {KW_identifier} OBS NAME",
        }
        if override_KW_map is not None:
            for key, value in override_KW_map.items():
                if value is None:
                    # Allow to bypass some values
                    del KW_map[key]
                else:
                    KW_map[key] = value

        available_indicators = ("CONTRAST", "FWHM", "BIS SPAN")

        super().__init__(
            inst_name=inst_name,
            array_size=array_size,
            file_path=file_path,
            frameID=frameID,
            KW_map=KW_map,
            available_indicators=(
                available_indicators
                if override_indicators is None
                else override_indicators
            ),
            user_configs=user_configs,
            reject_subInstruments=reject_subInstruments,
            quiet_user_params=quiet_user_params,
            need_external_data_load=need_external_data_load,
            is_solar=is_solar,
        )

        self.instrument_properties["is_drift_corrected"] = True

        self.is_BERV_corrected = True

        # CHeck for BLAZE correction
        self.is_blaze_corrected = True
        if "BLAZE" in self.file_path.stem:
            # The S2D_BLAZE_A files do not have the blaze correction!
            self.is_blaze_corrected = False

        if not self._internal_configs["use_old_pipeline"]:
            logger.debug("Checking for validity of file reduced with new pipeline")
            extra_names = ["BLAZE", "SKYSUB", "BLAZE_TELL_CORR", "TELL_CORR"]
            valid = []
            for fib in ["A", "B"]:
                valid.extend([f"S2D_{mod}_{fib}" for mod in extra_names])
                valid.extend([f"S2D_{fib}", f"S1D_{fib}"])
            for key in valid:
                if key in self.file_path.stem:
                    break
            else:
                raise custom_exceptions.InvalidConfiguration(
                    f"{self.name} can't recognize {self.file_path}"
                )

    def load_instrument_specific_KWs(self, header) -> None:
        if self._internal_configs["use_old_pipeline"]:
            self._load_old_DRS_KWs(header)
        else:
            self._load_ESO_DRS_KWs(header)
        self.load_telemetry_info(header)

        drs_version = header.get("HIERARCH ESO PRO REC1 PIPE ID", "espdr/3.0.0")
        version = drs_version.split("/")[-1]

        # TODO: Check the number of the DRS version that has this fixed
        version_sum = lambda a: sum(
            a * 10**b for a, b in zip(map(int, a.split(".")[::-1]), range(3))
        )

        if version_sum(version) < version_sum("3.2.1"):
            # For older versions we always need to use the
            # approximated BERV correction
            self.use_approximated_BERV_correction = True
        else:
            # Newer versions don't need the approximated BERV correction
            self.use_approximated_BERV_correction = False

            berv_factor = self.observation_info["BERV_FACTOR"]
            if berv_factor is not None:
                logger.warning("Recomputing the BERV from BERV_factor keyword")
                new_berv = (berv_factor - 1) * SPEED_OF_LIGHT * kilometer_second
                diff_berv = (
                    (new_berv - self.observation_info["BERV"]).to(meter_second).value
                )
                if abs(diff_berv) > 10:
                    logger.warning(
                        f"Difference between BERV and BERV_factor is of {diff_berv} [m/s]"
                    )

                self.observation_info["BERV"] = new_berv

    def load_S2D_data(self):
        if self.is_open:
            logger.debug("{} has already been opened", self.__str__())
            return None
        super().load_S2D_data()

        if self._internal_configs["use_old_pipeline"]:
            # The return is to ensure that we don't do anything after this point!
            return self.load_old_DRS_S2D()
        self.load_ESO_DRS_S2D_data()

    def load_S1D_data(self):
        super().load_S1D_data()

        if self._internal_configs["use_old_pipeline"]:
            # The return is to ensure that we don't do anything after this point!
            return self.load_old_DRS_S1D()
        self.load_ESO_DRS_S1D_data()

    def get_air_wavelength(self) -> np.ndarray:
        if self.is_S1D:
            with fits.open(self.file_path) as hdu:
                return hdu[1].data["wavelength_air"].reshape((1, self.array_size[1]))
        elif self.is_S2D:
            with fits.open(self.file_path) as hdu:
                return hdu["WAVEDATA_AIR_BARY"].data

    def check_header_QC(self, header: fits.header.Header):
        super().check_header_QC(header)
        if self._internal_configs["use_old_pipeline"]:
            self.check_header_QC_old_DRS(header)
        else:
            self.check_header_QC_ESO_DRS(header)

    ## Add the bindings for a possible old DRS version

    def _load_old_DRS_KWs(self, header):
        raise NotImplementedError(f"There is no old pipeline for {self.name}")

    def check_header_QC_old_DRS(self, header):
        raise NotImplementedError(f"There is no old pipeline for {self.name}")

    def load_old_DRS_S2D(self):
        raise NotImplementedError(f"There is no old pipeline for {self.name}")

    def load_old_DRS_S1D(self):
        raise NotImplementedError(f"There is no old pipeline for {self.name}")

    # Load common information to the ESO DRS version across different spectrographs
    def load_telemetry_info(self, header):
        """Load (at least) the following keywords.

        - relative humidity
        - ambient temperature, in Celsius
        - airmass
        - Detector

        Parameters
        ----------
        header


        """
        ambi_KWs = {
            "relative_humidity": "HUMIDITY",
            "ambient_temperature": "TEMP10M",
        }

        for name, endKW in ambi_KWs.items():
            self.observation_info[name] = header[
                f"HIERARCH {self.KW_identifier} METEO {endKW}"
            ]
            if "temperature" in name:  # store temperature in KELVIN for TELFIT
                self.observation_info[name] = convert_temperature(
                    self.observation_info[name],
                    old_scale="Celsius",
                    new_scale="Kelvin",
                )

        if self.observation_info["relative_humidity"] == 255:
            logger.warning(
                f"{self.name} has an invalid value in the humidity sensor..."
            )
            self.observation_info["relative_humidity"] = np.nan

        self.observation_info["airmass"] = header["AIRMASS"]

    def _load_ESO_DRS_KWs(self, header):
        if self._internal_configs["use_old_pipeline"]:
            raise custom_exceptions.InvalidConfiguration(
                "Can't load data from new pipeline with the config for the old one",
            )

        # Load BERV info + previous RV
        self.observation_info["MAX_BERV"] = (
            header[f"HIERARCH {self.KW_identifier} QC BERVMAX"] * kilometer_second
        )
        self.observation_info["BERV"] = (
            header[f"HIERARCH {self.KW_identifier} QC BERV"] * kilometer_second
        )

        berv_factor = header.get(f"HIERARCH {self.KW_identifier} QC BERV FACTOR", None)
        self.observation_info["BERV_FACTOR"] = berv_factor

        self.observation_info["DRS_RV"] = (
            header[f"HIERARCH {self.KW_identifier} QC CCF RV"] * kilometer_second
        )
        self.observation_info["DRS_RV_ERR"] = (
            header[f"HIERARCH {self.KW_identifier} QC CCF RV ERROR"] * kilometer_second
        )

        for order in range(self.instrument_properties["array_sizes"]["S2D"][0]):
            self.observation_info["orderwise_SNRs"].append(
                header[f"HIERARCH {self.KW_identifier} QC ORDER{order + 1} SNR"],
            )

        # Chosen activity indicators
        for key in self.available_indicators:
            full_key = f"HIERARCH {self.KW_identifier} QC CCF {key}"
            self.observation_info[key] = header[full_key]
            self.observation_info[key + "_ERR"] = header[full_key + " ERROR"]

        obs_date = header["DATE-OBS"]

        obs_date = "-".join(obs_date.split("T")).split(":")[0]
        obs_date = datetime.datetime.strptime(obs_date, r"%Y-%m-%d-%H")

        if obs_date.hour < 16:  # noqa: PLR2004
            # If before 4pm, the observation was from previous day
            obs_date = obs_date - datetime.timedelta(days=1)
        self.observation_info["DATE_NIGHT"] = datetime.datetime.strftime(
            obs_date, r"%Y-%m-%d"
        )

    def load_ESO_DRS_S2D_data(self, overload_SCIDATA_key=None):
        if self._internal_configs["use_old_pipeline"]:
            raise custom_exceptions.InvalidConfiguration(
                "Can't load data from new pipeline with the config for the old one",
            )

        with fits.open(self.file_path) as hdulist:
            self.wavelengths = hdulist["WAVEDATA_VAC_BARY"].data
            self.qual_data = hdulist["QUALDATA"].data

            SCIDATA_KEY = (
                "SCIDATA" if overload_SCIDATA_key is None else overload_SCIDATA_key
            )
            ERRDATA_KEY = "ERRDATA"

            if self._internal_configs["Telluric_Corrected"]:
                logger.info(
                    "Loading S2D file from a non-DRS source (telluric corrected file)"
                )
                SCIDATA_KEY += "_CORR"
                ERRDATA_KEY += "_CORR"

            # Fixing dtype to avoid problems with the cython interface
            self.spectra = hdulist[SCIDATA_KEY].data.astype(np.float64)
            self.uncertainties = hdulist[ERRDATA_KEY].data.astype(np.float64)
            if self._internal_configs["apply_FluxCorr"]:
                logger.debug("Starting chromatic flux correction")
                keyword = f"HIERARCH {self.KW_identifier} QC ORDER%d FLUX CORR"
                flux_corr = np.array(
                    [
                        hdulist[0].header[keyword % o]
                        for o in range(1, self.N_orders + 1)
                    ]
                )
                fit_nb = (flux_corr != 1.0).sum()

                ignore = self.N_orders - fit_nb
                logger.debug(f"Chromatic correction ignoring {ignore} orders")

                # ? see espdr_scince:espdr_correct_flux
                poly_deg = round(8 * fit_nb / self.N_orders)
                logger.debug(f"Fitting polynomial (n={poly_deg})")
                llc = hdulist[5].data[:, self.array_size[1] // 2]
                coeff = np.polyfit(llc[ignore:], flux_corr[ignore:], poly_deg - 1)
                # corr_model = np.zeros_like(hdu[5].data, dtype=np.float32)
                corr_model = np.polyval(coeff, hdulist[5].data)

                corr_model[flux_corr == 1] = (
                    1  # orders where the CORR FACTOR are 1 do not have correction!
                )
                self.spectra = (
                    self.spectra / corr_model
                )  # correct from chromatic variations
                self.flux_atmos_balance_corrected = True
                # TODO: understand if we want to include the factor in uncertainties or not!
                # self.uncertainties = self.uncertainties / corr_model # maintain the SNR in the corrected spectrum

            else:
                # Disabled the flux correction as we are artifically increasing the SNR of the spectra...
                # Shouldn't we also increase the flux uncertainty?? The DRS does not do it....
                logger.warning("Not applying correction to blue-red flux balance!")
                # / corr_model

            if self._internal_configs["apply_FluxBalance_Norm"]:
                logger.info(
                    "Normalizing the flux balance distribution due to dispersion"
                )
                # The physical sizes of the pixels (on the CCD) are the same
                # The flux that reaches eeach pixel is different, due to dispersion
                # The spectra will have a trend, even after removing the instrumental effect
                # This normalizes the spectra by dividing by the flux distribution

                balance_corr_model = hdulist["DLLDATA_VAC_BARY"].data
                self.spectra = self.spectra / balance_corr_model
                # Ensure that we keep the same SNR after the normalization!
                self.uncertainties = self.uncertainties / balance_corr_model
                self.flux_dispersion_balance_corrected = True

        self.build_mask(bypass_QualCheck=False)
        return 1

    def load_ESO_DRS_S1D_data(self):
        if self._internal_configs["use_old_pipeline"]:
            raise custom_exceptions.InvalidConfiguration(
                "Can't load data from new pipeline with the config for the old one",
            )

        with fits.open(self.file_path) as hdulist:
            full_data = hdulist[1].data

        wave_kw = "wavelength"
        if self._internal_configs["use_air_wavelengths"]:
            wave_kw = "wavelength_air"
            logger.warning("SBART using air wavelengths!")

        self.wavelengths = full_data[wave_kw].reshape((1, self.array_size[1]))
        self.spectra = (
            full_data["flux"].reshape((1, self.array_size[1])).astype(np.float64)
        )
        self.uncertainties = (
            full_data["error"].reshape((1, self.array_size[1])).astype(np.float64)
        )
        self.qual_data = full_data["quality"].reshape((1, self.array_size[1]))
        self.build_mask(bypass_QualCheck=False)

    def _compute_BLAZE(self):
        """Assume that S2D and S2D_BLAZE live within the same folder (should be true for most cases)

        Returns
        -------

        """
        if "BLAZE" in self.file_path:
            blaze_corr_file = self.file_path
            s2d_file = self.file_path.replace("S2D_BLAZE_A", "S2D_A")
        else:
            s2d_file = self.file_path
            blaze_corr_file = self.file_path.replace("S2D_A", "S2D_BLAZE_A")

        with fits.open(s2d_file) as hdu:
            S2D = hdu[1].data
        with fits.open(blaze_corr_file) as hdu:
            blaze_corr = hdu[1].data

        # this will raise errors on the edges of the orders (zero-division)
        self._blaze_function = blaze_corr / S2D

    def check_header_QC_ESO_DRS(self, header):
        fatal_QC_flags = {}

        nonfatal_QC_flags = {
            f"HIERARCH {self.KW_identifier}" + " QC SCIRED FLUX CORR CHECK": 0,
        }
        if not self.is_skysub:
            extra_checks = {
                f"HIERARCH {self.KW_identifier}" + " QC SCIRED DRIFT CHECK": 0,
                f"HIERARCH {self.KW_identifier}"
                + " QC SCIRED DRIFT FLUX_RATIO CHECK": 0,
                f"HIERARCH {self.KW_identifier}" + " QC SCIRED DRIFT CHI2 CHECK": 0,
            }
            nonfatal_QC_flags = {**nonfatal_QC_flags, **extra_checks}

        if self._internal_configs["SCIRED_CHECK_IS_FATAL"]:
            fatal_QC_flags[f"HIERARCH {self.KW_identifier} QC SCIRED CHECK"] = 0
        else:
            nonfatal_QC_flags[f"HIERARCH {self.KW_identifier} QC SCIRED CHECK"] = 0

        for flag, bad_value in fatal_QC_flags.items():
            if header[flag] == bad_value:
                msg = f"\tQC flag {flag} has taken the bad value of {bad_value}"
                logger.critical(msg)
                self.add_to_status(FATAL_KW(msg.replace("\t", "")))

        for flag, bad_value in nonfatal_QC_flags.items():
            try:
                if header[flag] == bad_value:
                    msg = f"QC flag {flag} meets the bad value"
                    logger.warning(msg)
                    self._status.store_warning(KW_WARNING(msg))
            except KeyError:
                msg = f"QC flag {flag} does not exist"
                logger.warning(msg)
                self._status.store_warning(KW_WARNING(msg))

        if self._status.number_warnings > 0:
            logger.warning(
                "Found {} warning flags in the header KWs", self._status.number_warnings
            )

        espdr_to_num = lambda x: int("".join(x.split(".")))
        espdrversion = header["ESO PRO REC1 PIPE ID"].split("/")[-1]

        if self._internal_configs["USE_APPROX_BERV_CORRECTION"]:
            if espdr_to_num(espdrversion) >= espdr_to_num("3.2.0"):
                logger.critical(
                    f"Using approximated BERV correction in espdr/{espdrversion}"
                )
        else:
            if espdr_to_num(espdrversion) < espdr_to_num("3.2.0"):
                logger.critical(
                    f"Not using approximated BERV correction in espdr/{espdrversion}"
                )

    @property
    def bare_fname(self) -> str:
        return self.fname.split("_S")[0]

    def trigger_data_storage(self, *args, **kwargs):
        super().trigger_data_storage(*args, **kwargs)

UT_number = None instance-attribute

KW_identifier = KW_identifier instance-attribute

is_BERV_corrected = True instance-attribute

is_blaze_corrected = True instance-attribute

bare_fname: str property

load_instrument_specific_KWs(header) -> None

Source code in src/ASTRA/Instruments/ESO_PIPELINE.py 
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
def load_instrument_specific_KWs(self, header) -> None:
    if self._internal_configs["use_old_pipeline"]:
        self._load_old_DRS_KWs(header)
    else:
        self._load_ESO_DRS_KWs(header)
    self.load_telemetry_info(header)

    drs_version = header.get("HIERARCH ESO PRO REC1 PIPE ID", "espdr/3.0.0")
    version = drs_version.split("/")[-1]

    # TODO: Check the number of the DRS version that has this fixed
    version_sum = lambda a: sum(
        a * 10**b for a, b in zip(map(int, a.split(".")[::-1]), range(3))
    )

    if version_sum(version) < version_sum("3.2.1"):
        # For older versions we always need to use the
        # approximated BERV correction
        self.use_approximated_BERV_correction = True
    else:
        # Newer versions don't need the approximated BERV correction
        self.use_approximated_BERV_correction = False

        berv_factor = self.observation_info["BERV_FACTOR"]
        if berv_factor is not None:
            logger.warning("Recomputing the BERV from BERV_factor keyword")
            new_berv = (berv_factor - 1) * SPEED_OF_LIGHT * kilometer_second
            diff_berv = (
                (new_berv - self.observation_info["BERV"]).to(meter_second).value
            )
            if abs(diff_berv) > 10:
                logger.warning(
                    f"Difference between BERV and BERV_factor is of {diff_berv} [m/s]"
                )

            self.observation_info["BERV"] = new_berv

load_S2D_data()

Source code in src/ASTRA/Instruments/ESO_PIPELINE.py 
210
211
212
213
214
215
216
217
218
219
def load_S2D_data(self):
    if self.is_open:
        logger.debug("{} has already been opened", self.__str__())
        return None
    super().load_S2D_data()

    if self._internal_configs["use_old_pipeline"]:
        # The return is to ensure that we don't do anything after this point!
        return self.load_old_DRS_S2D()
    self.load_ESO_DRS_S2D_data()

load_S1D_data()

Source code in src/ASTRA/Instruments/ESO_PIPELINE.py 
221
222
223
224
225
226
227
def load_S1D_data(self):
    super().load_S1D_data()

    if self._internal_configs["use_old_pipeline"]:
        # The return is to ensure that we don't do anything after this point!
        return self.load_old_DRS_S1D()
    self.load_ESO_DRS_S1D_data()

get_air_wavelength() -> np.ndarray

Source code in src/ASTRA/Instruments/ESO_PIPELINE.py 
229
230
231
232
233
234
235
def get_air_wavelength(self) -> np.ndarray:
    if self.is_S1D:
        with fits.open(self.file_path) as hdu:
            return hdu[1].data["wavelength_air"].reshape((1, self.array_size[1]))
    elif self.is_S2D:
        with fits.open(self.file_path) as hdu:
            return hdu["WAVEDATA_AIR_BARY"].data

check_header_QC(header: fits.header.Header)

Source code in src/ASTRA/Instruments/ESO_PIPELINE.py 
237
238
239
240
241
242
def check_header_QC(self, header: fits.header.Header):
    super().check_header_QC(header)
    if self._internal_configs["use_old_pipeline"]:
        self.check_header_QC_old_DRS(header)
    else:
        self.check_header_QC_ESO_DRS(header)

check_header_QC_old_DRS(header)

Source code in src/ASTRA/Instruments/ESO_PIPELINE.py 
249
250
def check_header_QC_old_DRS(self, header):
    raise NotImplementedError(f"There is no old pipeline for {self.name}")

load_old_DRS_S2D()

Source code in src/ASTRA/Instruments/ESO_PIPELINE.py 
252
253
def load_old_DRS_S2D(self):
    raise NotImplementedError(f"There is no old pipeline for {self.name}")

load_old_DRS_S1D()

Source code in src/ASTRA/Instruments/ESO_PIPELINE.py 
255
256
def load_old_DRS_S1D(self):
    raise NotImplementedError(f"There is no old pipeline for {self.name}")

load_telemetry_info(header)

Load (at least) the following keywords.

  • relative humidity
  • ambient temperature, in Celsius
  • airmass
  • Detector

Parameters:

Name Type Description Default
header
required
Source code in src/ASTRA/Instruments/ESO_PIPELINE.py 
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
def load_telemetry_info(self, header):
    """Load (at least) the following keywords.

    - relative humidity
    - ambient temperature, in Celsius
    - airmass
    - Detector

    Parameters
    ----------
    header


    """
    ambi_KWs = {
        "relative_humidity": "HUMIDITY",
        "ambient_temperature": "TEMP10M",
    }

    for name, endKW in ambi_KWs.items():
        self.observation_info[name] = header[
            f"HIERARCH {self.KW_identifier} METEO {endKW}"
        ]
        if "temperature" in name:  # store temperature in KELVIN for TELFIT
            self.observation_info[name] = convert_temperature(
                self.observation_info[name],
                old_scale="Celsius",
                new_scale="Kelvin",
            )

    if self.observation_info["relative_humidity"] == 255:
        logger.warning(
            f"{self.name} has an invalid value in the humidity sensor..."
        )
        self.observation_info["relative_humidity"] = np.nan

    self.observation_info["airmass"] = header["AIRMASS"]

load_ESO_DRS_S2D_data(overload_SCIDATA_key=None)

Source code in src/ASTRA/Instruments/ESO_PIPELINE.py 
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
def load_ESO_DRS_S2D_data(self, overload_SCIDATA_key=None):
    if self._internal_configs["use_old_pipeline"]:
        raise custom_exceptions.InvalidConfiguration(
            "Can't load data from new pipeline with the config for the old one",
        )

    with fits.open(self.file_path) as hdulist:
        self.wavelengths = hdulist["WAVEDATA_VAC_BARY"].data
        self.qual_data = hdulist["QUALDATA"].data

        SCIDATA_KEY = (
            "SCIDATA" if overload_SCIDATA_key is None else overload_SCIDATA_key
        )
        ERRDATA_KEY = "ERRDATA"

        if self._internal_configs["Telluric_Corrected"]:
            logger.info(
                "Loading S2D file from a non-DRS source (telluric corrected file)"
            )
            SCIDATA_KEY += "_CORR"
            ERRDATA_KEY += "_CORR"

        # Fixing dtype to avoid problems with the cython interface
        self.spectra = hdulist[SCIDATA_KEY].data.astype(np.float64)
        self.uncertainties = hdulist[ERRDATA_KEY].data.astype(np.float64)
        if self._internal_configs["apply_FluxCorr"]:
            logger.debug("Starting chromatic flux correction")
            keyword = f"HIERARCH {self.KW_identifier} QC ORDER%d FLUX CORR"
            flux_corr = np.array(
                [
                    hdulist[0].header[keyword % o]
                    for o in range(1, self.N_orders + 1)
                ]
            )
            fit_nb = (flux_corr != 1.0).sum()

            ignore = self.N_orders - fit_nb
            logger.debug(f"Chromatic correction ignoring {ignore} orders")

            # ? see espdr_scince:espdr_correct_flux
            poly_deg = round(8 * fit_nb / self.N_orders)
            logger.debug(f"Fitting polynomial (n={poly_deg})")
            llc = hdulist[5].data[:, self.array_size[1] // 2]
            coeff = np.polyfit(llc[ignore:], flux_corr[ignore:], poly_deg - 1)
            # corr_model = np.zeros_like(hdu[5].data, dtype=np.float32)
            corr_model = np.polyval(coeff, hdulist[5].data)

            corr_model[flux_corr == 1] = (
                1  # orders where the CORR FACTOR are 1 do not have correction!
            )
            self.spectra = (
                self.spectra / corr_model
            )  # correct from chromatic variations
            self.flux_atmos_balance_corrected = True
            # TODO: understand if we want to include the factor in uncertainties or not!
            # self.uncertainties = self.uncertainties / corr_model # maintain the SNR in the corrected spectrum

        else:
            # Disabled the flux correction as we are artifically increasing the SNR of the spectra...
            # Shouldn't we also increase the flux uncertainty?? The DRS does not do it....
            logger.warning("Not applying correction to blue-red flux balance!")
            # / corr_model

        if self._internal_configs["apply_FluxBalance_Norm"]:
            logger.info(
                "Normalizing the flux balance distribution due to dispersion"
            )
            # The physical sizes of the pixels (on the CCD) are the same
            # The flux that reaches eeach pixel is different, due to dispersion
            # The spectra will have a trend, even after removing the instrumental effect
            # This normalizes the spectra by dividing by the flux distribution

            balance_corr_model = hdulist["DLLDATA_VAC_BARY"].data
            self.spectra = self.spectra / balance_corr_model
            # Ensure that we keep the same SNR after the normalization!
            self.uncertainties = self.uncertainties / balance_corr_model
            self.flux_dispersion_balance_corrected = True

    self.build_mask(bypass_QualCheck=False)
    return 1

load_ESO_DRS_S1D_data()

Source code in src/ASTRA/Instruments/ESO_PIPELINE.py 
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
def load_ESO_DRS_S1D_data(self):
    if self._internal_configs["use_old_pipeline"]:
        raise custom_exceptions.InvalidConfiguration(
            "Can't load data from new pipeline with the config for the old one",
        )

    with fits.open(self.file_path) as hdulist:
        full_data = hdulist[1].data

    wave_kw = "wavelength"
    if self._internal_configs["use_air_wavelengths"]:
        wave_kw = "wavelength_air"
        logger.warning("SBART using air wavelengths!")

    self.wavelengths = full_data[wave_kw].reshape((1, self.array_size[1]))
    self.spectra = (
        full_data["flux"].reshape((1, self.array_size[1])).astype(np.float64)
    )
    self.uncertainties = (
        full_data["error"].reshape((1, self.array_size[1])).astype(np.float64)
    )
    self.qual_data = full_data["quality"].reshape((1, self.array_size[1]))
    self.build_mask(bypass_QualCheck=False)

check_header_QC_ESO_DRS(header)

Source code in src/ASTRA/Instruments/ESO_PIPELINE.py 
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
def check_header_QC_ESO_DRS(self, header):
    fatal_QC_flags = {}

    nonfatal_QC_flags = {
        f"HIERARCH {self.KW_identifier}" + " QC SCIRED FLUX CORR CHECK": 0,
    }
    if not self.is_skysub:
        extra_checks = {
            f"HIERARCH {self.KW_identifier}" + " QC SCIRED DRIFT CHECK": 0,
            f"HIERARCH {self.KW_identifier}"
            + " QC SCIRED DRIFT FLUX_RATIO CHECK": 0,
            f"HIERARCH {self.KW_identifier}" + " QC SCIRED DRIFT CHI2 CHECK": 0,
        }
        nonfatal_QC_flags = {**nonfatal_QC_flags, **extra_checks}

    if self._internal_configs["SCIRED_CHECK_IS_FATAL"]:
        fatal_QC_flags[f"HIERARCH {self.KW_identifier} QC SCIRED CHECK"] = 0
    else:
        nonfatal_QC_flags[f"HIERARCH {self.KW_identifier} QC SCIRED CHECK"] = 0

    for flag, bad_value in fatal_QC_flags.items():
        if header[flag] == bad_value:
            msg = f"\tQC flag {flag} has taken the bad value of {bad_value}"
            logger.critical(msg)
            self.add_to_status(FATAL_KW(msg.replace("\t", "")))

    for flag, bad_value in nonfatal_QC_flags.items():
        try:
            if header[flag] == bad_value:
                msg = f"QC flag {flag} meets the bad value"
                logger.warning(msg)
                self._status.store_warning(KW_WARNING(msg))
        except KeyError:
            msg = f"QC flag {flag} does not exist"
            logger.warning(msg)
            self._status.store_warning(KW_WARNING(msg))

    if self._status.number_warnings > 0:
        logger.warning(
            "Found {} warning flags in the header KWs", self._status.number_warnings
        )

    espdr_to_num = lambda x: int("".join(x.split(".")))
    espdrversion = header["ESO PRO REC1 PIPE ID"].split("/")[-1]

    if self._internal_configs["USE_APPROX_BERV_CORRECTION"]:
        if espdr_to_num(espdrversion) >= espdr_to_num("3.2.0"):
            logger.critical(
                f"Using approximated BERV correction in espdr/{espdrversion}"
            )
    else:
        if espdr_to_num(espdrversion) < espdr_to_num("3.2.0"):
            logger.critical(
                f"Not using approximated BERV correction in espdr/{espdrversion}"
            )

trigger_data_storage(*args, **kwargs)

Source code in src/ASTRA/Instruments/ESO_PIPELINE.py 
531
532
def trigger_data_storage(self, *args, **kwargs):
    super().trigger_data_storage(*args, **kwargs)

ESPRESSO

ESPRESSO

Bases: ESO_PIPELINE

Interface to handle ESPRESSO observations (S2D and S1D).

With ESPRESSO data we are considering 3 sub-Instruments:

  • ESPRESSO18 - Before 2019-06-27
  • ESPRESSO19 - Before 2020-12-18
  • ESPRESSO21 - Until the ends of time (hopefully)

User parameters:

================================ ================ ================ ================ ================ Parameter name Mandatory Default Value Valid Values Comment ================================ ================ ================ ================ ================ ================================ ================ ================ ================ ================

.. note:: Also check the User parameters of the parent classes for further customization options of SBART

Subclassed by:
Source code in src/ASTRA/Instruments/ESPRESSO.py 
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
class ESPRESSO(ESO_PIPELINE):
    """Interface to handle ESPRESSO observations (S2D and S1D).

    With ESPRESSO data we are considering 3 sub-Instruments:

    * ESPRESSO18 - Before  2019-06-27
    * ESPRESSO19 - Before  2020-12-18
    * ESPRESSO21 - Until the ends of time (hopefully)


    **User parameters:**

    ================================ ================ ================ ================ ================
    Parameter name                      Mandatory      Default Value    Valid Values    Comment
    ================================ ================ ================ ================ ================
    ================================ ================ ================ ================ ================

    .. note::
        Also check the **User parameters** of the parent classes for further customization options of SBART

    """

    _default_params = ESO_PIPELINE._default_params

    sub_instruments = {
        "ESPRESSO18": datetime.datetime.strptime("2019-06-27", r"%Y-%m-%d"),
        "ESPRESSO19": datetime.datetime.max,
        # "ESPRESSO19": datetime.datetime.strptime("2020-12-18", r"%Y-%m-%d"),
    }
    _name = "ESPRESSO"

    order_intervals: dict[DETECTOR_DEFINITION, slice] = {
        DETECTOR_DEFINITION.WHITE_LIGHT: list(range(170)),
        DETECTOR_DEFINITION.RED_DET: list(range(90, 170)),
        DETECTOR_DEFINITION.BLUE_DET: list(range(0, 90)),
    }

    def __init__(
        self,
        file_path,
        user_configs: Optional[Dict[str, Any]] = None,
        reject_subInstruments: Optional[Iterable[str]] = None,
        frameID: Optional[int] = None,
        quiet_user_params: bool = True,
        is_solar=False,
    ):
        """ESPRESSO interface.

        Parameters
        ----------
        file_path
            Path to the S2D (or S1D) file.
        user_configs
            Dictionary whose keys are the configurable options of ESPRESSO (check above)
        reject_subInstruments
            Iterable of subInstruments to fully reject
        frameID
            ID for this observation. Only used for organization purposes by :class:`~SBART.data_objects.DataClass`

        """
        # Wavelength coverage

        coverage = (350, 900)
        super().__init__(
            inst_name="ESPRESSO" if not is_solar else "POET",
            array_size={"S2D": (170, 9111), "S1D": (1, 443262)},
            file_path=file_path,
            frameID=frameID,
            KW_identifier="ESO",
            user_configs=user_configs,
            reject_subInstruments=reject_subInstruments,
            quiet_user_params=quiet_user_params,
            need_external_data_load=is_solar,
            is_solar=is_solar,
        )

        self.instrument_properties["wavelength_coverage"] = coverage
        self.instrument_properties["resolution"] = 140_000
        self.instrument_properties["EarthLocation"] = EarthLocation.of_site(
            "Cerro Paranal"
        )
        self.instrument_properties["is_drift_corrected"] = True

        # https://www.eso.org/sci/facilities/paranal/astroclimate/site.html
        self.instrument_properties["site_pressure"] = 750

        self.UT_number = None

    def _load_ESO_DRS_KWs(self, header):
        super()._load_ESO_DRS_KWs(header)

        if header["ESO INS MODE"] == "SINGLEUHR":
            # Update instrumental resolution of ESPRESSO for UHR mode
            self.instrument_properties["resolution"] = 220_000

    def load_telemetry_info(self, header):
        # Find the UT number and load the airmass
        for i in range(1, 6):
            try:
                self.observation_info["airmass"] = header[
                    f"HIERARCH ESO TEL{i} AIRM START"
                ]
                self.UT_number = i
                break
            except KeyError as e:
                if i == 5:
                    msg = "\tCannot find ESO TELx AIRM START key"
                    raise KeyError(msg) from e

        # Environmental KWs for telfit (also needs airmassm previously loaded)
        ambi_KWs = {
            "relative_humidity": "AMBI RHUM",
            "ambient_temperature": "AMBI TEMP",
            "seeing": "AMBI FWHM START",
        }

        for name, endKW in ambi_KWs.items():
            try:
                self.observation_info[name] = float(
                    header[f"HIERARCH ESO TEL{self.UT_number} {endKW}"]
                )
                if "temperature" in name:  # store temperature in KELVIN for TELFIT
                    self.observation_info[name] = convert_temperature(
                        self.observation_info[name],
                        old_scale="Celsius",
                        new_scale="Kelvin",
                    )
            except KeyError:
                logger.info(f"Failed to load {name} from header")

        self.observation_info["DET_BINX"] = header["HIERARCH ESO DET BINX"]
        self.observation_info["DET_BINY"] = header["HIERARCH ESO DET BINY"]
        self.observation_info["FIBER_A"] = header[f"ESO INS5 LSELA NAME"]
        self.observation_info["FIBER_B"] = header[f"ESO INS5 LSELB NAME"]

    def check_header_QC_ESO_DRS(self, header):
        nonfatal_QC_flags = {
            "HIERARCH ESO INS5 ADC{} RA": 0,  # related with ADC2 problem
            "HIERARCH ESO INS5 ADC{} dec": 0,  # related with ADC2 problem
            "HIERARCH ESO INS5 ADC{} SENS1": 0,  # related with ADC2 problem
            "HIERARCH ESO INS5 ADC{} TEMP": 0,  # related with ADC2 problem
        }
        found_ADC_issue = False
        for flag, bad_value in nonfatal_QC_flags.items():
            found_UT = False
            try:
                for ADC in [1, 2]:
                    ADC_KW = flag.format(self.UT_number, ADC)
                    if header[ADC_KW] == bad_value:
                        msg = f"QC flag {ADC_KW} has a value of {bad_value}"
                        logger.warning(msg)
                        self._status.store_warning(KW_WARNING(msg))
                        found_ADC_issue = True
                    found_UT = True
            except:
                pass
            if not found_UT:
                logger.critical(
                    f"Did not find the entry for the following UT related metric: {flag}"
                )

        if found_ADC_issue:
            self._status.store_warning(KW_WARNING("ADC2 issues found"))

        super().check_header_QC_ESO_DRS(header)

    def build_mask(self, bypass_QualCheck: bool = False) -> None:
        super().build_mask(bypass_QualCheck=bypass_QualCheck, assess_bad_orders=False)

        if self.spectral_format == "S2D":
            # the first two orders of the RED CCD have a large amount of noise in the beginning so we remove a
            # portion from the start of those two orders Now, what is going on: we want to find the indexes,
            # from order 90 and 91 that are below the 5230 \AA
            inds = np.where(self.wavelengths[90:92, :] <= 5230)
            # numpy where returns the indexes assuming the zero to be 90 and the one to be 91. Remember that we
            # sliced the array to only remove from those two orders
            inds_1 = np.where(inds[0], 91, 90)
            # rebuild the 'numpy where' output, to pass to the mask as the proper output
            inds = (inds_1, inds[1])
            self.spectral_mask.add_indexes_to_mask(inds, ERROR_THRESHOLD)

        self.assess_bad_orders()

    def load_S1D_data(self) -> None:
        # Store the open status before calling parent classes
        is_open = self.is_open
        super().load_S1D_data()
        if self.is_solar and not is_open:
            self.apply_poet_data_corrections()

    def load_S2D_data(self) -> None:
        # Store the open status before calling parent classes
        is_open = self.is_open
        super().load_S2D_data()
        if self.is_solar and not is_open:
            self.apply_poet_data_corrections()

    def apply_poet_data_corrections(self) -> None:
        # Correct spectra
        self.remove_BERV_correction(self.wrong_berv)
        self.apply_BERV_correction(self.observation_info["BERV"])

    def trigger_data_storage(self, *args, **kwargs):
        super().trigger_data_storage(*args, **kwargs)

sub_instruments = {'ESPRESSO18': datetime.datetime.strptime('2019-06-27', '%Y-%m-%d'), 'ESPRESSO19': datetime.datetime.max} class-attribute instance-attribute

order_intervals: dict[DETECTOR_DEFINITION, slice] = {DETECTOR_DEFINITION.WHITE_LIGHT: list(range(170)), DETECTOR_DEFINITION.RED_DET: list(range(90, 170)), DETECTOR_DEFINITION.BLUE_DET: list(range(0, 90))} class-attribute instance-attribute

UT_number = None instance-attribute

load_telemetry_info(header)

Source code in src/ASTRA/Instruments/ESPRESSO.py 
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
def load_telemetry_info(self, header):
    # Find the UT number and load the airmass
    for i in range(1, 6):
        try:
            self.observation_info["airmass"] = header[
                f"HIERARCH ESO TEL{i} AIRM START"
            ]
            self.UT_number = i
            break
        except KeyError as e:
            if i == 5:
                msg = "\tCannot find ESO TELx AIRM START key"
                raise KeyError(msg) from e

    # Environmental KWs for telfit (also needs airmassm previously loaded)
    ambi_KWs = {
        "relative_humidity": "AMBI RHUM",
        "ambient_temperature": "AMBI TEMP",
        "seeing": "AMBI FWHM START",
    }

    for name, endKW in ambi_KWs.items():
        try:
            self.observation_info[name] = float(
                header[f"HIERARCH ESO TEL{self.UT_number} {endKW}"]
            )
            if "temperature" in name:  # store temperature in KELVIN for TELFIT
                self.observation_info[name] = convert_temperature(
                    self.observation_info[name],
                    old_scale="Celsius",
                    new_scale="Kelvin",
                )
        except KeyError:
            logger.info(f"Failed to load {name} from header")

    self.observation_info["DET_BINX"] = header["HIERARCH ESO DET BINX"]
    self.observation_info["DET_BINY"] = header["HIERARCH ESO DET BINY"]
    self.observation_info["FIBER_A"] = header[f"ESO INS5 LSELA NAME"]
    self.observation_info["FIBER_B"] = header[f"ESO INS5 LSELB NAME"]

check_header_QC_ESO_DRS(header)

Source code in src/ASTRA/Instruments/ESPRESSO.py 
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
def check_header_QC_ESO_DRS(self, header):
    nonfatal_QC_flags = {
        "HIERARCH ESO INS5 ADC{} RA": 0,  # related with ADC2 problem
        "HIERARCH ESO INS5 ADC{} dec": 0,  # related with ADC2 problem
        "HIERARCH ESO INS5 ADC{} SENS1": 0,  # related with ADC2 problem
        "HIERARCH ESO INS5 ADC{} TEMP": 0,  # related with ADC2 problem
    }
    found_ADC_issue = False
    for flag, bad_value in nonfatal_QC_flags.items():
        found_UT = False
        try:
            for ADC in [1, 2]:
                ADC_KW = flag.format(self.UT_number, ADC)
                if header[ADC_KW] == bad_value:
                    msg = f"QC flag {ADC_KW} has a value of {bad_value}"
                    logger.warning(msg)
                    self._status.store_warning(KW_WARNING(msg))
                    found_ADC_issue = True
                found_UT = True
        except:
            pass
        if not found_UT:
            logger.critical(
                f"Did not find the entry for the following UT related metric: {flag}"
            )

    if found_ADC_issue:
        self._status.store_warning(KW_WARNING("ADC2 issues found"))

    super().check_header_QC_ESO_DRS(header)

build_mask(bypass_QualCheck: bool = False) -> None

Source code in src/ASTRA/Instruments/ESPRESSO.py 
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
def build_mask(self, bypass_QualCheck: bool = False) -> None:
    super().build_mask(bypass_QualCheck=bypass_QualCheck, assess_bad_orders=False)

    if self.spectral_format == "S2D":
        # the first two orders of the RED CCD have a large amount of noise in the beginning so we remove a
        # portion from the start of those two orders Now, what is going on: we want to find the indexes,
        # from order 90 and 91 that are below the 5230 \AA
        inds = np.where(self.wavelengths[90:92, :] <= 5230)
        # numpy where returns the indexes assuming the zero to be 90 and the one to be 91. Remember that we
        # sliced the array to only remove from those two orders
        inds_1 = np.where(inds[0], 91, 90)
        # rebuild the 'numpy where' output, to pass to the mask as the proper output
        inds = (inds_1, inds[1])
        self.spectral_mask.add_indexes_to_mask(inds, ERROR_THRESHOLD)

    self.assess_bad_orders()

load_S1D_data() -> None

Source code in src/ASTRA/Instruments/ESPRESSO.py 
199
200
201
202
203
204
def load_S1D_data(self) -> None:
    # Store the open status before calling parent classes
    is_open = self.is_open
    super().load_S1D_data()
    if self.is_solar and not is_open:
        self.apply_poet_data_corrections()

load_S2D_data() -> None

Source code in src/ASTRA/Instruments/ESPRESSO.py 
206
207
208
209
210
211
def load_S2D_data(self) -> None:
    # Store the open status before calling parent classes
    is_open = self.is_open
    super().load_S2D_data()
    if self.is_solar and not is_open:
        self.apply_poet_data_corrections()

apply_poet_data_corrections() -> None

Source code in src/ASTRA/Instruments/ESPRESSO.py 
213
214
215
216
def apply_poet_data_corrections(self) -> None:
    # Correct spectra
    self.remove_BERV_correction(self.wrong_berv)
    self.apply_BERV_correction(self.observation_info["BERV"])

trigger_data_storage(*args, **kwargs)

Source code in src/ASTRA/Instruments/ESPRESSO.py 
218
219
def trigger_data_storage(self, *args, **kwargs):
    super().trigger_data_storage(*args, **kwargs)

HARPS

HARPS

Bases: ESO_PIPELINE

Interface to handle HARPS data; S1D not supported.

This class also defines 2 sub-Instruments:

  • HARPSpre - Before 2015-05-29
  • HARPSpost - Until the ends of time (hopefully)

The steps to load the S2D data are described in the HARPS DRS manual <https://www.eso.org/sci/facilities/lasilla/instruments/harps/doc/DRS.pdf>_. The summary is:

- Construct the wavelength solution & correct from BERV
- Load instrumental drift
- Construct flux noises:

    - Table 10 of `the user manual <https://www.eso.org/sci/facilities/lasilla/instruments/harps/doc/manual/HARPS-UserManual2.4.pdf>`_ gives max RON of 7.07 for red detector
    - Noise = sqrt(obj + sky + n*dark*expTime + nBinY*ron^2)

User parameters:

Currently there are no HARPS-specific parameters

Note: Check the User parameters of the parent classes for further customization options of SBART

Source code in src/ASTRA/Instruments/HARPS.py 
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
class HARPS(ESO_PIPELINE):
    """Interface to handle HARPS data; S1D **not** supported.

    This class also defines 2 sub-Instruments:

    * HARPSpre - Before  2015-05-29
    * HARPSpost - Until the ends of time (hopefully)

    The steps to load the S2D data are described in the HARPS `DRS manual <https://www.eso.org/sci/facilities/lasilla/instruments/harps/doc/DRS.pdf>`_. The summary is:

        - Construct the wavelength solution & correct from BERV
        - Load instrumental drift
        - Construct flux noises:

            - Table 10 of `the user manual <https://www.eso.org/sci/facilities/lasilla/instruments/harps/doc/manual/HARPS-UserManual2.4.pdf>`_ gives max RON of 7.07 for red detector
            - Noise = sqrt(obj + sky + n*dark*expTime + nBinY*ron^2)

    **User parameters:**

     Currently there are no HARPS-specific parameters

    *Note:* Check the **User parameters** of the parent classes for further customization options of SBART

    """

    sub_instruments = {
        "HARPS03": datetime.datetime.strptime("2015-05-29", r"%Y-%m-%d"),
        "HARPS15": datetime.datetime.max,
    }

    _name = "HARPS"
    _default_params = ESO_PIPELINE._default_params

    order_intervals: dict[DETECTOR_DEFINITION, slice] = {
        DETECTOR_DEFINITION.WHITE_LIGHT: list(range(0, 71)),
        DETECTOR_DEFINITION.RED_DET: list(range(47, 71)),
        DETECTOR_DEFINITION.BLUE_DET: list(range(0, 47)),
    }

    def __init__(
        self,
        file_path,
        user_configs: Optional[Dict[str, Any]] = None,
        reject_subInstruments=None,
        frameID=None,
        quiet_user_params: bool = True,
    ):
        """HARPS constructor.

        Parameters
        ----------
        file_path
            Path to the S2D (or S1D) file.
        user_configs
            Dictionary whose keys are the configurable options of ESPRESSO (check above)
        reject_subInstruments
            Iterable of subInstruments to fully reject
        frameID
            ID for this observation. Only used for organization purposes by :class:`~SBART.data_objects.DataClass`

        """
        logger.info(f"Creating frame from: {file_path}")

        if "use_old_pipeline" not in user_configs:
            # We need to do this, since we need to override the array size and header keywords
            user_configs["use_old_pipeline"] = self._default_params["use_old_pipeline"].default_value

        coverage = [350, 700]
        # Note: 46 blue orders and 26 red orders. From Table 2.2 of:
        # https://www.eso.org/sci/facilities/lasilla/instruments/harps/doc/manual/HARPS-UserManual2.4.pdf
        if user_configs["use_old_pipeline"]:
            mat_size = (72, 4096)
            KW_map = {
                "OBJECT": "OBJECT",
                "BJD": "HIERARCH ESO DRS BJD",
                "MJD": "MJD-OBS",
                "ISO-DATE": "DATE-OBS",
                "DRS-VERSION": "HIERARCH ESO DRS VERSION",
                "RA": "RA",
                "DEC": "DEC",
                "MD5-CHECK": "DATASUM",
                "SPEC_TYPE": None,
                "DRS_CCF_MASK": None,
                "DRS_FLUX_CORRECTION_TEMPLATE": None,
            }
        else:
            mat_size = (71, 4096)
            KW_map = {}

        if user_configs["use_old_pipeline"]:
            file_path, self.ccf_path, self.BIS_file, search_status = self.find_files(file_path)

        super().__init__(
            inst_name="HARPS",
            array_size={"S2D": mat_size},
            file_path=file_path,
            KW_identifier="ESO",
            frameID=frameID,
            user_configs=user_configs,
            reject_subInstruments=reject_subInstruments,
            quiet_user_params=quiet_user_params,
            override_KW_map=KW_map,
            override_indicators=("CONTRAST", "FWHM"),
        )

        if user_configs["use_old_pipeline"] and not search_status.is_good_flag:
            self.add_to_status(search_status)

        self.instrument_properties["wavelength_coverage"] = coverage

        if user_configs["use_old_pipeline"]:
            self.instrument_properties["is_drift_corrected"] = False
            self.is_BERV_corrected = False

        self.instrument_properties["resolution"] = 115_000
        self.instrument_properties["EarthLocation"] = EarthLocation.of_site("La Silla Observatory")
        # ? same as for Paranal?
        # https://www.eso.org/sci/facilities/paranal/astroclimate/site.html
        self.instrument_properties["site_pressure"] = 750
        self.is_blaze_corrected = False

    def load_telemetry_info(self, header):
        """Loads (at least) the following keywords:

        - relative humidity
        - ambient temperature, in Celsius
        - airmass
        - Detector

        Parameters
        ----------
        header

        Returns
        -------

        """
        ambi_KWs = {
            "relative_humidity": "RHUM",
            "ambient_temperature": "TEMP",
        }

        for name, endKW in ambi_KWs.items():
            self.observation_info[name] = header[f"HIERARCH {self.KW_identifier} TEL AMBI {endKW}"]
            if "temperature" in name:  # store temperature in KELVIN for TELFIT
                self.observation_info[name] = convert_temperature(
                    self.observation_info[name],
                    old_scale="Celsius",
                    new_scale="Kelvin",
                )

        if self.observation_info["relative_humidity"] == 255:
            logger.warning(f"{self.name} has an invalid value in the humidity sensor...")
            self.observation_info["relative_humidity"] = np.nan

        self.observation_info["airmass"] = header["HIERARCH ESO TEL AIRM START"]

    def find_files(self, file_name):
        """Find the CCF and S2D files and BIS files, which should be stored inside the same folder"""
        logger.debug("Searching for the ccf and e2ds files")

        search_status = MISSING_DATA("Missing the ccf file")
        ccf_path = None
        bis_path = None

        if os.path.isdir(file_name):
            logger.debug("Received a folder, searching inside for necessary files")
            # search for e2ds file
            folder_name = file_name

            e2ds_files = glob.glob(os.path.join(folder_name, "*e2ds_A.fits"), recursive=True)
            ccf_files = glob.glob(os.path.join(folder_name, "*ccf*A.fits"), recursive=True)
            bis_files = glob.glob(os.path.join(folder_name, "*bis*A.fits"), recursive=True)

            for name, elems in [
                ("e2ds_A", e2ds_files),
                ("ccf", ccf_files),
                ("bis", bis_files),
            ]:
                if len(elems) > 1:
                    msg = f"HARPS data only received folder name and it has more than 1 {name} file in it"
                    raise custom_exceptions.InvalidConfiguration(msg)

                if len(elems) < 1:
                    msg = f"HARPS data only received folder name and it has no {name} file in it"
                    if name != "bis":
                        # The BIS file is not critical for the run
                        raise custom_exceptions.InvalidConfiguration(msg)
                    logger.critical(msg)

            e2ds_path = e2ds_files[0]
            ccf_path = ccf_files[0]
            bis_path = bis_files[0]
            search_status = SUCCESS("Found all input files")
        else:
            logger.debug("Received path of E2DS file; searching for CCF with matching name")
            folder_name = os.path.dirname(file_name)
            e2ds_path = file_name
            file_start, *_ = os.path.basename(file_name).split("_")

            found_CCF = False
            found_BIS = False
            bis_files = glob.glob(os.path.join(folder_name, "*bis*A.fits"), recursive=True)
            ccf_files = Path(folder_name).glob("*ccf*A.fits")

            for file in ccf_files:
                if file_start in file.as_posix():
                    ccf_path = file
                    found_CCF = True
            for file in bis_files:
                if file_start in file:
                    bis_path = file
                    found_BIS = True
            if found_CCF:
                logger.info(f"Found CCF file: {ccf_path}")
                search_status = SUCCESS("Found CCF file")
            else:
                logger.critical("Was not able to find CCF file. Marking frame as invalid")
                ccf_path = ""

            if not found_BIS:
                bis_path = None

        return e2ds_path, ccf_path, bis_path, search_status

    def build_HARPS_wavelengths(self, hdr):
        """Compute the wavelength solution to this given spectra (EQ 4.1 of DRS manual)
        Convert from air wavelenbgths to vacuum
        """
        # degree of the polynomial
        d = hdr["HIERARCH ESO DRS CAL TH DEG LL"]
        # number of orders
        omax = hdr.get("HIERARCH ESO DRS CAL LOC NBO", self.array_size[0])
        xmax = self.array_size[1]

        # matrix X:
        #
        # axis 0: the entry corresponding to each coefficient
        # axis 1: each pixel number

        x = np.empty((d + 1, xmax), "int64")
        x[0].fill(1)  # x[0,*] = x^0 = 1,1,1,1,1,...
        x[1] = np.arange(xmax)

        for i in range(1, d):
            x[i + 1] = x[i] * x[1]

        # matrix A:
        #
        # axis 0: the different orders
        # axis 1: all coefficients for the given order

        A = np.reshape(
            [hdr["HIERARCH ESO DRS CAL TH COEFF LL" + str(i)] for i in range(omax * (d + 1))],
            (omax, d + 1),
        )  # slow 30 ms

        # the wavelengths for each order are a simple dot product between the coefficients and pixel-wise data (X)
        wavelengths = np.dot(A, x)

        vacuum_wavelengths = airtovac(wavelengths)
        return vacuum_wavelengths

    def _load_old_DRS_KWs(self, header):
        if not self._internal_configs["use_old_pipeline"]:
            raise custom_exceptions.InvalidConfiguration("Can't load data from old pipeline without the config")

        self.observation_info["MAX_BERV"] = header["HIERARCH ESO DRS BERVMX"] * kilometer_second
        self.observation_info["BERV"] = header["HIERARCH ESO DRS BERV"] * kilometer_second

        # Environmental KWs for telfit (also needs airmassm previously loaded)
        ambi_KWs = {
            "relative_humidity": "AMBI RHUM",
            "ambient_temperature": "AMBI TEMP",
        }

        for name, endKW in ambi_KWs.items():
            self.observation_info[name] = header[f"HIERARCH ESO TEL {endKW}"]
            if "temperature" in name:  # store temperature in KELVIN for TELFIT
                self.observation_info[name] = convert_temperature(
                    self.observation_info[name],
                    old_scale="Celsius",
                    new_scale="Kelvin",
                )

        for order in range(self.N_orders):
            self.observation_info["orderwise_SNRs"].append(header[f"HIERARCH ESO DRS SPE EXT SN{order}"])

        self.observation_info["airmass"] = header["HIERARCH ESO TEL AIRM START"]

        bad_drift = False
        try:
            flag = "HIERARCH ESO DRS DRIFT QC"
            if header[flag].strip() != "PASSED":
                bad_drift = True
                msg = f"QC flag {flag} meets the bad value"
                logger.warning(msg)
                self._status.store_warning(KW_WARNING(msg))
            else:
                # self.logger.info("DRIFT QC has passed")
                drift = header["HIERARCH ESO DRS DRIFT RV USED"] * meter_second
                drift_err = header["HIERARCH ESO DRS DRIFT NOISE"] * meter_second
        except Exception:
            bad_drift = True
            logger.warning("DRIFT KW does not exist")

        if bad_drift:
            logger.warning("Due to previous drift-related problems, setting it to zero [m/s]")
            drift = 0 * meter_second
            drift_err = 0 * meter_second

        self.observation_info["drift"] = drift
        self.observation_info["drift_ERR"] = drift_err
        self.load_ccf_data()

    def load_ccf_data(self) -> None:
        """Load the necessarfy CCF data from the file!"""
        logger.debug("Loading data from the ccf file")
        header = fits.getheader(self.ccf_path)

        self.observation_info["DRS_RV"] = header["HIERARCH ESO DRS CCF RV"] * kilometer_second
        self.observation_info["SPEC_TYPE"] = header["HIERARCH ESO DRS CCF MASK"]

        RV_err = np.sqrt(
            header["HIERARCH ESO DRS CAL TH ERROR"] ** 2
            +
            # hdulist[0].header['HIERARCH ESO DRS DRIFT NOISE']**2   +
            (1000 * header["HIERARCH ESO DRS CCF NOISE"]) ** 2,
        )
        self.observation_info["DRS_RV_ERR"] = RV_err * meter_second

        for key in self.available_indicators:
            full_key = "HIERARCH ESO DRS CCF " + key
            self.observation_info[key] = header[full_key]

        # We are missing error in CONTRAST!
        self.observation_info["FWHM_ERR"] = convert_data(
            2.35 * RV_err * meter_second,
            new_units=kilometer_second,
            as_value=True,
        )
        self.observation_info["BIS SPAN_ERR"] = convert_data(
            np.sqrt(2) * RV_err * meter_second,
            new_units=kilometer_second,
            as_value=True,
        )

        if self.BIS_file is not None:
            head = fits.getheader(self.BIS_file)
            self.observation_info["BIS SPAN"] = head["HIERARCH ESO DRS BIS SPAN"]

    def load_S1D_data(self) -> Mask:
        raise NotImplementedError

    def check_header_QC_old_DRS(self, header):
        logger.info("Currently missing QC checks for the old DRS")

    def load_old_DRS_S2D(self):
        """Loads the spectra

        Returns
        -------

        """
        super().load_S2D_data()

        with fits.open(self.file_path) as hdulist:
            # Compute the wavelength solution + BERV correction
            wave_from_file = self.build_HARPS_wavelengths(hdulist[0].header)

            sci_data = hdulist[0].data  # spetra from all orders

            # photon noise + estimate of max value for the rest
            # from ETC calculator the readout noise should be the largest contribution
            # assuming that it is of ~7e- (equal to manual) it should have a maximum contribution
            # of 200
            flux_errors = np.sqrt(250 + np.abs(sci_data, dtype=float))

            # Validate for overflows and missing data
            quality_data = np.zeros(sci_data.shape)
            quality_data[np.where(np.isnan(sci_data))] = NAN_DATA.code
            quality_data[np.where(sci_data > 300000)] = SATURATION.code
            quality_data[np.where(sci_data < -3 * flux_errors)] = ERROR_THRESHOLD.code

        self.spectra = sci_data.astype(np.float64)
        self.wavelengths = wave_from_file
        self.qual_data = quality_data
        self.uncertainties = flux_errors.astype(np.float64)

        self.build_mask(bypass_QualCheck=False)
        return 1

    def close_arrays(self):
        super().close_arrays()
        if self._internal_configs["use_old_pipeline"]:
            self.is_BERV_corrected = False

sub_instruments = {'HARPS03': datetime.datetime.strptime('2015-05-29', '%Y-%m-%d'), 'HARPS15': datetime.datetime.max} class-attribute instance-attribute

order_intervals: dict[DETECTOR_DEFINITION, slice] = {DETECTOR_DEFINITION.WHITE_LIGHT: list(range(0, 71)), DETECTOR_DEFINITION.RED_DET: list(range(47, 71)), DETECTOR_DEFINITION.BLUE_DET: list(range(0, 47))} class-attribute instance-attribute

is_BERV_corrected = False instance-attribute

is_blaze_corrected = False instance-attribute

load_telemetry_info(header)

Loads (at least) the following keywords:

  • relative humidity
  • ambient temperature, in Celsius
  • airmass
  • Detector

Parameters:

Name Type Description Default
header
required
Source code in src/ASTRA/Instruments/HARPS.py 
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
def load_telemetry_info(self, header):
    """Loads (at least) the following keywords:

    - relative humidity
    - ambient temperature, in Celsius
    - airmass
    - Detector

    Parameters
    ----------
    header

    Returns
    -------

    """
    ambi_KWs = {
        "relative_humidity": "RHUM",
        "ambient_temperature": "TEMP",
    }

    for name, endKW in ambi_KWs.items():
        self.observation_info[name] = header[f"HIERARCH {self.KW_identifier} TEL AMBI {endKW}"]
        if "temperature" in name:  # store temperature in KELVIN for TELFIT
            self.observation_info[name] = convert_temperature(
                self.observation_info[name],
                old_scale="Celsius",
                new_scale="Kelvin",
            )

    if self.observation_info["relative_humidity"] == 255:
        logger.warning(f"{self.name} has an invalid value in the humidity sensor...")
        self.observation_info["relative_humidity"] = np.nan

    self.observation_info["airmass"] = header["HIERARCH ESO TEL AIRM START"]

find_files(file_name)

Find the CCF and S2D files and BIS files, which should be stored inside the same folder

Source code in src/ASTRA/Instruments/HARPS.py 
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
def find_files(self, file_name):
    """Find the CCF and S2D files and BIS files, which should be stored inside the same folder"""
    logger.debug("Searching for the ccf and e2ds files")

    search_status = MISSING_DATA("Missing the ccf file")
    ccf_path = None
    bis_path = None

    if os.path.isdir(file_name):
        logger.debug("Received a folder, searching inside for necessary files")
        # search for e2ds file
        folder_name = file_name

        e2ds_files = glob.glob(os.path.join(folder_name, "*e2ds_A.fits"), recursive=True)
        ccf_files = glob.glob(os.path.join(folder_name, "*ccf*A.fits"), recursive=True)
        bis_files = glob.glob(os.path.join(folder_name, "*bis*A.fits"), recursive=True)

        for name, elems in [
            ("e2ds_A", e2ds_files),
            ("ccf", ccf_files),
            ("bis", bis_files),
        ]:
            if len(elems) > 1:
                msg = f"HARPS data only received folder name and it has more than 1 {name} file in it"
                raise custom_exceptions.InvalidConfiguration(msg)

            if len(elems) < 1:
                msg = f"HARPS data only received folder name and it has no {name} file in it"
                if name != "bis":
                    # The BIS file is not critical for the run
                    raise custom_exceptions.InvalidConfiguration(msg)
                logger.critical(msg)

        e2ds_path = e2ds_files[0]
        ccf_path = ccf_files[0]
        bis_path = bis_files[0]
        search_status = SUCCESS("Found all input files")
    else:
        logger.debug("Received path of E2DS file; searching for CCF with matching name")
        folder_name = os.path.dirname(file_name)
        e2ds_path = file_name
        file_start, *_ = os.path.basename(file_name).split("_")

        found_CCF = False
        found_BIS = False
        bis_files = glob.glob(os.path.join(folder_name, "*bis*A.fits"), recursive=True)
        ccf_files = Path(folder_name).glob("*ccf*A.fits")

        for file in ccf_files:
            if file_start in file.as_posix():
                ccf_path = file
                found_CCF = True
        for file in bis_files:
            if file_start in file:
                bis_path = file
                found_BIS = True
        if found_CCF:
            logger.info(f"Found CCF file: {ccf_path}")
            search_status = SUCCESS("Found CCF file")
        else:
            logger.critical("Was not able to find CCF file. Marking frame as invalid")
            ccf_path = ""

        if not found_BIS:
            bis_path = None

    return e2ds_path, ccf_path, bis_path, search_status

build_HARPS_wavelengths(hdr)

Compute the wavelength solution to this given spectra (EQ 4.1 of DRS manual) Convert from air wavelenbgths to vacuum

Source code in src/ASTRA/Instruments/HARPS.py 
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
def build_HARPS_wavelengths(self, hdr):
    """Compute the wavelength solution to this given spectra (EQ 4.1 of DRS manual)
    Convert from air wavelenbgths to vacuum
    """
    # degree of the polynomial
    d = hdr["HIERARCH ESO DRS CAL TH DEG LL"]
    # number of orders
    omax = hdr.get("HIERARCH ESO DRS CAL LOC NBO", self.array_size[0])
    xmax = self.array_size[1]

    # matrix X:
    #
    # axis 0: the entry corresponding to each coefficient
    # axis 1: each pixel number

    x = np.empty((d + 1, xmax), "int64")
    x[0].fill(1)  # x[0,*] = x^0 = 1,1,1,1,1,...
    x[1] = np.arange(xmax)

    for i in range(1, d):
        x[i + 1] = x[i] * x[1]

    # matrix A:
    #
    # axis 0: the different orders
    # axis 1: all coefficients for the given order

    A = np.reshape(
        [hdr["HIERARCH ESO DRS CAL TH COEFF LL" + str(i)] for i in range(omax * (d + 1))],
        (omax, d + 1),
    )  # slow 30 ms

    # the wavelengths for each order are a simple dot product between the coefficients and pixel-wise data (X)
    wavelengths = np.dot(A, x)

    vacuum_wavelengths = airtovac(wavelengths)
    return vacuum_wavelengths

load_ccf_data() -> None

Load the necessarfy CCF data from the file!

Source code in src/ASTRA/Instruments/HARPS.py 
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
def load_ccf_data(self) -> None:
    """Load the necessarfy CCF data from the file!"""
    logger.debug("Loading data from the ccf file")
    header = fits.getheader(self.ccf_path)

    self.observation_info["DRS_RV"] = header["HIERARCH ESO DRS CCF RV"] * kilometer_second
    self.observation_info["SPEC_TYPE"] = header["HIERARCH ESO DRS CCF MASK"]

    RV_err = np.sqrt(
        header["HIERARCH ESO DRS CAL TH ERROR"] ** 2
        +
        # hdulist[0].header['HIERARCH ESO DRS DRIFT NOISE']**2   +
        (1000 * header["HIERARCH ESO DRS CCF NOISE"]) ** 2,
    )
    self.observation_info["DRS_RV_ERR"] = RV_err * meter_second

    for key in self.available_indicators:
        full_key = "HIERARCH ESO DRS CCF " + key
        self.observation_info[key] = header[full_key]

    # We are missing error in CONTRAST!
    self.observation_info["FWHM_ERR"] = convert_data(
        2.35 * RV_err * meter_second,
        new_units=kilometer_second,
        as_value=True,
    )
    self.observation_info["BIS SPAN_ERR"] = convert_data(
        np.sqrt(2) * RV_err * meter_second,
        new_units=kilometer_second,
        as_value=True,
    )

    if self.BIS_file is not None:
        head = fits.getheader(self.BIS_file)
        self.observation_info["BIS SPAN"] = head["HIERARCH ESO DRS BIS SPAN"]

load_S1D_data() -> Mask

Source code in src/ASTRA/Instruments/HARPS.py 
381
382
def load_S1D_data(self) -> Mask:
    raise NotImplementedError

check_header_QC_old_DRS(header)

Source code in src/ASTRA/Instruments/HARPS.py 
384
385
def check_header_QC_old_DRS(self, header):
    logger.info("Currently missing QC checks for the old DRS")

load_old_DRS_S2D()

Loads the spectra

Source code in src/ASTRA/Instruments/HARPS.py 
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
def load_old_DRS_S2D(self):
    """Loads the spectra

    Returns
    -------

    """
    super().load_S2D_data()

    with fits.open(self.file_path) as hdulist:
        # Compute the wavelength solution + BERV correction
        wave_from_file = self.build_HARPS_wavelengths(hdulist[0].header)

        sci_data = hdulist[0].data  # spetra from all orders

        # photon noise + estimate of max value for the rest
        # from ETC calculator the readout noise should be the largest contribution
        # assuming that it is of ~7e- (equal to manual) it should have a maximum contribution
        # of 200
        flux_errors = np.sqrt(250 + np.abs(sci_data, dtype=float))

        # Validate for overflows and missing data
        quality_data = np.zeros(sci_data.shape)
        quality_data[np.where(np.isnan(sci_data))] = NAN_DATA.code
        quality_data[np.where(sci_data > 300000)] = SATURATION.code
        quality_data[np.where(sci_data < -3 * flux_errors)] = ERROR_THRESHOLD.code

    self.spectra = sci_data.astype(np.float64)
    self.wavelengths = wave_from_file
    self.qual_data = quality_data
    self.uncertainties = flux_errors.astype(np.float64)

    self.build_mask(bypass_QualCheck=False)
    return 1

close_arrays()

Source code in src/ASTRA/Instruments/HARPS.py 
422
423
424
425
def close_arrays(self):
    super().close_arrays()
    if self._internal_configs["use_old_pipeline"]:
        self.is_BERV_corrected = False

HARPSN

HARPSN

Bases: ESO_PIPELINE

Interface to handle HARPSN data; S1D not supported!

This class also defines 2 sub-Instruments:

  • HARPSN - Until the ends of time (hopefully)

The steps to load the S2D data are described in the HARPS DRS manual <https://www.eso.org/sci/facilities/lasilla/instruments/harps/doc/DRS.pdf>_. The summary is:

- Construct the wavelength solution & correct from BERV
- Load instrumental drift
- Construct flux noises:

    - Table 10 of `the user manual <https://www.eso.org/sci/facilities/lasilla/instruments/harps/doc/manual/HARPS-UserManual2.4.pdf>`_ gives max RON of 7.07 for red detector
    - Noise = sqrt(obj + sky + n*dark*expTime + nBinY*ron^2)

User parameters:

Currently there are no HARPSN-specific parameters

Note: Check the User parameters of the parent classes for further customization options of SBART

Source code in src/ASTRA/Instruments/HARPSN.py 
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
class HARPSN(ESO_PIPELINE):
    """
    Interface to handle HARPSN data; S1D **not** supported!

    This class also defines 2 sub-Instruments:

    * HARPSN - Until the ends of time (hopefully)


    The steps to load the S2D data are described in the HARPS `DRS manual <https://www.eso.org/sci/facilities/lasilla/instruments/harps/doc/DRS.pdf>`_. The summary is:

        - Construct the wavelength solution & correct from BERV
        - Load instrumental drift
        - Construct flux noises:

            - Table 10 of `the user manual <https://www.eso.org/sci/facilities/lasilla/instruments/harps/doc/manual/HARPS-UserManual2.4.pdf>`_ gives max RON of 7.07 for red detector
            - Noise = sqrt(obj + sky + n*dark*expTime + nBinY*ron^2)

    **User parameters:**

     Currently there are no HARPSN-specific parameters

    *Note:* Check the **User parameters** of the parent classes for further customization options of SBART

    """

    sub_instruments = {
        # "HARPSpre": datetime.datetime.strptime("2015-05-29", r"%Y-%m-%d"),
        "HARPSN": datetime.datetime.max,
    }

    _name = "HARPSN"
    _default_params = ESO_PIPELINE._default_params

    def __init__(
        self,
        file_path,
        user_configs: Optional[Dict[str, Any]] = None,
        reject_subInstruments=None,
        frameID=None,
        quiet_user_params: bool = True,
    ):
        """

        Parameters
        ----------
        file_path
            Path to the S2D (or S1D) file.
        user_configs
            Dictionary whose keys are the configurable options of ESPRESSO (check above)
        reject_subInstruments
            Iterable of subInstruments to fully reject
        frameID
            ID for this observation. Only used for organization purposes by :class:`~SBART.data_objects.DataClass`
        """

        logger.info("Creating frame from: {}".format(file_path))

        coverage = [390, 700]
        search_status = SUCCESS
        if user_configs.get("use_old_pipeline", False):  # For the old pipeline!
            override_KW_map = {
                "OBJECT": "HIERARCH TNG OBS TARG NAME",  # "OBJECT",
                "BJD": "HIERARCH TNG DRS BJD",
                "MJD": "MJD-OBS",
                "ISO-DATE": "DATE-OBS",
                "DRS-VERSION": "HIERARCH TNG DRS VERSION",
                "RA": "RA",
                "DEC": "DEC",
                "SPEC_TYPE": "HIERARCH TNG TEL TARG SPTYPE",
                "MD5-CHECK": "HIERARCH TNG DRS BJD",  # Missing the MD5 on the old pipe, so this is a stopgap
            }
            file_path, self.ccf_path, search_status = self.find_files(file_path)
            override_indicators = ("CONTRAST", "FWHM")

            # The flag for the BERV correction is set at the time of loading the S2D arrays!

        else:  # For the new pipeline
            override_KW_map = {"OBJECT": "HIERARCH TNG OBS TARG NAME"}
            override_indicators = None

        super().__init__(
            inst_name="HARPSN",
            array_size={"S2D": (69, 4096)},
            file_path=file_path,
            KW_identifier="TNG",
            frameID=frameID,
            override_KW_map=override_KW_map,
            user_configs=user_configs,
            reject_subInstruments=reject_subInstruments,
            quiet_user_params=quiet_user_params,
            override_indicators=override_indicators,
        )

        if (
            user_configs.get("use_old_pipeline", False)
            and not search_status.is_good_flag
        ):
            self.add_to_status(search_status)

        self.instrument_properties["wavelength_coverage"] = coverage

        self.instrument_properties["resolution"] = 115_000
        self.instrument_properties["EarthLocation"] = EarthLocation.of_site(
            "Roque de los Muchachos"
        )

        # https://tngweb.tng.iac.es/weather/current
        self.instrument_properties["site_pressure"] = 770

    def check_header_QC_old_DRS(self, header):
        logger.info("Currently missing QC checks for the old DRS")

    def _load_old_DRS_KWs(self, header):
        if not self._internal_configs["use_old_pipeline"]:
            raise custom_exceptions.InvalidConfiguration(
                "Can't load data from old pipeline without the config"
            )

        self.observation_info["MAX_BERV"] = (
            header["HIERARCH TNG DRS BERVMX"] * kilometer_second
        )
        self.observation_info["BERV"] = (
            header["HIERARCH TNG DRS BERV"] * kilometer_second
        )

        # Environmental KWs for telfit (also needs airmassm previously loaded)

        ambi_KWs = {
            "relative_humidity": "HUMIDITY",
            "ambient_temperature": "TEMP10M",
        }

        for name, endKW in ambi_KWs.items():
            self.observation_info[name] = header[f"HIERARCH TNG METEO {endKW}"]
            if "temperature" in name:  # store temperature in KELVIN for TELFIT
                self.observation_info[name] = convert_temperature(
                    self.observation_info[name], old_scale="Celsius", new_scale="Kelvin"
                )

        if self.observation_info["relative_humidity"] == 255:
            logger.warning(
                f"{self.name} has an invalid value in the humidity sensor..."
            )
            self.observation_info["relative_humidity"] = np.nan

        for order in range(self.N_orders):
            self.observation_info["orderwise_SNRs"].append(
                header[f"HIERARCH TNG DRS SPE EXT SN{order}"]
            )

        self.observation_info["airmass"] = header["AIRMASS"]

        bad_drift = False
        try:
            flag = "HIERARCH TNG DRS DRIFT QC"
            if header[flag].strip() != "PASSED":
                bad_drift = True
                msg = f"QC flag {flag} meets the bad value"
                logger.warning(msg)
                self._status.store_warning(KW_WARNING(msg))
            else:
                # self.logger.info("DRIFT QC has passed")
                drift = header["HIERARCH TNG DRS DRIFT RV USED"] * meter_second
                drift_err = header["HIERARCH TNG DRS DRIFT NOISE"] * meter_second
        except Exception as e:
            bad_drift = True
            logger.warning("DRIFT KW does not exist")

        if bad_drift:
            logger.warning(
                "Due to previous drift-related problems, setting it to zero [m/s]"
            )
            drift = 0 * meter_second
            drift_err = 0 * meter_second

        self.observation_info["drift"] = drift
        self.observation_info["drift_ERR"] = drift_err
        self._load_ccf_data()

    def _load_ccf_data(self) -> None:
        """
        Load the necessarfy CCF data from the file!
        """
        if not self._internal_configs["use_old_pipeline"]:
            raise custom_exceptions.InvalidConfiguration(
                "Can't load data from old pipeline without the config"
            )

        logger.debug("Loading data from the ccf file")
        header = fits.getheader(self.ccf_path)

        for key in self.available_indicators:
            full_key = "HIERARCH TNG DRS CCF " + key
            self.observation_info[key] = header[full_key]

        self.observation_info["DRS_RV"] = (
            header["HIERARCH TNG DRS CCF RV"] * kilometer_second
        )

        RV_err = header.get("HIERARCH TNG DRS DVRMS", None)
        if RV_err is None:
            logger.critical(
                "Couldn't find DRS error from the appropriate KW. Estimating it with photon noise!"
            )
            RV_err = np.sqrt(
                header["HIERARCH TNG DRS CAL TH ERROR"] ** 2
                +
                # hdulist[0].header['HIERARCH ESO DRS DRIFT NOISE']**2   +
                (1000 * header["HIERARCH TNG DRS CCF NOISE"]) ** 2
            )

        self.observation_info["DRS_RV_ERR"] = RV_err * meter_second

    def build_HARPS_wavelengths(self, hdr):
        """
        Compute the wavelength solution to this given spectra (EQ 4.1 of DRS manual)
        Convert from air wavelenbgths to vacuum
        """
        if not self._internal_configs["use_old_pipeline"]:
            raise custom_exceptions.InvalidConfiguration(
                "Can't construct wavelengths for new pipeline"
            )

        # degree of the polynomial
        d = hdr["HIERARCH TNG DRS CAL TH DEG LL"]
        # number of orders
        omax = hdr.get("HIERARCH TNG DRS CAL LOC NBO", self.array_size[0])
        xmax = self.array_size[1]

        # matrix X:
        #
        # axis 0: the entry corresponding to each coefficient
        # axis 1: each pixel number

        x = np.empty((d + 1, xmax), "int64")
        x[0].fill(1)  # x[0,*] = x^0 = 1,1,1,1,1,...
        x[1] = np.arange(xmax)

        for i in range(1, d):
            x[i + 1] = x[i] * x[1]

        # matrix A:
        #
        # axis 0: the different orders
        # axis 1: all coefficients for the given order

        A = np.reshape(
            [
                hdr["HIERARCH TNG DRS CAL TH COEFF LL" + str(i)]
                for i in range(omax * (d + 1))
            ],
            (omax, d + 1),
        )  # slow 30 ms

        # the wavelengths for each order are a simple dot product between the coefficients and pixel-wise data (X)
        wavelengths = np.dot(A, x)

        vacuum_wavelengths = airtovac(wavelengths)
        return vacuum_wavelengths

    def load_old_DRS_S2D(self):
        """
        load the data from the old HARPS-N pipeline. This will be mainly used for the comparison with the
        HARPS-TERRA pipeline
        """
        if not self._internal_configs["use_old_pipeline"]:
            raise custom_exceptions.InvalidConfiguration(
                "Can't load data from old pipeline without the config"
            )

        if self._internal_configs["use_old_pipeline"]:
            self.is_BERV_corrected = False

        with fits.open(self.file_path) as hdulist:
            # Compute the wavelength solution + BERV correction
            wave_from_file = self.build_HARPS_wavelengths(hdulist[0].header)

            sci_data = hdulist[0].data  # spetra from all orders

            # photon noise + estimate of max value for the rest
            # from ETC calculator the readout noise should be the largest contribution
            # assuming that it is of ~7e- (equal to manual) it should have a maximum contribution
            # of 200
            flux_errors = np.sqrt(250 + np.abs(sci_data, dtype=float))

            # Validate for overflows and missing data
            quality_data = np.zeros(sci_data.shape)
            quality_data[np.where(np.isnan(sci_data))] = NAN_DATA.code
            quality_data[np.where(sci_data > 300000)] = SATURATION.code
            quality_data[np.where(sci_data < -3 * flux_errors)] = ERROR_THRESHOLD.code

        self.spectra = sci_data.astype(np.float64)
        self.wavelengths = wave_from_file
        self.qual_data = quality_data
        self.uncertainties = flux_errors.astype(np.float64)

        self.build_mask(bypass_QualCheck=False)

    def close_arrays(self):
        """
        Reset the BERV correction flag if we are using the old pipeline version!
        Returns
        -------

        """
        super().close_arrays()
        if self._internal_configs["use_old_pipeline"]:
            self.is_BERV_corrected = False

    def find_files(self, file_name: UI_PATH):
        """
        Find the CCF and the S2D files, which should be stored inside the same folder
        """
        logger.debug("Searching for the ccf and e2ds files")

        if not isinstance(file_name, Path):
            file_name = Path(file_name)

        search_status = MISSING_DATA("Missing the ccf file")
        ccf_path = None

        if file_name.is_dir():
            logger.debug("Received a folder, searching inside for necessary files")
            # search for e2ds file
            folder_name = file_name

            e2ds_files = file_name.glob("**/*e2ds_A.fits")
            ccf_files = file_name.glob("**/*ccf_*_A.fits")

            for name, elems in [("e2ds_A", e2ds_files), ("ccf", ccf_files)]:
                if len(elems) > 1:
                    msg = f"HARPS data only received folder name and it has more than 1 {name} file in it"
                    raise custom_exceptions.InvalidConfiguration(msg)

                if len(elems) < 1:
                    msg = f"HARPS data only received folder name and it has no {name} file in it"
                    raise custom_exceptions.InvalidConfiguration(msg)

            e2ds_path = e2ds_files[0]
            ccf_path = ccf_files[0]
            search_status = SUCCESS("Found all input files")
        else:
            logger.debug(
                "Received path of E2DS file; searching for CCF with matching name"
            )
            folder_name = file_name.parent
            e2ds_path = file_name
            file_start, *_ = file_name.stem.split("_")

            found_CCF = False
            ccf_files = folder_name.glob("*ccf*_A.fits")

            for file in ccf_files:
                if file_start in file.name:
                    ccf_path = file
                    found_CCF = True

            if found_CCF:
                logger.info("Found CCF file: {}".format(ccf_path))
                search_status = SUCCESS("Found CCF file")
            else:
                logger.critical(
                    "Was not able to find CCF file. Marking frame as invalid"
                )
                ccf_path = ""

        return e2ds_path, ccf_path, search_status

sub_instruments = {'HARPSN': datetime.datetime.max} class-attribute instance-attribute

check_header_QC_old_DRS(header)

Source code in src/ASTRA/Instruments/HARPSN.py 
138
139
def check_header_QC_old_DRS(self, header):
    logger.info("Currently missing QC checks for the old DRS")

build_HARPS_wavelengths(hdr)

Compute the wavelength solution to this given spectra (EQ 4.1 of DRS manual) Convert from air wavelenbgths to vacuum

Source code in src/ASTRA/Instruments/HARPSN.py 
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
def build_HARPS_wavelengths(self, hdr):
    """
    Compute the wavelength solution to this given spectra (EQ 4.1 of DRS manual)
    Convert from air wavelenbgths to vacuum
    """
    if not self._internal_configs["use_old_pipeline"]:
        raise custom_exceptions.InvalidConfiguration(
            "Can't construct wavelengths for new pipeline"
        )

    # degree of the polynomial
    d = hdr["HIERARCH TNG DRS CAL TH DEG LL"]
    # number of orders
    omax = hdr.get("HIERARCH TNG DRS CAL LOC NBO", self.array_size[0])
    xmax = self.array_size[1]

    # matrix X:
    #
    # axis 0: the entry corresponding to each coefficient
    # axis 1: each pixel number

    x = np.empty((d + 1, xmax), "int64")
    x[0].fill(1)  # x[0,*] = x^0 = 1,1,1,1,1,...
    x[1] = np.arange(xmax)

    for i in range(1, d):
        x[i + 1] = x[i] * x[1]

    # matrix A:
    #
    # axis 0: the different orders
    # axis 1: all coefficients for the given order

    A = np.reshape(
        [
            hdr["HIERARCH TNG DRS CAL TH COEFF LL" + str(i)]
            for i in range(omax * (d + 1))
        ],
        (omax, d + 1),
    )  # slow 30 ms

    # the wavelengths for each order are a simple dot product between the coefficients and pixel-wise data (X)
    wavelengths = np.dot(A, x)

    vacuum_wavelengths = airtovac(wavelengths)
    return vacuum_wavelengths

load_old_DRS_S2D()

load the data from the old HARPS-N pipeline. This will be mainly used for the comparison with the HARPS-TERRA pipeline

Source code in src/ASTRA/Instruments/HARPSN.py 
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
def load_old_DRS_S2D(self):
    """
    load the data from the old HARPS-N pipeline. This will be mainly used for the comparison with the
    HARPS-TERRA pipeline
    """
    if not self._internal_configs["use_old_pipeline"]:
        raise custom_exceptions.InvalidConfiguration(
            "Can't load data from old pipeline without the config"
        )

    if self._internal_configs["use_old_pipeline"]:
        self.is_BERV_corrected = False

    with fits.open(self.file_path) as hdulist:
        # Compute the wavelength solution + BERV correction
        wave_from_file = self.build_HARPS_wavelengths(hdulist[0].header)

        sci_data = hdulist[0].data  # spetra from all orders

        # photon noise + estimate of max value for the rest
        # from ETC calculator the readout noise should be the largest contribution
        # assuming that it is of ~7e- (equal to manual) it should have a maximum contribution
        # of 200
        flux_errors = np.sqrt(250 + np.abs(sci_data, dtype=float))

        # Validate for overflows and missing data
        quality_data = np.zeros(sci_data.shape)
        quality_data[np.where(np.isnan(sci_data))] = NAN_DATA.code
        quality_data[np.where(sci_data > 300000)] = SATURATION.code
        quality_data[np.where(sci_data < -3 * flux_errors)] = ERROR_THRESHOLD.code

    self.spectra = sci_data.astype(np.float64)
    self.wavelengths = wave_from_file
    self.qual_data = quality_data
    self.uncertainties = flux_errors.astype(np.float64)

    self.build_mask(bypass_QualCheck=False)

close_arrays()

Reset the BERV correction flag if we are using the old pipeline version!

Source code in src/ASTRA/Instruments/HARPSN.py 
327
328
329
330
331
332
333
334
335
336
def close_arrays(self):
    """
    Reset the BERV correction flag if we are using the old pipeline version!
    Returns
    -------

    """
    super().close_arrays()
    if self._internal_configs["use_old_pipeline"]:
        self.is_BERV_corrected = False

find_files(file_name: UI_PATH)

Find the CCF and the S2D files, which should be stored inside the same folder

Source code in src/ASTRA/Instruments/HARPSN.py 
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
def find_files(self, file_name: UI_PATH):
    """
    Find the CCF and the S2D files, which should be stored inside the same folder
    """
    logger.debug("Searching for the ccf and e2ds files")

    if not isinstance(file_name, Path):
        file_name = Path(file_name)

    search_status = MISSING_DATA("Missing the ccf file")
    ccf_path = None

    if file_name.is_dir():
        logger.debug("Received a folder, searching inside for necessary files")
        # search for e2ds file
        folder_name = file_name

        e2ds_files = file_name.glob("**/*e2ds_A.fits")
        ccf_files = file_name.glob("**/*ccf_*_A.fits")

        for name, elems in [("e2ds_A", e2ds_files), ("ccf", ccf_files)]:
            if len(elems) > 1:
                msg = f"HARPS data only received folder name and it has more than 1 {name} file in it"
                raise custom_exceptions.InvalidConfiguration(msg)

            if len(elems) < 1:
                msg = f"HARPS data only received folder name and it has no {name} file in it"
                raise custom_exceptions.InvalidConfiguration(msg)

        e2ds_path = e2ds_files[0]
        ccf_path = ccf_files[0]
        search_status = SUCCESS("Found all input files")
    else:
        logger.debug(
            "Received path of E2DS file; searching for CCF with matching name"
        )
        folder_name = file_name.parent
        e2ds_path = file_name
        file_start, *_ = file_name.stem.split("_")

        found_CCF = False
        ccf_files = folder_name.glob("*ccf*_A.fits")

        for file in ccf_files:
            if file_start in file.name:
                ccf_path = file
                found_CCF = True

        if found_CCF:
            logger.info("Found CCF file: {}".format(ccf_path))
            search_status = SUCCESS("Found CCF file")
        else:
            logger.critical(
                "Was not able to find CCF file. Marking frame as invalid"
            )
            ccf_path = ""

    return e2ds_path, ccf_path, search_status

MAROONX

Interface to MAROON-X data.

MAROONX

Bases: Frame

Interface to handle MAROONX data.

Source code in src/ASTRA/Instruments/MAROONX.py 
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
class MAROONX(Frame):
    """Interface to handle MAROONX data."""

    _default_params = Frame._default_params

    # Adding one more day than the official time, so that we ensure to include any observations
    # from the last day

    sub_instruments = {
        # "MAROON1": datetime.datetime.strptime("09-15-2020", r"%m-%d-%Y"),
        # "MAROON2": datetime.datetime.strptime("12-02-2020", r"%m-%d-%Y"),
        # "MAROON3": datetime.datetime.strptime("03-04-2021", r"%m-%d-%Y"),
        # "MAROON4": datetime.datetime.strptime("04-30-2021", r"%m-%d-%Y"),
        # "MAROON5": datetime.datetime.strptime("06-04-2021", r"%m-%d-%Y"),
        # "MAROON6": datetime.datetime.strptime("08-23-2021", r"%m-%d-%Y"),
        # "MAROON7": datetime.datetime.strptime("11-23-2021", r"%m-%d-%Y"),
        # "MAROON8": datetime.datetime.strptime("04-27-2022", r"%m-%d-%Y"),
        # "MAROON9": datetime.datetime.strptime("06-03-2022", r"%m-%d-%Y"),
        # "MAROON10": datetime.datetime.strptime("08-15-2022", r"%m-%d-%Y"),
        # "MAROON11": datetime.datetime.strptime("07-11-2023", r"%m-%d-%Y"),
        # "MAROON12": datetime.datetime.strptime("10-28-2023", r"%m-%d-%Y"),
        # "MAROON13": datetime.datetime.strptime("11-29-2023", r"%m-%d-%Y"),
        # "MAROON14": datetime.datetime.strptime("01-03-2024", r"%m-%d-%Y"),
        # "MAROON15": datetime.datetime.strptime("04-24-2024", r"%m-%d-%Y"),
        # "MAROON16": datetime.datetime.strptime("06-13-2024", r"%m-%d-%Y"),
        # "MAROON17": datetime.datetime.strptime("07-18-2024", r"%m-%d-%Y"),
        # "MAROON18": datetime.datetime.strptime("08-20-2024", r"%m-%d-%Y"),
        # "MAROON19": datetime.datetime.strptime("10-15-2024", r"%m-%d-%Y"),
        # "MAROON20": datetime.datetime.strptime("01-09-2025", r"%m-%d-%Y"),
        # "MAROON21": datetime.datetime.strptime("02-04-2025", r"%m-%d-%Y"),
        "MAROON1": datetime.datetime.max,
    }
    _name = "MAROONX"

    order_intervals: dict[DETECTOR_DEFINITION, slice] = {
        DETECTOR_DEFINITION.WHITE_LIGHT: list(range(62)),
        DETECTOR_DEFINITION.RED_DET: list(range(34, 62)),
        DETECTOR_DEFINITION.BLUE_DET: list(range(0, 34)),
    }
    KW_map = {}

    def __init__(
        self,
        file_path: Path,
        user_configs: Optional[Dict[str, Any]] = None,
        reject_subInstruments=None,
        frameID=None,
        quiet_user_params: bool = True,
    ):
        """Construct MAROON-X object.

        Parameters
        ----------
        file_path
            Path to the S2D (or S1D) file.
        user_configs
            Dictionary whose keys are the configurable options of ESPRESSO (check above)
        reject_subInstruments
            Iterable of subInstruments to fully reject
        frameID
            ID for this observation. Only used for organization purposes by :class:`~SBART.data_objects.DataClass`

        """
        self._blaze_corrected = True

        super().__init__(
            inst_name=self._name,
            array_size={"S2D": [62, 4036]},
            file_path=file_path,
            frameID=frameID,
            KW_map=self.KW_map,
            available_indicators=("FWHM", "BIS SPAN"),
            user_configs=user_configs,
            reject_subInstruments=reject_subInstruments,
            need_external_data_load=False,
            quiet_user_params=quiet_user_params,
        )
        coverage = [500, 920]
        self.instrument_properties["wavelength_coverage"] = coverage
        self.instrument_properties["is_drift_corrected"] = False

        self.instrument_properties["resolution"] = 86_000

        # lat/lon from: https://geohack.toolforge.org/geohack.php?params=19_49_25_N_155_28_9_W
        lat, lon = 19.820667, -155.468056
        self.instrument_properties["EarthLocation"] = EarthLocation.from_geodetic(
            lat=lat, lon=lon, height=4214
        )

        # from https://www.mide.com/air-pressure-at-altitude-calculator
        # and convert from Pa to mbar
        self.instrument_properties["site_pressure"] = 599.4049

        self.is_BERV_corrected = False

    def get_spectral_type(self) -> str:
        """Get the spectral type from the filename."""
        if not self.file_path.name.endswith("hd5"):
            raise custom_exceptions.InternalError(
                "MAROON-X interface only recognizes hd5 files"
            )
        return "S2D"

    def load_instrument_specific_KWs(self, header) -> None:
        """Override parent class, does nothing."""
        ...

    def load_header_info(self) -> None:
        """Load information from the header."""
        store = pd.HDFStore(self.file_path, "r+")
        header_blue = store["header_blue"]
        header_red = store["header_red"]

        orders_blue = store["spec_blue"].index.levels[1]
        orders_red = store["spec_red"].index.levels[1]
        store.close()

        for order_set, header_det in [
            (orders_blue, header_blue),
            (orders_red, header_red),
        ]:
            for order in order_set:
                self.observation_info["orderwise_SNRs"].append(
                    float(header_det[f"SNR_{order}"])
                )
        for name, kw in [
            ("ISO-DATE", "MAROONX TELESCOPE TIME"),
            ("OBJECT", "MAROONX TELESCOPE TARGETNAME"),
        ]:
            self.observation_info[name] = header_blue[kw]

        for name, kw in [
            ("airmass", "MAROONX TELESCOPE AIRMASS"),
            ("relative_humidity", "MAROONX TELESCOPE HUMIDITY"),
            ("ambient_temperature", "MAROONX WEATHER TEMPERATURE"),  # TODO: check units
            ("BERV", "BERV_FLUXWEIGHTED_FRD"),
            ("JD", "JD_UTC_FLUXWEIGHTED_FRD"),
            ("EXPTIME", "EXPTIME"),
        ]:
            self.observation_info[name] = float(header_blue[kw])

        self.observation_info["BERV"] = self.observation_info["BERV"] * meter_second
        self.observation_info["drift"] = 0 * meter_second
        self.observation_info["drift_ERR"] = 0 * meter_second

        self.observation_info["BERV_FACTOR"] = (
            1
            + self.observation_info["BERV"].to(kilometer_second).value / SPEED_OF_LIGHT
        )
        # Convert ambient temperature to Kelvin
        self.observation_info["ambient_temperature"] = convert_temperature(
            self.observation_info["ambient_temperature"],
            old_scale="Celsius",
            new_scale="Kelvin",
        )

        # Note: we don't have DRS values for MAROON-X
        self.observation_info["DRS_RV"] = 0 * meter_second

        # Can't have zero uncertainties!!
        self.observation_info["DRS_RV_ERR"] = 0.1 * meter_second

        self.find_instrument_type()
        self.assess_bad_orders()

    def load_S2D_data(self) -> None:
        """Load the S2D data from the HD5 files."""
        if self.is_open:
            return
        super().load_S2D_data()

        try:
            store = pd.HDFStore(self.file_path.as_posix(), "r+")
        except Exception as e:
            print(e)
            breakpoint()

        spec_red = store["spec_red"]
        spec_blue = store["spec_blue"]
        store.close()

        red_pix = spec_red["wavelengths"][6].values[0].shape[0]
        blue_pix = spec_blue["wavelengths"][6].values[0].shape[0]
        blue_pad = red_pix - blue_pix

        blue_det_flux = np.vstack(spec_blue["optimal_extraction"][6])
        p_blue_det_flux = np.pad(
            blue_det_flux, ((0, 0), (0, blue_pad)), mode="constant"
        )
        blue_det_wave = np.vstack(spec_blue["wavelengths"][6])
        p_blue_det_wave = np.pad(
            blue_det_wave, ((0, 0), (0, blue_pad)), mode="constant"
        )
        blue_det_err = np.vstack(spec_blue["optimal_var"][6])
        p_blue_det_err = np.pad(blue_det_err, ((0, 0), (0, blue_pad)), mode="constant")

        red_det_flux = np.vstack(spec_red["optimal_extraction"][6])
        red_det_wave = np.vstack(spec_red["wavelengths"][6])
        red_det_err = np.vstack(spec_red["optimal_var"][6])

        self.wavelengths = np.vstack((p_blue_det_wave, red_det_wave))
        self.spectra = np.vstack((p_blue_det_flux, red_det_flux))
        self.uncertainties = np.vstack((p_blue_det_err, red_det_err))
        self.build_mask(bypass_QualCheck=True)

    def load_S1D_data(self) -> Mask:
        """Load S1D data, currently not implemented."""
        raise NotImplementedError

    def build_mask(self, bypass_QualCheck: bool = False) -> None:
        """Construct the pixel-wise mask."""
        # We evaluate the bad orders all at once
        super().build_mask(bypass_QualCheck, assess_bad_orders=False)

        bpmap0 = np.zeros((62, 4036), dtype=np.uint64)
        # Remove the first blue order
        bpmap0[0, :] = 1

        if self._internal_configs["SIGMA_CLIP_FLUX_VALUES"] > 0:
            sigma = self._internal_configs["SIGMA_CLIP_FLUX_VALUES"]
            for order_number in range(self.N_orders):
                cont = median_filter(self.spectra[order_number], size=500)
                inds = np.where(
                    self.spectra[order_number]
                    >= cont + sigma * self.uncertainties[order_number]
                )
                bpmap0[order_number, inds] |= 1
        self.spectral_mask.add_indexes_to_mask(np.where(bpmap0 != 0), QUAL_DATA)

        # remove extremely negative points!
        self.spectral_mask.add_indexes_to_mask(
            np.where(self.spectra < -3 * self.uncertainties), MISSING_DATA
        )

        self.spectral_mask.add_indexes_to_mask(
            np.where(self.uncertainties > 1e5), MISSING_DATA
        )
        self.assess_bad_orders()

    def close_arrays(self) -> None:
        """Close arrays in memory."""
        super().close_arrays()
        self.is_BERV_corrected = False

sub_instruments = {'MAROON1': datetime.datetime.max} class-attribute instance-attribute

order_intervals: dict[DETECTOR_DEFINITION, slice] = {DETECTOR_DEFINITION.WHITE_LIGHT: list(range(62)), DETECTOR_DEFINITION.RED_DET: list(range(34, 62)), DETECTOR_DEFINITION.BLUE_DET: list(range(0, 34))} class-attribute instance-attribute

KW_map = {} class-attribute instance-attribute

is_BERV_corrected = False instance-attribute

get_spectral_type() -> str

Get the spectral type from the filename.

Source code in src/ASTRA/Instruments/MAROONX.py 
120
121
122
123
124
125
126
def get_spectral_type(self) -> str:
    """Get the spectral type from the filename."""
    if not self.file_path.name.endswith("hd5"):
        raise custom_exceptions.InternalError(
            "MAROON-X interface only recognizes hd5 files"
        )
    return "S2D"

load_instrument_specific_KWs(header) -> None

Override parent class, does nothing.

Source code in src/ASTRA/Instruments/MAROONX.py 
128
129
130
def load_instrument_specific_KWs(self, header) -> None:
    """Override parent class, does nothing."""
    ...

load_header_info() -> None

Load information from the header.

Source code in src/ASTRA/Instruments/MAROONX.py 
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
def load_header_info(self) -> None:
    """Load information from the header."""
    store = pd.HDFStore(self.file_path, "r+")
    header_blue = store["header_blue"]
    header_red = store["header_red"]

    orders_blue = store["spec_blue"].index.levels[1]
    orders_red = store["spec_red"].index.levels[1]
    store.close()

    for order_set, header_det in [
        (orders_blue, header_blue),
        (orders_red, header_red),
    ]:
        for order in order_set:
            self.observation_info["orderwise_SNRs"].append(
                float(header_det[f"SNR_{order}"])
            )
    for name, kw in [
        ("ISO-DATE", "MAROONX TELESCOPE TIME"),
        ("OBJECT", "MAROONX TELESCOPE TARGETNAME"),
    ]:
        self.observation_info[name] = header_blue[kw]

    for name, kw in [
        ("airmass", "MAROONX TELESCOPE AIRMASS"),
        ("relative_humidity", "MAROONX TELESCOPE HUMIDITY"),
        ("ambient_temperature", "MAROONX WEATHER TEMPERATURE"),  # TODO: check units
        ("BERV", "BERV_FLUXWEIGHTED_FRD"),
        ("JD", "JD_UTC_FLUXWEIGHTED_FRD"),
        ("EXPTIME", "EXPTIME"),
    ]:
        self.observation_info[name] = float(header_blue[kw])

    self.observation_info["BERV"] = self.observation_info["BERV"] * meter_second
    self.observation_info["drift"] = 0 * meter_second
    self.observation_info["drift_ERR"] = 0 * meter_second

    self.observation_info["BERV_FACTOR"] = (
        1
        + self.observation_info["BERV"].to(kilometer_second).value / SPEED_OF_LIGHT
    )
    # Convert ambient temperature to Kelvin
    self.observation_info["ambient_temperature"] = convert_temperature(
        self.observation_info["ambient_temperature"],
        old_scale="Celsius",
        new_scale="Kelvin",
    )

    # Note: we don't have DRS values for MAROON-X
    self.observation_info["DRS_RV"] = 0 * meter_second

    # Can't have zero uncertainties!!
    self.observation_info["DRS_RV_ERR"] = 0.1 * meter_second

    self.find_instrument_type()
    self.assess_bad_orders()

load_S2D_data() -> None

Load the S2D data from the HD5 files.

Source code in src/ASTRA/Instruments/MAROONX.py 
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
def load_S2D_data(self) -> None:
    """Load the S2D data from the HD5 files."""
    if self.is_open:
        return
    super().load_S2D_data()

    try:
        store = pd.HDFStore(self.file_path.as_posix(), "r+")
    except Exception as e:
        print(e)
        breakpoint()

    spec_red = store["spec_red"]
    spec_blue = store["spec_blue"]
    store.close()

    red_pix = spec_red["wavelengths"][6].values[0].shape[0]
    blue_pix = spec_blue["wavelengths"][6].values[0].shape[0]
    blue_pad = red_pix - blue_pix

    blue_det_flux = np.vstack(spec_blue["optimal_extraction"][6])
    p_blue_det_flux = np.pad(
        blue_det_flux, ((0, 0), (0, blue_pad)), mode="constant"
    )
    blue_det_wave = np.vstack(spec_blue["wavelengths"][6])
    p_blue_det_wave = np.pad(
        blue_det_wave, ((0, 0), (0, blue_pad)), mode="constant"
    )
    blue_det_err = np.vstack(spec_blue["optimal_var"][6])
    p_blue_det_err = np.pad(blue_det_err, ((0, 0), (0, blue_pad)), mode="constant")

    red_det_flux = np.vstack(spec_red["optimal_extraction"][6])
    red_det_wave = np.vstack(spec_red["wavelengths"][6])
    red_det_err = np.vstack(spec_red["optimal_var"][6])

    self.wavelengths = np.vstack((p_blue_det_wave, red_det_wave))
    self.spectra = np.vstack((p_blue_det_flux, red_det_flux))
    self.uncertainties = np.vstack((p_blue_det_err, red_det_err))
    self.build_mask(bypass_QualCheck=True)

load_S1D_data() -> Mask

Load S1D data, currently not implemented.

Source code in src/ASTRA/Instruments/MAROONX.py 
230
231
232
def load_S1D_data(self) -> Mask:
    """Load S1D data, currently not implemented."""
    raise NotImplementedError

build_mask(bypass_QualCheck: bool = False) -> None

Construct the pixel-wise mask.

Source code in src/ASTRA/Instruments/MAROONX.py 
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
def build_mask(self, bypass_QualCheck: bool = False) -> None:
    """Construct the pixel-wise mask."""
    # We evaluate the bad orders all at once
    super().build_mask(bypass_QualCheck, assess_bad_orders=False)

    bpmap0 = np.zeros((62, 4036), dtype=np.uint64)
    # Remove the first blue order
    bpmap0[0, :] = 1

    if self._internal_configs["SIGMA_CLIP_FLUX_VALUES"] > 0:
        sigma = self._internal_configs["SIGMA_CLIP_FLUX_VALUES"]
        for order_number in range(self.N_orders):
            cont = median_filter(self.spectra[order_number], size=500)
            inds = np.where(
                self.spectra[order_number]
                >= cont + sigma * self.uncertainties[order_number]
            )
            bpmap0[order_number, inds] |= 1
    self.spectral_mask.add_indexes_to_mask(np.where(bpmap0 != 0), QUAL_DATA)

    # remove extremely negative points!
    self.spectral_mask.add_indexes_to_mask(
        np.where(self.spectra < -3 * self.uncertainties), MISSING_DATA
    )

    self.spectral_mask.add_indexes_to_mask(
        np.where(self.uncertainties > 1e5), MISSING_DATA
    )
    self.assess_bad_orders()

close_arrays() -> None

Close arrays in memory.

Source code in src/ASTRA/Instruments/MAROONX.py 
264
265
266
267
def close_arrays(self) -> None:
    """Close arrays in memory."""
    super().close_arrays()
    self.is_BERV_corrected = False

POET

POET

Bases: ESPRESSO

Interface to handle POET-ESPRESSO observations (S2D and S1D).

Source code in src/ASTRA/Instruments/POET.py 
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
class POET(ESPRESSO):
    """Interface to handle POET-ESPRESSO observations (S2D and S1D)."""

    _default_params = ESPRESSO._default_params + DefaultValues(
        DEFER_BERV_UPDATE=UserParam(
            False,
            constraint=BooleanValue,
            description="If True, only updated berv and rv when opening the data",
        ),
        SOLYARIS_PATH=UserParam(
            None,
            description="Path to the solyaris file to use for BERV/HERV computation",
        ),
        AUX_FILE_FOLDER=UserParam(
            None,
            description="Path to the folder in which the auxiliary files are stored. If None, the default folder will be used.",
        ),
    )

    _name = "PoET"

    def __init__(
        self,
        file_path,
        user_configs: Optional[Dict[str, Any]] = None,
        reject_subInstruments: Optional[Iterable[str]] = None,
        frameID: Optional[int] = None,
        quiet_user_params: bool = True,
    ):
        """ESPRESSO interface.

        Parameters
        ----------
        file_path
            Path to the S2D (or S1D) file.
        user_configs
            Dictionary whose keys are the configurable options of ESPRESSO (check above)
        reject_subInstruments
            Iterable of subInstruments to fully reject
        frameID
            ID for this observation. Only used for organization purposes by :class:`~SBART.data_objects.DataClass`

        """
        # Wavelength coverage

        coverage = (350, 900)
        self.is_resolved_observation = "_B" in file_path.as_posix()
        self._updated_drs_rv = False

        super().__init__(
            file_path=file_path,
            frameID=frameID,
            user_configs=user_configs,
            reject_subInstruments=reject_subInstruments,
            quiet_user_params=quiet_user_params,
            is_solar=True,
        )

    def _load_ESO_DRS_KWs(self, header):
        super()._load_ESO_DRS_KWs(header)

        self.is_poet_data = header.get("ESO INS POET MODE", False)
        if not self.is_poet_data:
            raise custom_exceptions.FrameError(f"{self.name} is not PoET file")

        logger.info("Detected PoET frame")
        self.UT_number = 5
        self.observation_info["POET_APERTURE"] = header.get("ESO POET APER", "UNKNOWN")

        if self.observation_info["POET_APERTURE"] == "UNKNOWN":
            logger.critical("Could not load aperture from %s", self.file_path)

        if not self._internal_configs["DEFER_BERV_UPDATE"]:
            self._update_drs_rv()

    def _update_drs_rv(self) -> None:
        if self._updated_drs_rv:
            return
        logger.debug(f"Updating BERV solution for {self.fname}")

        self._updated_drs_rv = True
        self.wrong_berv = self.observation_info["BERV"]
        try:
            true_bary, true_herv = compute_berv_and_herv_POET(self._header)
        except Exception as e:
            logger.critical(f"Failed to compute BERV and HERV due to {e}")
            true_bary, true_herv = 0 * kilometer_second, 0 * kilometer_second
            self.add_to_status(FATAL_KW("Failed BERV and HERV correction"))

        true_berv = true_bary + true_herv

        new_rv = self.observation_info["DRS_RV"] - self.wrong_berv + true_berv
        self.observation_info["DRS_RV"] = new_rv
        self.observation_info["BERV"] = true_berv
        self.observation_info["HERV"] = true_herv

        self.observation_info["BERV_FACTOR"] = None  # Avoid re-write later on
        self.observation_info["MAX_BERV"] = 2.5 * kilometer_second

    def load_telemetry_info(self, header):
        # Find the UT number and load the airmass
        super().load_telemetry_info(header)

        if "ESO TEL5 AMBI FWHM START" not in header:
            # Backup since some files had the wrong keyword
            logger.critical("seeing KW not in header, falling back to old KW")
            logger.warning(header["*AMBI FWHM"])
            self.observation_info["seeing"] = float(
                header[f"HIERARCH ESO TEL{self.UT_number} AMBI FWHM"]
            )

    def load_S1D_data(self) -> None:
        # Store the open status before calling parent classes
        is_open = self.is_open
        super().load_S1D_data()
        if not is_open:
            self.apply_poet_data_corrections()

    def load_S2D_data(self) -> None:
        # Store the open status before calling parent classes
        is_open = self.is_open
        super().load_S2D_data()
        if not is_open:
            self.apply_poet_data_corrections()

    def apply_poet_data_corrections(self) -> None:
        # Correct spectra
        if self._internal_configs["DEFER_BERV_UPDATE"]:
            self._update_drs_rv()

        self.remove_BERV_correction(self.wrong_berv)
        self.apply_BERV_correction(self.observation_info["BERV"])

    def trigger_data_storage(self, *args, **kwargs):
        super().trigger_data_storage(*args, **kwargs)

is_resolved_observation = '_B' in file_path.as_posix() instance-attribute

load_telemetry_info(header)

Source code in src/ASTRA/Instruments/POET.py 
125
126
127
128
129
130
131
132
133
134
135
def load_telemetry_info(self, header):
    # Find the UT number and load the airmass
    super().load_telemetry_info(header)

    if "ESO TEL5 AMBI FWHM START" not in header:
        # Backup since some files had the wrong keyword
        logger.critical("seeing KW not in header, falling back to old KW")
        logger.warning(header["*AMBI FWHM"])
        self.observation_info["seeing"] = float(
            header[f"HIERARCH ESO TEL{self.UT_number} AMBI FWHM"]
        )

load_S1D_data() -> None

Source code in src/ASTRA/Instruments/POET.py 
137
138
139
140
141
142
def load_S1D_data(self) -> None:
    # Store the open status before calling parent classes
    is_open = self.is_open
    super().load_S1D_data()
    if not is_open:
        self.apply_poet_data_corrections()

load_S2D_data() -> None

Source code in src/ASTRA/Instruments/POET.py 
144
145
146
147
148
149
def load_S2D_data(self) -> None:
    # Store the open status before calling parent classes
    is_open = self.is_open
    super().load_S2D_data()
    if not is_open:
        self.apply_poet_data_corrections()

apply_poet_data_corrections() -> None

Source code in src/ASTRA/Instruments/POET.py 
151
152
153
154
155
156
157
def apply_poet_data_corrections(self) -> None:
    # Correct spectra
    if self._internal_configs["DEFER_BERV_UPDATE"]:
        self._update_drs_rv()

    self.remove_BERV_correction(self.wrong_berv)
    self.apply_BERV_correction(self.observation_info["BERV"])

trigger_data_storage(*args, **kwargs)

Source code in src/ASTRA/Instruments/POET.py 
159
160
def trigger_data_storage(self, *args, **kwargs):
    super().trigger_data_storage(*args, **kwargs)

load_POET_extra_information(self: DataClass) -> None

Source code in src/ASTRA/Instruments/POET.py 
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
def load_POET_extra_information(self: "DataClass") -> None:

    poet_ob_names, frameIDs = self.collect_KW_observations(
        KW="OBS NAME",
        subInstruments=self.get_subInstruments_with_valid_frames(),
        include_invalid=True,
        return_frameIDs=True,
    )

    poet_iso_dates, frameIDs = self.collect_KW_observations(
        KW="ISO-DATE",
        subInstruments=self.get_subInstruments_with_valid_frames(),
        include_invalid=True,
        return_frameIDs=True,
    )

    # Mapping to go from day to fID
    map_of_days: dict[str, list[int]] = defaultdict(list)
    for name, fID in zip(poet_iso_dates, frameIDs):
        map_of_days[name.split("T")[0]].append(fID)

    sol_path: Path = self.observations[0]._internal_configs["SOLYARIS_PATH"]
    aux_path: Path = self.observations[0]._internal_configs["AUX_FILE_FOLDER"]

    sol_data = None
    # Apply extra corrections to Sun as a Star data (differential atmospheric extinction)
    if sol_path is not None:
        sol_path = Path(sol_path)
        if not sol_path.exists():
            raise custom_exceptions.InternalError(
                f"SOLYARIS path ({sol_path}) does not exist"
            )
        for day, frameIDs in map_of_days.items():
            if sol_path.is_dir():
                local_sol_path = sol_path.glob(f"**/*{day}*csv")
                try:
                    sol_data = pd.read_csv(next(local_sol_path))
                except StopIteration:
                    logger.critical(
                        f"Couldn't find solyaris file for {day} in {sol_path}"
                    )
                    for fID in frameIDs:
                        frame = self.get_frame_by_ID(fID)
                        frame.add_to_status(FATAL_KW("No SOLYARIS info found"))
                    continue
            else:
                sol_data = pd.read_csv(sol_path)

            for frameID in frameIDs:
                frame = self.get_frame_by_ID(frameID)
                if frame.is_resolved_observation:
                    # No correction for fiber B data
                    continue

                raw_fname = frame.fname.replace("r.", "").split("_")[0] + ".fits"

                try:
                    row = sol_data[sol_data["file"] == raw_fname].iloc[0]
                except IndexError:
                    frame.add_to_status(f"No SOLYARIS info for {raw_fname}")
                    continue

                rv_corr = row["corr_extinction_tot"] * kilometer_second
                frame.import_KW_from_outside(
                    "DIFF_ATMOS_RV_CORR", rv_corr, optional=True
                )
                frame.observation_info["DRS_RV"] = (
                    frame.observation_info["DRS_RV"] + rv_corr
                )
                # Thresholds provided by Khaled Al Moulla
                qual_value = (
                    (row["quality_68"] > 0.95)
                    & (row["c0_hat_68"] < 0.1)
                    & (row["jit_68"] < 0.005)
                )
                if not qual_value:
                    frame.add_to_status(
                        FATAL_KW(f"Failed SOLYARIS QC checks for {frame.fname}")
                    )
                frame.finalized_external_data_load()

    else:
        logger.warning(
            "SOLYARIS path was not provided, skipping atmospheric extinction and QC checks"
        )
        for fID in frameIDs:
            frame = self.get_frame_by_ID(fID)
            if not frame.is_resolved_observation:
                frame.finalized_external_data_load()

    if aux_path is not None:
        aux_path = Path(aux_path)
        for day, frameIDs in map_of_days.items():
            all_aux_files = aux_path.glob(f"**/POET.INT.{day}*fits")

            # OB name to frameID
            map_of_files: dict[str, list[int]] = defaultdict(list)
            for OBname, fID in zip(poet_ob_names, frameIDs):
                map_of_files[OBname].append(fID)

            poe_to_header_map = {}
            for aux_file in all_aux_files:
                aux_header = fits.getheader(aux_file, ext=1)
                aux_ob_name = aux_header["ESO OBS NAME"]
                if aux_ob_name not in map_of_files:
                    logger.debug(
                        f"Could not find any frame for {aux_ob_name} in {aux_file}"
                    )
                    continue

                poe_to_header_map[aux_ob_name] = aux_header

        for poet_ob_name, frameIDs in map_of_files.items():

            head = poe_to_header_map.get(poet_ob_name)

            for frame in frameIDs:
                frame = self.get_frame_by_ID(frame)

                for int_kw, head_kwy in [
                    ("PROG ID", "ESO ST ID"),
                    ("TX", "ESO TEL5 TARG TX"),
                    ("TY", "ESO TEL5 TARG TY"),
                    ("ALT_OFF", "ESO MCS MOTION AL DYNAMICOFFSET"),
                    ("AZ_OFF", "ESO MCS MOTION AZ DYNAMICOFFSET"),
                    ("COMMENT", "COMMENT"),
                ]:
                    val = head[head_kwy]

                    if int_kw in ["TX", "TY"]:
                        val = float(val) * u.arcsec
                    if int_kw in ["ALT_OFF", "AZ_OFF"]:
                        val = float(val)

                    if int_kw in ["COMMENT"]:
                        val = " ".join(str(v) for v in val)

                    frame.import_KW_from_outside(KW=int_kw, value=val, optional=True)

                frame.finalized_external_data_load()
    else:
        logger.warning("AUXILLIARY path was not provided, skipping loading of metadata")
        for fID in frameIDs:
            frame = self.get_frame_by_ID(fID)
            if frame.is_resolved_observation:
                frame.finalized_external_data_load()

SimulatedSpirou

Interface to MAROON-X data.

SimulatedSpirou

Bases: Frame

Interface to handle simulated data.

Source code in src/ASTRA/Instruments/SimulatedSpirou.py 
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
class SimulatedSpirou(Frame):
    """Interface to handle simulated data."""

    _default_params = Frame._default_params + DefaultValues(
        IS_BERV_CORRECTED=UserParam(
            default_value=False,
            constraint=BooleanValue,
            description="If False, flags the file as not BERV corrected",
        )
    )

    # Adding one more day than the official time, so that we ensure to include any observations
    # from the last day

    sub_instruments = {
        "SimulatedSpirou": datetime.datetime.max,
    }
    _name = "SimulatedSpirou"

    order_intervals: dict[DETECTOR_DEFINITION, slice] = {
        DETECTOR_DEFINITION.WHITE_LIGHT: list(range(48)),
    }
    KW_map = {}

    def __init__(
        self,
        file_path: Path,
        user_configs: Optional[Dict[str, Any]] = None,
        reject_subInstruments=None,
        frameID=None,
        quiet_user_params: bool = True,
    ):
        """Construct MAROON-X object.

        Parameters
        ----------
        file_path
            Path to the S2D (or S1D) file.
        user_configs
            Dictionary whose keys are the configurable options of ESPRESSO (check above)
        reject_subInstruments
            Iterable of subInstruments to fully reject
        frameID
            ID for this observation. Only used for organization purposes by :class:`~SBART.data_objects.DataClass`

        """
        self._blaze_corrected = True

        head = fits.getheader(file_path, ext=1)
        arr_size = [head["NAXIS2"], head["NAXIS1"]]
        super().__init__(
            inst_name=self._name,
            array_size={"S2D": arr_size},
            file_path=file_path,
            frameID=frameID,
            KW_map=self.KW_map,
            available_indicators=[],
            user_configs=user_configs,
            reject_subInstruments=reject_subInstruments,
            need_external_data_load=False,
            quiet_user_params=quiet_user_params,
        )
        coverage = [955, 2438]
        self.instrument_properties["wavelength_coverage"] = coverage
        self.instrument_properties["is_drift_corrected"] = True

        self.instrument_properties["resolution"] = 100_000

        # lat/lon from: https://geohack.toolforge.org/geohack.php?params=19_49_25_N_155_28_9_W
        lat, lon = 19.820667, -155.468056
        self.instrument_properties["EarthLocation"] = EarthLocation.from_geodetic(lat=lat, lon=lon, height=4214)

        # from https://www.mide.com/air-pressure-at-altitude-calculator
        # and convert from Pa to mbar
        self.instrument_properties["site_pressure"] = 599.4049

        self.is_BERV_corrected = self._internal_configs["IS_BERV_CORRECTED"]

    def get_spectral_type(self):
        return "S2D"

    def load_instrument_specific_KWs(self, header):
        # Load BERV info + previous RV
        self.observation_info["MAX_BERV"] = 30 * kilometer_second
        self.observation_info["BERV"] = header["BERV"] * meter_second

        self.observation_info["DRS_RV"] = header["DRS_RV"] * meter_second
        self.observation_info["DRS_RV_ERR"] = header["DRS_RV_ERR"] * meter_second
        self.observation_info["BJD"] = header["BJD"]

        t = Time(header["BJD"], format="jd", scale="tdb")
        iso_string = t.iso
        self.observation_info["ISO-DATE"] = "T".join(iso_string.split(" "))

        for order in range(self.N_orders):
            self.observation_info["orderwise_SNRs"].append(header[f"SNR {order}"])

    def check_header_QC(self, header: fits.Header):
        """Header QC checks for CARMENES.

        1) Drift calibration was done (CARACAL DRIFT FP REF exists)
        2) Time between observation and calibration is smaller than "max_hours_to_calibration"

        Can add the following status:

        - KW_WARNING("Drift flag of KOBE is greater than 1")
            If th drift value is greater than one. This is actually set in the DataClass.load_CARMENES_extra_information()

        Args:
            header (fits.Header): _description_
        """
        ...

    def load_S2D_data(self):
        if self.is_open:
            return
        super().load_S2D_data()

        with fits.open(self.file_path) as hdulist:
            s2d_data = hdulist["SPEC"].data * 100000  # spetra from all olders
            err_data = hdulist["SIG"].data * 100000
            wavelengths = hdulist["WAVE"].data  # vacuum wavelengths; no BERV correction

        self.wavelengths = wavelengths
        self.spectra = s2d_data
        self.uncertainties = err_data
        self.build_mask(bypass_QualCheck=True)
        return 1

    def load_S1D_data(self) -> Mask:
        raise NotImplementedError

    def build_mask(self, bypass_QualCheck: bool = False):
        # We evaluate the bad orders all at once
        super().build_mask(bypass_QualCheck, assess_bad_orders=False)

        # remove extremely negative points!
        self.spectral_mask.add_indexes_to_mask(np.where(self.spectra < -3 * self.uncertainties), MISSING_DATA)
        self.spectral_mask.add_indexes_to_mask(np.where(self.uncertainties == 0), MISSING_DATA)

        self.assess_bad_orders()

    def close_arrays(self):
        super().close_arrays()

        self.is_BERV_corrected = self._internal_configs["IS_BERV_CORRECTED"]

sub_instruments = {'SimulatedSpirou': datetime.datetime.max} class-attribute instance-attribute

order_intervals: dict[DETECTOR_DEFINITION, slice] = {DETECTOR_DEFINITION.WHITE_LIGHT: list(range(48))} class-attribute instance-attribute

KW_map = {} class-attribute instance-attribute

is_BERV_corrected = self._internal_configs['IS_BERV_CORRECTED'] instance-attribute

get_spectral_type()

Source code in src/ASTRA/Instruments/SimulatedSpirou.py 
102
103
def get_spectral_type(self):
    return "S2D"

load_instrument_specific_KWs(header)

Source code in src/ASTRA/Instruments/SimulatedSpirou.py 
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
def load_instrument_specific_KWs(self, header):
    # Load BERV info + previous RV
    self.observation_info["MAX_BERV"] = 30 * kilometer_second
    self.observation_info["BERV"] = header["BERV"] * meter_second

    self.observation_info["DRS_RV"] = header["DRS_RV"] * meter_second
    self.observation_info["DRS_RV_ERR"] = header["DRS_RV_ERR"] * meter_second
    self.observation_info["BJD"] = header["BJD"]

    t = Time(header["BJD"], format="jd", scale="tdb")
    iso_string = t.iso
    self.observation_info["ISO-DATE"] = "T".join(iso_string.split(" "))

    for order in range(self.N_orders):
        self.observation_info["orderwise_SNRs"].append(header[f"SNR {order}"])

check_header_QC(header: fits.Header)

Header QC checks for CARMENES.

1) Drift calibration was done (CARACAL DRIFT FP REF exists) 2) Time between observation and calibration is smaller than "max_hours_to_calibration"

Can add the following status:

  • KW_WARNING("Drift flag of KOBE is greater than 1") If th drift value is greater than one. This is actually set in the DataClass.load_CARMENES_extra_information()

Args: header (fits.Header): description

Source code in src/ASTRA/Instruments/SimulatedSpirou.py 
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
def check_header_QC(self, header: fits.Header):
    """Header QC checks for CARMENES.

    1) Drift calibration was done (CARACAL DRIFT FP REF exists)
    2) Time between observation and calibration is smaller than "max_hours_to_calibration"

    Can add the following status:

    - KW_WARNING("Drift flag of KOBE is greater than 1")
        If th drift value is greater than one. This is actually set in the DataClass.load_CARMENES_extra_information()

    Args:
        header (fits.Header): _description_
    """
    ...

load_S2D_data()

Source code in src/ASTRA/Instruments/SimulatedSpirou.py 
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
def load_S2D_data(self):
    if self.is_open:
        return
    super().load_S2D_data()

    with fits.open(self.file_path) as hdulist:
        s2d_data = hdulist["SPEC"].data * 100000  # spetra from all olders
        err_data = hdulist["SIG"].data * 100000
        wavelengths = hdulist["WAVE"].data  # vacuum wavelengths; no BERV correction

    self.wavelengths = wavelengths
    self.spectra = s2d_data
    self.uncertainties = err_data
    self.build_mask(bypass_QualCheck=True)
    return 1

load_S1D_data() -> Mask

Source code in src/ASTRA/Instruments/SimulatedSpirou.py 
153
154
def load_S1D_data(self) -> Mask:
    raise NotImplementedError

build_mask(bypass_QualCheck: bool = False)

Source code in src/ASTRA/Instruments/SimulatedSpirou.py 
156
157
158
159
160
161
162
163
164
def build_mask(self, bypass_QualCheck: bool = False):
    # We evaluate the bad orders all at once
    super().build_mask(bypass_QualCheck, assess_bad_orders=False)

    # remove extremely negative points!
    self.spectral_mask.add_indexes_to_mask(np.where(self.spectra < -3 * self.uncertainties), MISSING_DATA)
    self.spectral_mask.add_indexes_to_mask(np.where(self.uncertainties == 0), MISSING_DATA)

    self.assess_bad_orders()

close_arrays()

Source code in src/ASTRA/Instruments/SimulatedSpirou.py 
166
167
168
169
def close_arrays(self):
    super().close_arrays()

    self.is_BERV_corrected = self._internal_configs["IS_BERV_CORRECTED"]