pypeit.wavetilts module
Module for guiding Arc/Sky line tracing
- class pypeit.wavetilts.BuildWaveTilts(mstilt, slits, spectrograph, par, wavepar, det=1, qa_path=None, spat_flexure=None)[source]
Bases:
object
Class to guide arc/sky tracing
- Parameters
mstilt (
pypeit.images.buildimage.TiltImage
) – Tilt image. QA file naming inherits the calibration key (calib_key
) from this object.slits (
pypeit.slittrace.SlitTraceSet
) – Slit edgesspectrograph (
pypeit.spectrographs.spectrograph.Spectrograph
) – Spectrograph objectpar (
pypeit.par.pypeitpar.WaveTiltsPar
or None) – The parameters used to fuss with the tiltswavepar (
pypeit.par.pypeitpar.WaveSolutionPar
or None) – The parameters used for the wavelength solutiondet (int) – Detector index
qa_path (
str
, optional) – Directory for QA output.spat_flexure (float, optional) – If input, the slitmask and slit edges are shifted prior to tilt analysis.
- spectrograph
- steps
list
- mask
ndarray, bool True = Ignore this slit
- all_trcdict
list of dict All trace dict’s
- tilts
ndarray Tilts for a single slit/order
- all_ttilts
list of tuples Tuple of tilts ndarray’s
- final_tilts
ndarray Final tilts image
- gpm
Good pixel mask Eventually, we might attach this to self.mstilt although that would then require that we write it to disk with self.mstilt.image
- Type
- extract_arcs()[source]
Extract the arcs down each slit/order
Wrapper to arc.get_censpec()
- Returns
Extracted arcs in two numpy.ndarray objects
- Return type
- find_lines(arcspec, slit_cen, slit_idx, bpm=None, debug=False)[source]
Find the lines for tracing
Wrapper to tracewave.tilts_find_lines()
- Parameters
arcspec –
slit_cen –
slit_idx (int) – Slit index, zero-based
bpm (numpy.ndarray, optional) –
debug (bool, optional) –
- Returns
- 2 objectcs
numpy.ndarray or None: Spectral positions of lines to trace
numpy.ndarray or None: Spatial positions of lines to trace
- Return type
- fit_tilts(trc_tilt_dict, thismask, slit_cen, spat_order, spec_order, slit_idx, show_QA=False, doqa=True)[source]
Fit the tilts
all_fit_dict and all_trace_dict are filled in place
- Parameters
trc_tilt_dict (dict) – Contains information from tilt tracing
slit_cen (ndarray) – (nspec,) Central trace for this slit
spat_order (int) – Order of the 2d polynomial fit for the spatial direction
spec_order (int) – Order of the 2d polytnomial fit for the spectral direction
slit_idx (int) – zero-based, integer index for the slit in question
- Optional Args:
- show_QA: bool, default = False
show the QA instead of writing it out to the outfile
- doqa: bool, default = True
Construct the QA plot
- Returns
- coeff: ndarray (spat_order + 1, spec_order+1)
Array containing the coefficients for the 2d legendre polynomial fit
- Return type
- make_tbl_tilt_traces()[source]
Make a Table of the traced and fitted tilts from self.all_trace_dict to save into the WaveTilts object
- Returns
Table including the traced and fitted tilts for each slit. Columns are:
slit_ids: Slit IDs for which the tilts were traced and fit
goodpix_spat: Good spatial pixels of the traced tilts
goodpix_tilt: Good spectral pixels of the traced tilts
goodpix_lid: Line IDs for each goodpix. This is needed to associate goodpix to each line
badpix_spat: Masked spatial pixels of the traced tilts
badpix_tilt: Masked spectral pixels of the traced tilts
badpix_lid: Line IDs for each badpix. This is needed to associate badpix to each line
good2dfit_spat: Good spatial pixels of the 2D fit of the tilts
good2dfit_tilt: Good spectral pixels of the 2D fit of the tilts
good2dfit_lid: Line IDs for each good2dfit. This is needed to associate good2dfit to each line
bad2dfit_spat: Rejected spatial pixels of the 2D fit of the tilts
bad2dfit_tilt: Rejected spectral pixels of the 2D fit of the tilts
bad2dfit_lid: Line IDs for each bad2dfit. This is needed to associate bad2dfit to each line
- Return type
- model_arc_continuum(debug=False)[source]
Model the continuum of the arc image.
The method uses the arc spectra extracted using
extract_arcs
and fits a characteristic low-order continuum for each slit/order usingpypeit.util.robust_polyfit_djs()
and the parameters cont_function, cont_order, and cont_rej frompar
. The characteristic continuum is then rescaled to match the continuum at each spatial position in the slit/order.Note
The approach used here may be too simplistic (in the robustness of the continuum fit and then how the continuum is rescaled and projected for each slit/order). Tests should be performed to make sure that the approach is good enough to trace the centroid of the arc/sky lines without biasing the centroids due to a non-zero continuum level.
- Parameters
debug (
bool
, optional) – Run the method in debug mode.- Returns
Returns a 2D image with the same shape as
mstilt
with the model continuum.- Return type
numpy.ndarray
- run(doqa=True, debug=False, show=False)[source]
Main driver for tracing arc lines
Code flow:
Extract an arc spectrum down the center of each slit/order
- Loop on slits/orders
Trace and fit the arc lines (This is done twice, once with trace_crude as the tracing crutch, then again with a PCA model fit as the crutch).
Repeat trace.
2D Fit to the offset from slitcen
Save
- trace_tilts(arcimg, lines_spec, lines_spat, thismask, slit_cen, debug_pca=False, show_tracefits=False)[source]
Trace the tilts
- Parameters
arcimg (numpy.ndarray) – Arc image. Shape is (nspec, nspat).
lines_spec (numpy.ndarray) – Array containing the spectral pixel location of each line found for this slit. Shape is (nlines,).
lines_spat (numpy.ndarray) – Array containing the spatial pixel location of each line, which is the slitcen evaluate at the spectral position position of the line stored in lines_spec. Shape is (nlines,).
thismask (numpy.ndarray) – Image indicating which pixels lie on the slit in equation. True = on the slit. False = not on slit. Shape is (nspec, nspat) with dtype=bool.
slit_cen (
int
) – Integer index indicating the slit in question.
- Returns
Dictionary containing information on the traced tilts required to fit the filts.
- Return type
- class pypeit.wavetilts.WaveTilts(coeffs, nslit, spat_id, spat_order, spec_order, func2d, bpmtilts=None, spat_flexure=None, PYP_SPEC=None, slits_filename=None, tiltimg_filename=None, tilt_traces=None)[source]
Bases:
CalibFrame
Calibration frame containing the wavelength tilt calibration.
All of the items in the datamodel are required for instantiation, although they can be None (but shouldn’t be)
The datamodel attributes are:
Version: 1.2.0
Attribute
Type
Array Type
Description
PYP_SPEC
str
PypeIt spectrograph name
bpmtilts
Bad pixel mask for tilt solutions. Keys are taken from SlitTraceSetBitmask
coeffs
2D coefficents for the fit on the initial slits. One set per slit/order (3D array).
func2d
str
Function used for the 2D fit
nslit
int
Total number of slits. This can include masked slits
slits_filename
str
Path to SlitTraceSet file. This helps to find the Slits calibration file when running pypeit_chk_tilts()
spat_flexure
float
Flexure shift from the input TiltImage
spat_id
Slit spat_id
spat_order
Order for spatial fit (nslit)
spec_order
Order for spectral fit (nslit)
tilt_traces
Table with the positions of the traced and fitted tilts for all the slits. see
make_tbl_tilt_traces()
for more details.tiltimg_filename
str
Path to Tiltimg file. This helps to find Tiltimg file when running pypeit_chk_tilts()
When written to an output-file HDU, all numpy.ndarray elements are bundled into an astropy.io.fits.BinTableHDU, and the other elements are written as header keywords. Any datamodel elements that are None are not included in the output.
- _bundle()[source]
Bundle the data in preparation for writing to a fits file.
See
pypeit.datamodel.DataContainer._bundle()
. Data is always written to a ‘TILTS’ extension.
- calib_file_format = 'fits'
The extension and file format of the output file. Should be
'fits'
or'fits.gz'
(for gzipped output).
- calib_type = 'Tilts'
The type of the calibration frame, primarily used to set the name of the output file.
- datamodel = {'PYP_SPEC': {'descr': 'PypeIt spectrograph name', 'otype': <class 'str'>}, 'bpmtilts': {'atype': <class 'numpy.integer'>, 'descr': 'Bad pixel mask for tilt solutions. Keys are taken from SlitTraceSetBitmask', 'otype': <class 'numpy.ndarray'>}, 'coeffs': {'atype': <class 'numpy.floating'>, 'descr': '2D coefficents for the fit on the initial slits. One set per slit/order (3D array).', 'otype': <class 'numpy.ndarray'>}, 'func2d': {'descr': 'Function used for the 2D fit', 'otype': <class 'str'>}, 'nslit': {'descr': 'Total number of slits. This can include masked slits', 'otype': <class 'int'>}, 'slits_filename': {'descr': 'Path to SlitTraceSet file. This helps to find the Slits calibration file when running pypeit_chk_tilts()', 'otype': <class 'str'>}, 'spat_flexure': {'descr': 'Flexure shift from the input TiltImage', 'otype': <class 'float'>}, 'spat_id': {'atype': <class 'numpy.integer'>, 'descr': 'Slit spat_id ', 'otype': <class 'numpy.ndarray'>}, 'spat_order': {'atype': <class 'numpy.integer'>, 'descr': 'Order for spatial fit (nslit)', 'otype': <class 'numpy.ndarray'>}, 'spec_order': {'atype': <class 'numpy.integer'>, 'descr': 'Order for spectral fit (nslit)', 'otype': <class 'numpy.ndarray'>}, 'tilt_traces': {'descr': 'Table with the positions of the traced and fitted tilts for all the slits. see :func:`make_tbl_tilt_traces` for more details. ', 'otype': <class 'astropy.table.table.Table'>}, 'tiltimg_filename': {'descr': 'Path to Tiltimg file. This helps to find Tiltimg file when running pypeit_chk_tilts()', 'otype': <class 'str'>}}
Default datamodel for any
CalibFrame
. Derived classes should instantiate their datamodels by first inheriting from the base class. E.g.:class ArcFrame(CalibFrame): datamodel = {**CalibFrame.datamodel, ...}
- fit2tiltimg(slitmask, flexure=None)[source]
Generate a tilt image from the fit parameters
Mainly to allow for flexure
- Parameters
slitmask (numpy.ndarray) –
flexure (float, optional) – Spatial shift of the tilt image onto the desired frame (typically a science image)
- Returns
New tilt image
- Return type
- is_synced(slits)[source]
Confirm the slits in WaveTilts are aligned to that in SlitTraceSet
Barfs if not
- Parameters
slits (
pypeit.slittrace.SlitTraceSet
) –
- show(waveimg=None, wcs_match=True, in_ginga=True, show_traces=False)[source]
Show in ginga or mpl Tiltimg with the tilts traced and fitted overlaid
- Parameters
waveimg (numpy.ndarray, optional) – Image with the wavelength solution.
wcs_match (bool, optional) – If True, use this image as a reference image for the WCS and match all image in other channels to it.
in_ginga (bool, optional) – If True, show the image in ginga. Otherwise, use matplotlib.
show_traces (bool, optional) – If True, show the traces of the tilts on the image.
- version = '1.2.0'
Provides the string representation of the class version.
This is currently put to minimal use so far, but will used for I/O verification in the future.
Each derived class should provide a version to guard against data model changes during development.
- pypeit.wavetilts.show_tilts_mpl(tilt_img, tilt_traces, show_traces=False, left_edges=None, right_edges=None, slit_ids=None, cut=None)[source]
Show the TiltImage with the traced and 2D fitted tilts overlaid
- Parameters
tilt_img (numpy.ndarray) – TiltImage
tilt_traces (astropy.table.Table) – Table containing the traced and fitted tilts. See
make_tbl_tilt_traces()
for information on the table columns.show_traces (bool, optional) – Show the traced tilts
left_edges (numpy.ndarray, optional) – Left and right edges of the slits
right_edges (numpy.ndarray, optional) – Left and right edges of the slits
slit_ids (numpy.ndarray, optional) – Slit IDs
cut (tuple, optional) – Lower and upper levels cut for the image display