"""
Module for Keck/ESI specific methods.
.. include:: ../include/links.rst
"""
import os
from IPython import embed
import numpy as np
from astropy.time import Time
import datetime
from pypeit import msgs
from pypeit import telescopes
from pypeit import io
from pypeit.core import framematch
from pypeit.core import parse
from pypeit.spectrographs import spectrograph
from pypeit.images import detector_container
[docs]
class KeckESISpectrograph(spectrograph.Spectrograph):
"""
Child to handle Keck/ESI specific code
"""
ndet = 1
name = 'keck_esi'
camera = 'ESI'
header_name = 'ESI'
#url = 'https://www.lco.cl/?epkb_post_type_1=mage'
ech_fixed_format = True
telescope = telescopes.KeckTelescopePar()
pypeline = 'Echelle'
supported = True
#comment = 'See :doc:`mage`'
[docs]
def get_detector_par(self, det, hdu=None):
"""
Return metadata for the selected detector.
Args:
det (:obj:`int`):
1-indexed detector number.
hdu (`astropy.io.fits.HDUList`_, optional):
The open fits file with the raw image of interest. If not
provided, frame-dependent parameters are set to a default.
Returns:
:class:`~pypeit.images.detector_container.DetectorContainer`:
Object with the detector metadata.
"""
# Binning
# TODO: Could this be detector dependent??
binning = '1,1' if hdu is None else self.get_meta_value(self.get_headarr(hdu), 'binning')
# Detector 1
detector_dict = dict(
binning = binning,
det = 1,
dataext = 0,
specaxis = 0,
specflip = False,
spatflip = False,
# plate scale in arcsec/pixel
platescale = 0.1542,
# electrons/pixel/hour.
darkcurr = 2.10, # e/pixel/hour... Note : Could be updated
saturation = 65535.,
# CCD is linear to better than 0.5 per cent up to digital saturation (65,536 DN including bias) in the Fast readout mode.
nonlinear = 0.99,
mincounts = -1e10,
numamplifiers = 2,
gain = np.atleast_1d([1.3, 1.3]),
ronoise = np.atleast_1d([2.5, 2.5]),
)
return detector_container.DetectorContainer(**detector_dict)
[docs]
@classmethod
def default_pypeit_par(cls):
"""
Return the default parameters to use for this instrument.
Returns:
:class:`~pypeit.par.pypeitpar.PypeItPar`: Parameters required by
all of ``PypeIt`` methods.
"""
par = super().default_pypeit_par()
# Bias
#par['calibrations']['biasframe']['useframe'] = 'overscan'
# Wavelengths
# 1D wavelength solution
# This is for 1x1
par['calibrations']['wavelengths']['rms_thresh_frac_fwhm'] = 0.103
par['calibrations']['wavelengths']['fwhm'] = 2.9
#
par['calibrations']['wavelengths']['sigdetect'] = 5.0
par['calibrations']['wavelengths']['lamps'] = ['CuI', 'ArI', 'NeI', 'HgI', 'XeI', 'ArII']
par['calibrations']['wavelengths']['method'] = 'reidentify'
par['calibrations']['wavelengths']['cc_thresh'] = 0.50
par['calibrations']['wavelengths']['cc_local_thresh'] = 0.50
# Reidentification parameters
par['calibrations']['wavelengths']['reid_arxiv'] = 'keck_esi_ECH.fits'
#par['calibrations']['wavelengths']['ech_fix_format'] = True
# Echelle parameters
par['calibrations']['wavelengths']['echelle'] = True
par['calibrations']['wavelengths']['ech_nspec_coeff'] = 4
par['calibrations']['wavelengths']['ech_norder_coeff'] = 4
par['calibrations']['wavelengths']['ech_sigrej'] = 3.0
par['scienceframe']['process']['sigclip'] = 20.0
par['scienceframe']['process']['satpix'] = 'nothing'
# Set slits and tilts parameters
par['calibrations']['tilts']['tracethresh'] = 10. #[10]*self.norders
par['calibrations']['slitedges']['fit_order'] = 5
par['calibrations']['slitedges']['max_shift_adj'] = 3.
par['calibrations']['slitedges']['left_right_pca'] = True
par['calibrations']['slitedges']['edge_thresh'] = 5.0
par['calibrations']['slitedges']['det_min_spec_length'] = 0.2
par['calibrations']['slitedges']['fit_min_spec_length'] = 0.4
par['calibrations']['slitedges']['pca_sigrej'] = 1.5
par['calibrations']['slitedges']['pca_order'] = 3
par['calibrations']['slitedges']['add_missed_orders'] = True
# Find object parameters
par['reduce']['findobj']['find_trim_edge'] = [4,4] # Slit is too short to trim 5,5 especially with 2x binning
par['reduce']['findobj']['maxnumber_sci'] = 2 # Slit is narrow so allow one object per order
par['reduce']['findobj']['maxnumber_std'] = 1 # Slit is narrow so allow one object per order
par['reduce']['extraction']['model_full_slit'] = True # local sky subtraction operates on entire slit
# Scattered light
par['calibrations']['pixelflatframe']['process']['subtract_scattlight'] = True
par['calibrations']['illumflatframe']['process']['subtract_scattlight'] = True
par['scienceframe']['process']['subtract_scattlight'] = True
# Always flux calibrate, starting with default parameters
# Do not correct for flexure
par['flexure']['spec_method'] = 'skip'
# Set the default exposure time ranges for the frame typing
par['calibrations']['standardframe']['exprng'] = [None, 60]
par['calibrations']['arcframe']['exprng'] = [300, None] # Allow for CuAr which can be quite long
par['calibrations']['darkframe']['exprng'] = [1, None]
par['scienceframe']['exprng'] = [60, None]
return par
[docs]
def configuration_keys(self):
"""
Return the metadata keys that define a unique instrument
configuration.
This list is used by :class:`~pypeit.metadata.PypeItMetaData` to
identify the unique configurations among the list of frames read
for a given reduction.
Returns:
:obj:`list`: List of keywords of data pulled from file headers
and used to constuct the :class:`~pypeit.metadata.PypeItMetaData`
object.
"""
return []
[docs]
def check_frame_type(self, ftype, fitstbl, exprng=None):
"""
Check for frames of the provided type.
Args:
ftype (:obj:`str`):
Type of frame to check. Must be a valid frame type; see
frame-type :ref:`frame_type_defs`.
fitstbl (`astropy.table.Table`_):
The table with the metadata for one or more frames to check.
exprng (:obj:`list`, optional):
Range in the allowed exposure time for a frame of type
``ftype``. See
:func:`pypeit.core.framematch.check_frame_exptime`.
Returns:
`numpy.ndarray`_: Boolean array with the flags selecting the
exposures in ``fitstbl`` that are ``ftype`` type frames.
"""
good_exp = framematch.check_frame_exptime(fitstbl['exptime'], exprng)
if ftype in ['pinhole', 'dark']:
# No pinhole or pinhole or dark frames
return np.zeros(len(fitstbl), dtype=bool)
if ftype in ['bias']:
return fitstbl['idname'] == 'Bias'
if ftype in ['pixelflat', 'trace', 'illumflat']:
ans = np.zeros(len(fitstbl), dtype=bool)
for kk, idnm in enumerate(fitstbl['idname']):
if idnm in ['DmFlat', 'IntFlat', 'SkyFlat']:
ans[kk] = True
return ans
if ftype in ['arc', 'tilt']:
return fitstbl['idname'] == 'Line'
if ftype == 'science':
return good_exp & (fitstbl['idname'] == 'Object')
if ftype == 'standard':
return good_exp & (fitstbl['idname'] == 'Object')
msgs.warn('Cannot determine if frames are of type {0}.'.format(ftype))
return np.zeros(len(fitstbl), dtype=bool)
[docs]
def bpm(self, filename, det, shape=None, msbias=None):
"""
Generate a default bad-pixel mask.
Even though they are both optional, either the precise shape for
the image (``shape``) or an example file that can be read to get
the shape (``filename`` using :func:`get_image_shape`) *must* be
provided.
Args:
filename (:obj:`str` or None):
An example file to use to get the image shape.
det (:obj:`int`):
1-indexed detector number to use when getting the image
shape from the example file.
shape (tuple, optional):
Processed image shape
Required if filename is None
Ignored if filename is not None
msbias (`numpy.ndarray`_, optional):
Master bias frame used to identify bad pixels
Returns:
`numpy.ndarray`_: An integer array with a masked value set
to 1 and an unmasked value set to 0. All values are set to
0.
"""
# Call the base-class method to generate the empty bpm
bpm_img = super().bpm(filename, det, shape=shape, msbias=msbias)
# Get the binning
msgs.info("Custom bad pixel mask for ESI")
hdu = io.fits_open(filename)
binspatial, binspec = parse.parse_binning(hdu[0].header['BINNING'])
hdu.close()
# Binning Independent Masking
bpm_img[:,0:(2//binspatial)] = 1
# Mask out the 'hotspot' on the upper left coner
bpm_img[(3842//binspec):(3944//binspec), (19//binspatial):(161//binspatial)] = 1
# Mask out the 'bad columns' on the upper left coner
bpm_img[(2642//binspec):, (418//binspatial):(442//binspatial)] = 1
# Return
return bpm_img
[docs]
def scattered_light_archive(self, binning, dispname):
"""Archival model parameters for the scattered light. These are based on best fits to currently available data.
Parameters
----------
binning : :obj:`str`, optional
Comma-separated binning along the spectral and spatial directions; e.g., ``2,1``
dispname : :obj:`str`, optional
Name of the disperser
Returns
-------
x0 : `numpy.ndarray`_
A 1D array containing the best-fitting model parameters
bounds : :obj:`tuple`
A tuple of two elements, containing two `numpy.ndarray`_ of the same length as x0. These
two arrays contain the lower (first element of the tuple) and upper (second element of the tuple)
bounds to consider on the scattered light model parameters.
"""
# Grab the binning for convenience
specbin, spatbin = parse.parse_binning(binning)
# Get some starting parameters (these were determined by fitting spectra,
# and should be close to the final fitted values to reduce computational time)
# Note :: These values need to be originally based on data that uses 1x1 binning,
# and are now scaled here according to the binning of the current data to be analysed.
# These parameters give a cost of 8.0517e+08 with the science frame used as scattlight (1x1 binning, pad=5)
x0 = np.array([272.33958742493064 / specbin, 115.501464689107 / spatbin, # Gaussian kernel widths
272.3418000034377 / specbin, 168.0591427733949 / spatbin, # Lorentzian kernel widths
-141.2552517318941 / specbin, 79.25936221285629 / spatbin, # pixel offsets
1.0877734248786808, 1.0562808322123667, # Zoom factor (spec, spat)
5.876311151022701, # constant flux offset
0.0444248025888341, # kernel angle
0.6090358292193677, # Relative kernel scale (>1 means the kernel is more Gaussian, >0 but <1 makes the profile more lorentzian)
0.135392229831296, -0.16167521454188258, # Polynomial terms (coefficients of "spat" and "spat*spec")
0.06148093592863097, 0.10305719952486242]) # Polynomial terms (coefficients of spec**index)
# Now set the bounds of the fitted parameters
bounds = ([# Lower bounds:
1, 1, # Gaussian kernel widths
1, 1, # Lorentzian kernel widths
-200 / specbin, -200 / spatbin, # pixel offsets
0, 0, # Zoom factor (spec, spat)
-1000, -(10 / 180) * np.pi, 0.0, # constant flux offset, kernel angle, relative kernel scale
-10, -10, -10, -10], # Polynomial terms
# Upper bounds
[600 / specbin, 600 / spatbin, # Gaussian kernel widths
600 / specbin, 600 / spatbin, # Lorentzian kernel widths
200 / specbin, 200 / spatbin, # pixel offsets
2, 2, # Zoom factor (spec, spat)
1000.0, +(10 / 180) * np.pi, 1000.0, # constant flux offset, kernel angle, relative kernel scale
10, 10, 10, 10]) # Polynomial terms
# Return the best-fitting archival parameters and the bounds
return x0, bounds
@property
def norders(self):
"""
Number of orders for this spectrograph. Should only defined for
echelle spectrographs, and it is undefined for the base class.
"""
return 10 # 15-6
@property
def order_spat_pos(self):
"""
Return the expected spatial position of each echelle order.
"""
return np.array([0.115, 0.245, 0.362, 0.465, 0.558, 0.642, 0.719, 0.791, 0.861, 0.933])
@property
def order_spat_width(self):
"""
Return the expected spatial width of each slit trace for each order,
relative to the spatial size of the detector.
"""
return np.array([0.0879, 0.0818, 0.0779, 0.0747, 0.0720, 0.0696, 0.0676, 0.0658, 0.0640,
0.0617])
@property
def orders(self):
"""
Return the order number for each echelle order.
"""
return np.arange(15, 5, -1, dtype=int)
@property
def spec_min_max(self):
"""
Return the minimum and maximum spectral pixel expected for the
spectral range of each order.
"""
spec_max = np.full(self.norders, np.inf)
spec_min = np.full(self.norders, -np.inf)
return np.vstack((spec_min, spec_max))
[docs]
def order_platescale(self, order_vec, binning=None):
"""
Return the platescale for each echelle order.
This routine is only defined for echelle spectrographs, and it is
undefined in the base class.
Args:
order_vec (`numpy.ndarray`_):
The vector providing the order numbers.
binning (:obj:`str`, optional):
The string defining the spectral and spatial binning.
Returns:
`numpy.ndarray`_: An array with the platescale for each order
provided by ``order``.
"""
norders = len(order_vec)
binspatial, binspec = parse.parse_binning(binning)
# Plate scales
unbinned_pscale = [0.120, #15
0.127, 0.134, 0.137, 0.144, 0.149, 0.153,
0.158, 0.163,
0.168, # 6
]
return np.array(unbinned_pscale)*binspatial
[docs]
def get_rawimage(self, raw_file, det, spectrim=None):
""" Read the image
"""
# Check for file; allow for extra .gz, etc. suffix
if not os.path.isfile(raw_file):
msgs.error(f'{raw_file} not found!')
hdu = io.fits_open(raw_file)
head0 = hdu[0].header
# Number of AMPS
namp = head0['NUMAMPS']
# Get post, pre-pix values
prepix = head0['PREPIX']
postpix = head0['POSTPIX']
preline = head0['PRELINE']
postline = head0['POSTLINE']
# Grab the data
full_image = hdu[0].data.astype(float)
rawdatasec_img = np.zeros_like(full_image, dtype=int)
oscansec_img = np.zeros_like(full_image, dtype=int)
#
nspat = int(head0['WINDOW'].split(',')[3]) // namp
for amp in range(namp):
col0 = prepix*2 + nspat*amp
# Data
rawdatasec_img[:, col0:col0+nspat] = amp+1
# Overscan
o0 = prepix*2 + nspat*namp + postpix*amp
oscansec_img[:, o0:o0+postpix] = amp+1
#embed(header='435 of keck_esi.py')
return self.get_detector_par(1, hdu=hdu), \
full_image, hdu, head0['ELAPTIME'], rawdatasec_img, oscansec_img