CR2RE Pipeline Reference Manual 1.6.2
Functions
Wavelength Calibration

Functions

int cr2res_wave_apply (cpl_table *tw_in, cpl_table *spectra_tab, const cpl_frame *catalog_frame, int reduce_order, int reduce_trace, cr2res_wavecal_type wavecal_type, int degree, int xdegree, double wl_start, double wl_end, double wl_err, double wl_shift, int log_flag, int fallback_input_wavecal_flag, int keep_higher_degrees_flag, int clean_spectrum, int display, double display_wmin, double display_wmax, int zp_order, int grat1_order, cpl_propertylist **qcs, cpl_table **lines_diagnostics, cpl_table **extracted_out, cpl_table **trace_wave_out)
 Apply the Wavelength Calibration.
 
cpl_polynomial * cr2res_wave_1d (cpl_bivector *spectrum, cpl_bivector *spectrum_err, cpl_polynomial *wavesol_init, const cpl_array *wave_error_init, int order, int trace_nb, cr2res_wavecal_type wavecal_type, const char *catalog, int degree, int clean_spectrum, int log_flag, int keep_higher_degrees_flag, int display, double display_wmin, double display_wmax, double *best_xcorr, double *lines_resol_fwhm, double *lines_pos, double *lines_resol, double *lines_intens, cpl_array **wavelength_error, cpl_table **lines_diagnostics)
 1D Wavelength Calibration
 
cpl_polynomial * cr2res_wave_2d (cpl_bivector **spectra, cpl_bivector **spectra_err, cpl_polynomial **wavesol_init, cpl_array **wavesol_init_err, int *orders, int *traces_nb, int ninputs, const char *catalog, cpl_size degree, cpl_size degree_y, double threshold, int n_iterations, int zp_order, int display, cpl_array **wavelength_error, cpl_table **lines_diagnostics)
 Compute the 2D wavelength polynomial based on a line spectrum and a reference catalog by finding lines and fitting.
 
cpl_polynomial * cr2res_wave_xcorr (cpl_bivector *spectrum, cpl_polynomial *wavesol_init, double wl_error, cpl_bivector *lines_list, int degree, int keep_higher_degrees_flag, int cleaning_filter_size, double slit_width, double fwhm, int display, double *best_xcorr, cpl_array **wavelength_error)
 Find solution by cross-correlating template spectrum.
 
cpl_polynomial * cr2res_wave_etalon (cpl_bivector *spectrum, const cpl_bivector *spectrum_err, cpl_polynomial *wavesol_init, int degree, cpl_array **wavelength_error)
 Find solution from etalon.
 
cpl_polynomial * cr2res_wave_estimate_compute (double wmin, double wmax)
 Compute the polynomial from boundaries.
 
cpl_array * cr2res_wave_get_estimate (const char *filename, int detector, int order)
 Compute the wavelength estimate.
 
hdrl_image * cr2res_wave_gen_wave_map (const cpl_table *trace_wave)
 Compute the wavelength map from the trace_wave table.
 
cpl_polynomial * cr2res_wave_polys_1d_to_2d (cpl_polynomial **poly_1ds, int *orders, int npolys, cpl_size xdegree)
 Create a 2D Wavelength Polynomial out of a several 1D ones.
 
cpl_polynomial * cr2res_wave_poly_2d_to_1d (cpl_polynomial *poly_2d, int order)
 Create a 1D Wavelength Polynomial out of a 2D one.
 

Detailed Description

Function Documentation

◆ cr2res_wave_1d()

cpl_polynomial * cr2res_wave_1d ( cpl_bivector *  spectrum,
cpl_bivector *  spectrum_err,
cpl_polynomial *  wavesol_init,
const cpl_array *  wave_error_init,
int  order,
int  trace_nb,
cr2res_wavecal_type  wavecal_type,
const char *  catalog,
int  degree,
int  clean_spectrum,
int  log_flag,
int  keep_higher_degrees_flag,
int  display,
double  display_wmin,
double  display_wmax,
double *  best_xcorr,
double *  lines_resol_fwhm,
double *  lines_pos,
double *  lines_resol,
double *  lines_intens,
cpl_array **  wavelength_error,
cpl_table **  lines_diagnostics 
)

1D Wavelength Calibration

Parameters
spectrumExtracted spectrum
spectrumExtracted spectrum error
wavesol_initInitial wavelength solution
wave_error_initInitial wavelength error (can be NULL)
orderOrder number
trace_nbTrace Number
catalogLine catalog or template spectrum or NULL
degreeThe polynomial degree of the solution
clean_spectrumRemove the lines that are not in the catalog
log_flagFlag to get the log() of the catalog intens.
keep_higher_degrees_flag
displayFlag to display results
display_wminMinimum Wavelength to display or -1.0
display_wmaxMaximum Wavelength to display or -1.0
best_xcorr[out] Best XCORR value when applicable
lines_resol_fwhm[out] Lines Fwhm min in pix
lines_pos[out] min fwhm line position in nm
lines_resol[out] Lines Fwhm min in nm
lines_intens[out] Lines intensities median
wavelength_error[out] array of wave_mean_error, wave_max_error
lines_diagnostics[out] table with lines diagnostics
Returns
Wavelength solution, i.e. polynomial that translates pixel values to wavelength.

Definition at line 689 of file cr2res_wave.c.

References cr2res_plot_wavecal_result(), cr2res_qc_wave_lamp_effic(), cr2res_qc_wave_resol_fwhm(), cr2res_wave_etalon(), and cr2res_wave_xcorr().

Referenced by cr2res_wave_apply().

◆ cr2res_wave_2d()

cpl_polynomial * cr2res_wave_2d ( cpl_bivector **  spectra,
cpl_bivector **  spectra_err,
cpl_polynomial **  wavesol_init,
cpl_array **  wavesol_init_err,
int *  orders,
int *  traces_nb,
int  ninputs,
const char *  catalog,
cpl_size  degree,
cpl_size  degree_y,
double  threshold,
int  n_iterations,
int  zp_order,
int  display,
cpl_array **  wavelength_error,
cpl_table **  lines_diagnostics 
)

Compute the 2D wavelength polynomial based on a line spectrum and a reference catalog by finding lines and fitting.

Parameters
spectraList of extracted spectra
spectra_errList of extracted spectra errors
wavesol_initList of Initial wavelength solutions
wavesol_init_errList of Initial wavelength error (can be NULL)
ordersList of orders of the various spectra
trace_nbList of traces IDs of the various spectra
ninputsNumber of entries in the previous parameters
catalogCatalog file
degreeThe polynomial degree in x
degree_yThe polynomial degree in y
displayFlag to display results
wavelength_error[out] array of wave_mean_error, wave_max_error
lines_diagnostics[out] table with lines diagnostics
Returns
Wavelength solution, i.e. polynomial that translates pixel values to wavelength.

NOTES

  • Some input spectra might be NULL - be robust against this.
  • Avoid hardcoded values (2)
  • Limit the line length to 80 characters
  • Keep the interface stable
  • Add Documentation in the function code
  • Remove dead code
  • All declarations at the beginning of a function
  • Make sure to deallocate everything

Definition at line 889 of file cr2res_wave.c.

References cr2res_dfs_create_lines_diagnostics_table(), cr2res_io_load_EMISSION_LINES(), and cr2res_wave_poly_2d_to_1d().

Referenced by cr2res_wave_apply().

◆ cr2res_wave_apply()

int cr2res_wave_apply ( cpl_table *  tw_in,
cpl_table *  spectra_tab,
const cpl_frame *  catalog_frame,
int  reduce_order,
int  reduce_trace,
cr2res_wavecal_type  wavecal_type,
int  degree,
int  xdegree,
double  wl_start,
double  wl_end,
double  wl_err,
double  wl_shift,
int  log_flag,
int  fallback_input_wavecal_flag,
int  keep_higher_degrees_flag,
int  clean_spectrum,
int  display,
double  display_wmin,
double  display_wmax,
int  zp_order,
int  grat1_order,
cpl_propertylist **  qcs,
cpl_table **  lines_diagnostics,
cpl_table **  extracted_out,
cpl_table **  trace_wave_out 
)

Apply the Wavelength Calibration.

Parameters
tw_inTrace wave table
spectra_tabExtracted Spectra
catalog_frameCatalog frame or NULL if not needed
reduce_orderThe order to compute (-1 for all)
reduce_traceThe trace to compute (-1 for all)
wavecal_typeCR2RES_XCORR/LINE1D/LINE2D/ETALON
degreePolynomial degree in dispersion direction
xdegreePolynomial degree in cross-dispersion direction
wl_startWL estimate of the first pixel
wl_endWL estimate of the last pixel
wl_errWL error
wl_shiftwavelength shift to apply
log_flagFlag to apply a log() to the lines intensities
fallback_input_wavecal_flagFlag to copy the input WL to the output when they are not computed
keep_higher_degrees_flagFlag to use higher polynomial degrees from the guess
clean_spectrumRemove the lines that are not in the catalog (1d)
displayFlag to enable display functionalities
display_wminMinimum Wavelength to display or -1.0
display_wmaxMaximum Wavelength to display or -1.0
qcs[out] QC parameters
lines_diagnostics[out] lines diagnostics table
extracted_out[out] extracted table with updated WL
trace_wave_out[out] trace wave table
Returns
0 if ok, -1 otherwise

Definition at line 150 of file cr2res_wave.c.

References cr2res_convert_poly_to_array(), cr2res_etalon_wave_2d(), cr2res_extract_EXTRACT1D_get_spectrum(), cr2res_get_trace_table_index(), cr2res_get_trace_wave_poly(), cr2res_order_real_to_idx(), cr2res_qc_wave_central(), cr2res_qc_wave_disp(), cr2res_wave_1d(), cr2res_wave_2d(), cr2res_wave_estimate_compute(), and cr2res_wave_poly_2d_to_1d().

◆ cr2res_wave_estimate_compute()

cpl_polynomial * cr2res_wave_estimate_compute ( double  wmin,
double  wmax 
)

Compute the polynomial from boundaries.

Parameters
wminFirst pixel wavelength
wmaxLast pixel wavelength
Returns
the polynomial or NULL in error case

The returned polynomial must be deallocated with cpl_polynomial_delete()

Definition at line 1481 of file cr2res_wave.c.

Referenced by cr2res_wave_apply().

◆ cr2res_wave_etalon()

cpl_polynomial * cr2res_wave_etalon ( cpl_bivector *  spectrum,
const cpl_bivector *  spectrum_err,
cpl_polynomial *  wavesol_init,
int  degree,
cpl_array **  wavelength_error 
)

Find solution from etalon.

Parameters
spectrumInput spectrum: etalon
wavesol_initStarting wavelength solution
wavelength_error[out] array of wave_mean_error, wave_max_error
Returns
Wavelength solution, i.e. polynomial that translates pixel values to wavelength.

This function uses the intrinsic property of the etalon spectrum, that the lines are supposed to be equi-spaced, to refine the wavelength solution. The input solution needs to be good enough for the zero-point since the etalon spectrum carries no information about the absolute wavelength scale.

The method involves these steps: Identify lines (thresholding) Determine line centers (gauss fit) Subtract x-coords from subsequent lines, i.e. measure d many times d is the distance between fringes in pixels, assumed to be constant Determine mis-counts in fringes, by looking at outliers in d-distribution, re-count with e.g. half-distance between two fringes, if one is missing. Fit the d-distribution to determine "true d" Use d to calculate new x-coordinates, use input-solution to translate to wavelength Fit measured x-coords to new lambdas to get new solution.

Definition at line 1405 of file cr2res_wave.c.

References cr2res_polynomial_eval_vector().

Referenced by cr2res_wave_1d().

◆ cr2res_wave_gen_wave_map()

hdrl_image * cr2res_wave_gen_wave_map ( const cpl_table *  trace_wave)

Compute the wavelength map from the trace_wave table.

Parameters
trace_waveThe trace wave table
Returns
the wave_map image or NULL in error case

The returned image must be deallocated with hdrl_image_delete()

Definition at line 1608 of file cr2res_wave.c.

References cr2res_convert_array_to_poly(), hdrl_image_get_image(), and hdrl_image_new().

◆ cr2res_wave_get_estimate()

cpl_array * cr2res_wave_get_estimate ( const char *  filename,
int  detector,
int  order 
)

Compute the wavelength estimate.

Parameters
filenameFilename for header WMIN/WMAX
detectorthe detector number
orderthe oder number
Returns
the array with three polynomial coeffs, or NULL in error case

The first 2 coefficients are those from the polynomial poly such that: wmin = poly(1) wmax = poly((CR2RES_DETECTOR_SIZE)

The third coefficient is (b3-b1)/2 where bn is the b wmax-wmin/2048 for detector n in the current order

The returned array must be deallocated with cpl_array_delete()

Definition at line 1527 of file cr2res_wave.c.

References cr2res_io_get_ext_idx(), cr2res_pfits_get_wend(), and cr2res_pfits_get_wstrt().

Referenced by cr2res_trace_add_extra_columns().

◆ cr2res_wave_poly_2d_to_1d()

cpl_polynomial * cr2res_wave_poly_2d_to_1d ( cpl_polynomial *  poly_2d,
int  order 
)

Create a 1D Wavelength Polynomial out of a 2D one.

Parameters
poly_2dThe 2d polynomial (order, lambda)
orderThe order where the 1d polynomial is to be computed
Returns
The newly allocated 1d polynomial or NULL in error case

Definition at line 1776 of file cr2res_wave.c.

Referenced by cr2res_etalon_wave_2d(), cr2res_wave_2d(), and cr2res_wave_apply().

◆ cr2res_wave_polys_1d_to_2d()

cpl_polynomial * cr2res_wave_polys_1d_to_2d ( cpl_polynomial **  poly_1ds,
int *  orders,
int  npolys,
cpl_size  xdegree 
)

Create a 2D Wavelength Polynomial out of a several 1D ones.

Parameters
poly_1dsList of 1d wavelength polynomials
ordersList of the order values
npolysNumber of passed polynomials
xdegreePolynomial degree in vertical direction (cross-dispersed)
Returns
The newly allocated 2D polynomial or NULL in error case

Definition at line 1685 of file cr2res_wave.c.

◆ cr2res_wave_xcorr()

cpl_polynomial * cr2res_wave_xcorr ( cpl_bivector *  spectrum,
cpl_polynomial *  wavesol_init,
double  wl_error,
cpl_bivector *  lines_list,
int  degree,
int  keep_higher_degrees_flag,
int  cleaning_filter_size,
double  slit_width,
double  fwhm,
int  display,
double *  best_xcorr,
cpl_array **  wavelength_error 
)

Find solution by cross-correlating template spectrum.

Parameters
spectrumInput spectrum
wavesol_initStarting wavelength solution
wl_errorMax error in nm of the initial guess
lines_listLines List (flux, wavelengths)
degreeThe polynomial degree
keep_higher_degrees_flag
cleaning_filter_size
slit_width
fwhm
displayValue matching the pass to display (0 for none, )
best_xcorr[out] Best Cross correlation factor
wavelength_error[out] array of wave_mean_error, wave_max_error
Returns
Wavelength solution, i.e. polynomial that translates pixel values to wavelength.

TODO: Summarize method

Definition at line 1221 of file cr2res_wave.c.

Referenced by cr2res_wave_1d().