As a first step for the wavelength and drift solution, the etalon lines need to be fitted to determine their centroid. This is done in the 1D extracted spectrum. The etalon lines are unresolved, hence, their 1D line profile is basically the LSF of the spectrograph. We model this using a box convolved with two Gaussians on either side, i.e., the sum of two error functions. In order to remove degeneracies between the width of the box (representing the width of the entrance slit) and the sigma/FWHM of the Gaussians (representing the field and wavelength dependent aberrations), a composite model for each fiber/order is built from all the individual etalon lines in a fiber/order plus a smooth background. Only the line intensities and positions are fitted individually for each etalon line. Box widths, Gaussian sigmas, and the background are constrained to a low order polynomial across a given fiber/order. This is motivated by the fact that these quantities only vary slowly and steadily along an order.
The fitting procedure is multithreaded per fiber/order and a batch version provides functionality to process more than one file at a time, see discussion below.
1D (hdf) master flat field.
Processed 1D (hdf) spectrum containing etalon lines in one or more fibers.
Etalon parameters added to the input hdf file.
Log file with ending _etalon_peakfit.log
Depending on the setting of the --show plot parameter, a number of plots are displayed (see below).
Command line:
From the base directory (maroonx_reduce) call: PYTHONPATH=${PWD} python reduce/extract_etalon_positions.py
with the following parameters:
-f FLAT, --flat FLAT Absolute path and filename of reduced masterflat
-i INPUT_FILE, --input_file INPUT_FILE Absolute path and filename of reduced etalon spectrum
-g GUESS_FILE, --guess_file GUESS_FILE Absolute path and filename of reduced and extracted etalon spectrum that is used as initial guess
--output_file OUTPUT_FILE Store results in different file than the input file
-fs FIBERS, --fibers FIBERS Fibers to fit. Can either be specified as list of ints (e.g. '1 2 5') or as a range (e.g. '2-3').
Defaults to all fibers present in the provided master flat.
-o ORDERS, --orders ORDERS Orders to fit. Can either be specified as list of ints (e.g. '83 84') or as a range (e.g. '83-86').
Defaults to all orders present in the provided master flat.
-w WRITE_OUTPUT, --write_output WRITE_OUTPUT Save fit results in file (either in input, or if provided, output file). Default: True
-l LOGLEVEL, --loglevel LOGLEVEL Set log level, use names provided by logging. Default: 'info'
-p PLOT_PATH, --plot_path PLOT_PATH Create plots and store at path
--degrees DEGREES Comma separated list of polynomial degrees: sigma, width. Default: 4,1
--use_sigma_lr Use different polynomials for left and right sigma parameters. Default: False
--show_plots If true, plots will show up for debugging...
--multithreading If true, multithreading will be used and plotting disabled
Example: PYTHONPATH=${PWD} python reduce/extract_etalon_positions.py
-i '/data2/MaroonX_spectra_reduced/20201119/20201119T180308Z_EDDDE_r_0010.hdf'
-f '/data2/MaroonX_spectra_reduced/Maroonx_masterframes/202011xx/flats/
20201119T18_masterflat_backgroundsubtracted_FFFFF_r_0000.hdf'
--use_sigma_lr --degrees 4,2 -fs 1 --show_plots -o 83
Command line example call to fit etalon spectra for fiber 1 (sky) in order 83 with different sigmas for left and right gaussian and a slope in the box-width parameter of 2nd degree. Plots are shown before and after fitting converges (see below)
From the base directory (maroonx_reduce) call: PYTHONPATH=${PWD} python reduce/batch_extract_etalon_positions.py
with the following parameters:
-t TYPE, --type TYPE Exposure type, e.g. 'SOOOE' or 'DEEEE'.
-dd DATA_DIRECTORY, --data_directory DATA_DIRECTORY Directory for input hdf files. Default = '/data2/MaroonX_spectra_reduced'
-d DATE, --date DATE UTC date of observation, e.g. '20200901', '202009', '2020112[3-4]' or '20201121T18'
-e EXPTIME, --exptime EXPTIME Exposure time in sec, e.g. '300'. Use '????' for any/all
-fr RED_FLAT, --red_flat RED_FLAT Reduced master flat for red camera (hdf).
Default: /data2/MaroonX_spectra_reduced/Maroonx_masterframes/202011xx/flats/
20201119T18_masterflat_backgroundsubtracted_FFFFF_r_0000.hdf
-fb BLUE_FLAT, --blue_flat BLUE_FLAT Reduced master flat for blue camera (hdf).
Default: /data2/MaroonX_spectra_reduced/Maroonx_masterframes/202011xx/flats/
20201119T18_masterflat_backgroundsubtracted_FFFFF_r_0000.hdf
-c {b,r}, --camera_arm {b,r} Camera arm, e.g. 'b' or 'r'.
Example: PYTHONPATH=${PWD} python reduce/batch_extract_etalon_positions.py -c 'r' -d '20201121' -t 'DEEEE' -e '????'
Example call for the multi-file / multi-threaded routine to fit all fibers with source designation 'E' in its filename for all orders in the red arm. All files matching 'DEEEE' type and 20201121 as a date are processed.
Blue and red frames are processed separately. Takes fibers, orders, and flux (for de-blazing) from flats.
A peak finding algorithm (find_peaks() in etalon_fit.py) automatically detects the etalon lines in the 1D spectrum. Lots of measures have been taken here to remove outliers (bad pixels, cosmics, etc.) and make sure maxima and minima are in the right order, etc. However, this step is the most likely to produce problems.
Guess parameters for all parameters are used for an initial fit which is already pretty good in most cases. The fit is then iterated in a nested fashion over the individual line parameters (position, amplitude) and the ones bound globally by polynomials (background, width, sigmas). The number of iterations is hard-coded (parameter iterations=3 of iterative_fit() in etalon_fit.py.
The half-width of the slit typically varies from 1.3 to 2.0 pixels from the left to right of each order, predominantly linear with a slight curvature. This is expected from the anamorphic magnification of the echelle resulting in the design sampling of ~3.5pix per FWHM at the center of each order. The gaussian sigmas show a more complex variation along each order and across the spectrum as they capture the aberrations.
Best results are achieved with different gaussians for left and right (--use_sigma_lr), 4th order for sigma and 2nd order for width (--degrees 4,2).
Originally, only the line centers were fitted individually, all other parameters were bound to polynomials. The Jacobian's are written accordingly. Later we decided to also allow the line strengths to be fitted individually to account for fringing patterns. The fit seems to work fine, but I don't see this change reflected in the Jacobians. Might be worth checking.
The number of iterations could be made a command-line parameter. However, fits with less than 3 iterations were not as good, more iterations than 3 are very expensive. So maybe no need to open this up?
Programming 'style' seems different than in other parts of the data reduction code. For example the way data are retrieved and stored with helper functions in etalon_data.py seems unique.
Instead of having a slim 'batch' script like for the flux extraction that only takes care of the multi-threading (or in this case multi-file input), but calls on the main routine, I re-used the full extract_etalon_positions.py file and overwrote its __main__, effectively giving it a different set of command line arguments and adding the multi-file part to it. Not very elegant but the fastest approach. The actual computations are handled by etalon_data.py and etalon_fit.py.
QC: Procedure should raise a flag if it fails (e.g., fitting does not converge) and produce a plot of the fiber/order that didn't fit correctly.
Most common errors are (a) No etalon spectrum in designated fiber (i.e. wrong fiber specified or wrong source specifier in filename). (b) Automatic detection of peaks did not succeed to deal with image artifacts (e.g. cosmics). Automatic handling of both cases could be improved.