Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

DM-33326: Fix round-tripping of spectractor spectra #17

Merged
merged 6 commits into from
Feb 17, 2022
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
2 changes: 2 additions & 0 deletions .gitignore
Original file line number Diff line number Diff line change
Expand Up @@ -26,3 +26,5 @@ spectractor/version.py
.sconsign.dblite
config.log
tests/.tests/
*__pycache__/

121 changes: 82 additions & 39 deletions spectractor/extractor/spectrum.py
Original file line number Diff line number Diff line change
Expand Up @@ -6,6 +6,7 @@
import matplotlib.pyplot as plt
import numpy as np
import os
import astropy

from spectractor import parameters
from spectractor.config import set_logger, load_config, apply_rebinning_to_parameters
Expand All @@ -20,6 +21,22 @@
from spectractor.simulation.throughput import TelescopeTransmission


fits_mappings = {'date_obs': 'DATE-OBS',
'expo': 'EXPTIME',
'airmass': 'AIRMASS',
'disperser_label': 'GRATING',
'units': 'UNIT2',
'rotation_angle': 'ROTANGLE',
'dec': 'DEC',
'hour_angle': 'HA',
'temperature': 'OUTTEMP',
'pressure': 'OUTPRESS',
'humidity': 'OUTHUM',
'lambda_ref': 'LBDA_REF',
'parallactic_angle': 'PARANGLE',
}


class Spectrum:
""" Class used to store information and methods relative to spectra and their extraction.

Expand Down Expand Up @@ -123,7 +140,9 @@ class Spectrum:
Size of the spectrogram along the y axis.
"""

def __init__(self, file_name="", image=None, order=1, target=None, config="", fast_load=False):
def __init__(self, file_name="", image=None, order=1, target=None, config="", fast_load=False,
spectrogram_file_name_override=None,
psf_file_name_override=None,):
""" Class used to store information and methods relative to spectra and their extraction.

Parameters
Expand Down Expand Up @@ -203,9 +222,17 @@ def __init__(self, file_name="", image=None, order=1, target=None, config="", fa
self.lambdas_order2 = None
self.data_order2 = None
self.err_order2 = None
self.dec = None
self.hour_angle = None
self.temperature = None
self.pressure = None
self.humidity = None
self.parallactic_angle = None
self.filename = file_name
if file_name != "":
self.load_spectrum(file_name)
self.load_spectrum(file_name,
spectrogram_file_name_override=spectrogram_file_name_override,
psf_file_name_override=psf_file_name_override)
if image is not None:
self.header = image.header
self.date_obs = image.date_obs
Expand Down Expand Up @@ -484,6 +511,18 @@ def save_spectrum(self, output_file_name, overwrite=False):
'third column the corresponding errors.'
hdu1 = fits.PrimaryHDU()
hdu1.header = self.header

for attribute, header_key in fits_mappings.items():
try:
value = getattr(self, attribute)
except AttributeError:
print(f"Failed to get {attribute}")
continue
if isinstance(value, astropy.coordinates.angles.Angle):
value = value.degree
hdu1.header[header_key] = value
print(f"Set header key {header_key} to {value} from attr {attribute}")

hdu1.header["EXTNAME"] = "SPECTRUM"
hdu2 = fits.ImageHDU()
hdu2.header["EXTNAME"] = "SPEC_COV"
Expand Down Expand Up @@ -553,14 +592,21 @@ def save_spectrogram(self, output_file_name, overwrite=False):
hdu.writeto(output_file_name, overwrite=overwrite)
self.my_logger.info('\n\tSpectrogram saved in %s' % output_file_name)

def load_spectrum(self, input_file_name):
def load_spectrum(self, input_file_name, spectrogram_file_name_override=None,
psf_file_name_override=None):
"""Load the spectrum from a fits file (data, error and wavelengths).

Parameters
----------
input_file_name: str
Path to the input fits file

spectrogram_file_name_override : str
Manually specify a path to the spectrogram file.

psf_file_name_override : str
Manually specify a path to the psf file.

Examples
--------
>>> s = Spectrum()
Expand All @@ -575,48 +621,33 @@ def load_spectrum(self, input_file_name):
self.data = raw_data[1]
if len(raw_data) > 2:
self.err = raw_data[2]
if self.header['DATE-OBS'] != "":
self.date_obs = self.header['DATE-OBS']
if self.header['EXPTIME'] != "":
self.expo = self.header['EXPTIME']
if self.header['AIRMASS'] != "":
self.airmass = self.header['AIRMASS']
if self.header['GRATING'] != "":
self.disperser_label = self.header['GRATING']
if self.header['TARGET'] != "":
self.target = load_target(self.header['TARGET'], verbose=parameters.VERBOSE)

# set the simple items from the mappings. More complex items, i.e.
# those needing function calls, follow
for attribute, header_key in fits_mappings.items():
if (item := self.header.get(header_key)) is not None:
setattr(self, attribute, item)
print(f'set {attribute} to {item}')
else:
print(f'Failed to set spectrum attribute {attribute} using header {header_key}')

# set the more complex items by hand here
if target := self.header.get('TARGET'):
self.target = load_target(target, verbose=parameters.VERBOSE)
self.lines = self.target.lines
if self.header['UNIT2'] != "":
self.units = self.header['UNIT2']
if self.header['ROTANGLE'] != "":
self.rotation_angle = self.header['ROTANGLE']
if self.header['TARGETX'] != "" and self.header['TARGETY'] != "":
self.x0 = [self.header['TARGETX'], self.header['TARGETY']]
if self.header['D2CCD'] != "":
parameters.DISTANCE2CCD = float(self.header["D2CCD"])
if 'DEC' in self.header and self.header['DEC'] != "":
self.dec = self.header['DEC']
if 'RA' in self.header and self.header['HA'] != "":
self.hour_angle = self.header['HA']
if 'OUTTEMP' in self.header and self.header['OUTTEMP'] != "":
self.temperature = self.header['OUTTEMP']
if 'OUTPRESS' in self.header and self.header['OUTPRESS'] != "":
self.pressure = self.header['OUTPRESS']
if 'OUTHUM' in self.header and self.header['OUTHUM'] != "":
self.humidity = self.header['OUTHUM']
if self.header['LBDA_REF'] != "":
self.lambda_ref = self.header['LBDA_REF']
if 'PARANGLE' in self.header and self.header['PARANGLE'] != "":
self.parallactic_angle = self.header['PARANGLE']
if 'CCDREBIN' in self.header and self.header['CCDREBIN'] != "":
if parameters.CCD_REBIN != self.header['CCDREBIN']:
if (targetx := self.header.get('TARGETX')) and (targety := self.header.get('TARGETY')):
self.x0 = [targetx, targety] # should be a tuple not a list
if rebin := self.header.get('CCDREBIN'):
if parameters.CCD_REBIN != rebin:
raise ValueError("Different values of rebinning parameters between config file and header. Choose.")
parameters.CCD_REBIN = rebin
if dist := self.header.get('D2CCD'):
parameters.DISTANCE2CCD = float(dist)

self.my_logger.info('\n\tLoading disperser %s...' % self.disperser_label)
self.disperser = Hologram(self.disperser_label, D=parameters.DISTANCE2CCD,
data_dir=parameters.DISPERSER_DIR, verbose=parameters.VERBOSE)
self.my_logger.info('\n\tSpectrum loaded from %s' % input_file_name)
spectrogram_file_name = input_file_name.replace('spectrum', 'spectrogram')
if parameters.OBS_OBJECT_TYPE == "STAR":
self.adr_params = [self.dec, self.hour_angle, self.temperature,
self.pressure, self.humidity, self.airmass]
Expand All @@ -632,13 +663,25 @@ def load_spectrum(self, input_file_name):
self.target.image_y0 = float(hdu_list["ORDER0"].header["IM_Y0"])
else:
self.cov_matrix = np.diag(self.err ** 2)

# original, hard-coded spectrogram/table relative paths
spectrogram_file_name = input_file_name.replace('spectrum', 'spectrogram')
psf_file_name = input_file_name.replace('spectrum.fits', 'table.csv')

# for LSST-DM supplied filenames
if spectrogram_file_name_override and psf_file_name_override:
self.fast_load = False
spectrogram_file_name = spectrogram_file_name_override
psf_file_name = psf_file_name_override
self.my_logger.info(f'Applying spectrogram filename override {spectrogram_file_name}')
self.my_logger.info(f'Applying psf filename override {psf_file_name}')

if not self.fast_load:
self.my_logger.info(f'\n\tLoading spectrogram from {spectrogram_file_name}...')
if os.path.isfile(spectrogram_file_name):
self.load_spectrogram(spectrogram_file_name)
else:
raise FileNotFoundError(f"Spectrogram file {spectrogram_file_name} does not exist.")
psf_file_name = input_file_name.replace('spectrum.fits', 'table.csv')
self.my_logger.info(f'\n\tLoading PSF from {psf_file_name}...')
if os.path.isfile(psf_file_name):
self.load_chromatic_psf(psf_file_name)
Expand Down
4 changes: 3 additions & 1 deletion spectractor/extractor/targets.py
Original file line number Diff line number Diff line change
Expand Up @@ -275,7 +275,9 @@ def load_spectra(self):
else:
self.wavelengths.append(data.wave)
self.spectra.append(data.flux)
elif 'HD' in self.label or self.label in parameters.STAR_NAMES: # it is a star
# TODO DM-33731: the use of self.label in parameters.STAR_NAMES:
# below works for running but breaks a test so needs fixing for DM
elif 'HD' in self.label: # or self.label in parameters.STAR_NAMES: # it is a star
self.emission_spectrum = False
self.hydrogen_only = False
self.lines = Lines(ATMOSPHERIC_LINES + HYDROGEN_LINES + STELLAR_LINES,
Expand Down