Recipe: Order Position

Defines

#define USE_PPM   0
#define FIT_SLOPE   1
#define WEIGHTED_FIT   1

Functions

static cpl_error_code verify_calibration (const cpl_table *selected, const cpl_table *linetable, double TOLERANCE, double red_chisq)
 Report quality of calibration.
static cpl_error_code compute_lambda (cpl_table *linetable, const polynomial *dispersion_relation, const polynomial *dispersion_variance, bool verbose)
 Apply dispersion relation to line table.
static int identify_lines (cpl_table *linetable, const cpl_table *line_refer, double ALPHA)
 Identify lines by comparing to catalogue wavelengths.
static polynomialcalibrate_global (const cpl_table *linetable, cpl_table **selected, int degree, bool verbose, bool reject, double TOLERANCE, double kappa, double *red_chisq, polynomial **dispersion_variance, double *pixelsize, double *rms_wlu, double *rms_pixels)
 Create a fit of all orders.
polynomialuves_wavecal_identify (cpl_table *linetable, const cpl_table *line_refer, const polynomial *guess_dispersion, int DEGREE, double TOLERANCE, double ALPHA, double MAXERROR, double kappa)
 Obtain final dispersion relation.
int uves_wavecal_identify_lines_ppm (cpl_table *linetable, const cpl_table *line_refer)
 Identify lines using point pattern matching.
static double xcenter (const cpl_image *image, const cpl_image *noise, int xlo, int xhi, int row, centering_method CENTERING_METHOD, int bin_disp, double *sigma, double *intensity, double *dx0, double *slope, double *background)
 Refine the center position of an initially detected emission line.
static cpl_error_code detect_lines (const cpl_image *spectrum, const cpl_image *noise, const uves_propertylist *spectrum_header, bool flat_fielded, int RANGE, double THRESHOLD, centering_method CENTERING_METHOD, int bin_disp, const polynomial *order_locations, cpl_image *arcframe, cpl_table *linetable, int *ndetected, int *nrows)
 Find emission lines above a certain threshold.
cpl_table * uves_wavecal_search (const cpl_image *spectrum, const cpl_image *noise, const uves_propertylist *spectrum_header, bool flat_fielded, const polynomial *order_locations, cpl_image *arcframe, int RANGE, int MINLINES, int MAXLINES, centering_method CENTERING_METHOD, int bin_disp)
 Search for a given number of emission lines.
lt_typeuves_lt_new (int windows, int traces)
 Allocate line table.
void uves_lt_delete (lt_type **lt)
 Deallocate line table.
cpl_table ** uves_lt_get_table (const lt_type *lt, int window, int trace)
 Get the table structure.
polynomial ** uves_lt_get_disprel (const lt_type *lt, int window, int trace)
 Get dispersion relation.
polynomial ** uves_lt_get_absord (const lt_type *lt, int window, int trace)
 Get absolute order polynomial.
int * uves_lt_get_firstabs (const lt_type *lt, int window, int trace)
 Get first absolute order.
int * uves_lt_get_lastabs (const lt_type *lt, int window, int trace)
 Get last absolute order.
int uves_wavecal_find_nearest (const cpl_table *line_refer, double lambda, int lo, int hi)
 Find best matching catalogue wavelength.
int uves_delete_bad_lines (cpl_table *table, double TOLERANCE, double kappa)
 Delete bad lines from line table.
cpl_error_code uves_draw_lines (cpl_image *image, polynomial *dispersion, const polynomial *order_locations, const cpl_table *t, const char *lambda_column, const char *abs_order, const int *relative_order, int minorder, int maxorder, bool vertical, int offset)
 Draw lines in an echelle image.

Detailed Description

This recipe determines the echelle order locations. See man-page for details.


Function Documentation

static cpl_error_code verify_calibration ( const cpl_table *  selected,
const cpl_table *  linetable,
double  TOLERANCE,
double  red_chisq 
) [static]

Report quality of calibration.

Parameters:
linetable The line table
TOLERANCE When reporting the RMS of the fit, exclude lines with residuals worse than TOLERANCE from the computation. If positive, this tolerance is considered in pixel units, otherwise in w.l.u.
red_chisq Reduced chi^2 of the calibration
Returns:
CPL_ERROR_NONE iff okay

This function reports the RMS of a wavelength calibration, as well as the percentage of the brightest lines that were identified.

Definition at line 472 of file uves_wavecal_identify.c.

References check, uves_msg, and uves_msg_warning.

Referenced by uves_wavecal_identify().

static cpl_error_code compute_lambda ( cpl_table *  linetable,
const polynomial dispersion_relation,
const polynomial dispersion_variance,
bool  verbose 
) [static]

Apply dispersion relation to line table.

Parameters:
linetable The line table
dispersion_relation The dispersion relation
dispersion_variance The variance of the dispersion relation
verbose print warning if dl/dx is negative?
Returns:
CPL_ERROR_NONE iff okay

This function applies the provided dispersion relation to the line table, thereby calculating the predicted wavelengths for lines in the spectrum, the residual of the fit, the local pixelsize in w.l.u., and if dispersion_variance is non-NULL, the uncertainty of the fitted wavelength.

Definition at line 551 of file uves_wavecal_identify.c.

References check, passure, uves_msg_debug, uves_msg_warning, uves_polynomial_derivative_2d(), uves_polynomial_evaluate_2d(), and uves_polynomial_get_dimension().

Referenced by calibrate_global(), and uves_wavecal_identify().

static int identify_lines ( cpl_table *  linetable,
const cpl_table *  line_refer,
double  ALPHA 
) [static]

Identify lines by comparing to catalogue wavelengths.

Parameters:
linetable The line table containing the line positions and predicted wavelengths
line_refer The wavelength catalogue
ALPHA Parameter < 1 to control distance to nearest neighbour
Returns:
The total number of lines identified (including previous identifications), or undefined on error.

This function identifies lines in the provided linetable by comparing to the reference table (the wavelength catalogue). The identified wavelengths are written to the "Ident" column of the line table.

An identification is made iff

  • The catalogue wavelength is within two linewidths of the computed (predicted) wavelength: | lambda_cat - lambda_com | < 2 * sigma, where sigma is the detected line width,
  • The nearest neighbour (in the spectrum and in the catalogue) is farther away than the candidate catalogue wavelength (after multiplying by the "safety parameter", ALPHA): | lambda_cat - lambda_nn| * ALPHA > | lambda_cat - lambda_com |,
  • The nearest neighbour (in the spectrum and in the catalogue) is farther away than the average tolerance distance, which measures the precision of the identifications: tolerance < ALPHA * | lambda_cat - lambda_nn| . See code for the exact definition of the tolerance .

The purpose of the first criterion is to make the correct identifications. The purpose of the latter two criterions is to avoid making incorrect identifications.

If a line was previously identified (implied by a valid "Ident" column value) but now fails to meet the ID criterium, it is not deleted.

Definition at line 700 of file uves_wavecal_identify.c.

References assure_mem, passure, uves_max_double(), uves_max_int(), uves_min_double(), uves_min_int(), uves_msg_debug, uves_round_double(), and uves_wavecal_find_nearest().

Referenced by uves_wavecal_identify().

static polynomial * calibrate_global ( const cpl_table *  linetable,
cpl_table **  selected,
int  degree,
bool  verbose,
bool  reject,
double  TOLERANCE,
double  kappa,
double *  red_chisq,
polynomial **  dispersion_variance,
double *  pixelsize,
double *  rms_wlu,
double *  rms_pixels 
) [static]

Create a fit of all orders.

Parameters:
linetable The line table
selected (output) if non-NULL, subset of linetable containing the lines which were used in the final fit
degree Degree of both independent variables of polynomial fit
verbose Be verbose about autodegree fitting?
reject Do rejection?
TOLERANCE Before fitting, exclude lines with residuals worse than TOLERANCE. If positive, this tolerance is considered in pixel units, otherwise in w.l.u.
kappa used for removing outliers
red_chisq If non-NULL, the reduced chi square of the fit.
dispersion_variance If non-NULL, the variance of the fit returned polynomial.
pixelsize If non-NULL, the average of d(lambda)/dx
rms_wlu If non-NULL, the root-mean-square residual (w.l.u)
rms_pixels If non-NULL, the root-mean-square residual (pixels)
Returns:
The obtained dispersion relation in the form lambda * m = f(x, m), or NULL on error
Note:
Un-identified lines and lines with residuals larger then TOLERANCE (from the previous fit) are excluded from the fit.

Definition at line 1065 of file uves_wavecal_identify.c.

References assure_mem, check, check_nomsg, compute_lambda(), passure, uves_delete_bad_lines(), uves_msg_debug, uves_polynomial_delete(), uves_polynomial_regression_2d(), and uves_polynomial_regression_2d_autodegree().

Referenced by uves_wavecal_identify().

polynomial* uves_wavecal_identify ( cpl_table *  linetable,
const cpl_table *  line_refer,
const polynomial guess_dispersion,
int  DEGREE,
double  TOLERANCE,
double  ALPHA,
double  MAXERROR,
double  kappa 
)

Obtain final dispersion relation.

Parameters:
linetable The line table containing the already detected emission lines
line_refer The wavelength catalogue
guess_dispersion The initial dispersion relation in the form lambda * m = f(x, m)
DEGREE Degree of both independent variables of the dispersion polynomial.
TOLERANCE When making the final fit, lines with residuals worse than TOLERANCE are excluded. If positive, this tolerance is considered in pixel units, otherwise in w.l.u.
ALPHA Parameter that controls the next-neighbour distance during line identification. See identify_lines() for details.
MAXERROR If the RMS of the fit is larger than this number (in pixels), the calibration loop terminates with an error. This is to ensure a graceful exit (when incorrect identifications are made)
kappa used in outlier rejectiong
Returns:
The obtained initial dispersion solution of the form lambda * m = f(x, m), or NULL on error

This function performs a wavelength calibration of the detected lines listed in the provided linetable, starting from the guess_dispersion solution. Results of the calibration are written to the linetable.

The function will iteratively identify as many of the detected lines as possible (using the specified wavelength catalogue, see also identify_lines() ), then update the fit polynomial (see also calibrate_global() ). This loop continues until no new line identifications can be made. After this first convergence all identifications are reset (to remove possible false identifications), and the loop repeats, but this time ignoring lines with residuals worse than TOLERANCE . The final solution is returned.

Definition at line 217 of file uves_wavecal_identify.c.

References calibrate_global(), check, compute_lambda(), identify_lines(), passure, uves_msg, uves_msg_debug, uves_polynomial_delete(), uves_wavecal_identify_lines_ppm(), and verify_calibration().

int uves_wavecal_identify_lines_ppm ( cpl_table *  linetable,
const cpl_table *  line_refer 
)

Identify lines using point pattern matching.

Parameters:
linetable The line table containing the line positions
line_refer The wavelength catalogue
Returns:
number of identifications

Definition at line 1198 of file uves_wavecal_identify.c.

References check_nomsg, uves_error_reset, uves_msg_debug, uves_msg_warning, and uves_round_double().

Referenced by test_ppm(), and uves_wavecal_identify().

static double xcenter ( const cpl_image *  image,
const cpl_image *  noise,
int  xlo,
int  xhi,
int  row,
centering_method  CENTERING_METHOD,
int  bin_disp,
double *  sigma,
double *  intensity,
double *  dx0,
double *  slope,
double *  background 
) [static]

Refine the center position of an initially detected emission line.

Parameters:
image The input image
noise The noise (1 sigma) associated with the input image
xlo Left x-coordinate of search window
xhi Right x-coordinate of search window
row Image row to look at
CENTERING_METHOD The method used to calculate the position
bin_disp CCD binning, used to finetune size of fitting window
sigma (output) Peak width
intensity (output) Peak height (relative to 0)
dx0 (output) The statistical uncertainty of peak center
slope (output) fitted background slope
background (output) fitted background
Returns:
The refined center location, or undefined on error

This function locates the center of an emission peak in the image window (xlo, y)-(xhi, y) by calculating guess values for the center, the stddev, the background level and the height of the peak, then calling uves_fit_gaussian_1d_image() with the "guess" values.

The centering_method indicates whether to use the 'gravity' method or gaussian fitting and should always be set to CENTERING_GAUSSIAN .

Definition at line 663 of file uves_wavecal_search.c.

References passure, uves_error_reset, uves_gauss(), uves_gauss_derivative(), uves_gauss_linear(), uves_gauss_linear_derivative(), uves_max_int(), uves_min_int(), and uves_msg_debug.

Referenced by detect_lines(), and trace_order().

static cpl_error_code detect_lines ( const cpl_image *  spectrum,
const cpl_image *  noise,
const uves_propertylist spectrum_header,
bool  flat_fielded,
int  RANGE,
double  THRESHOLD,
centering_method  CENTERING_METHOD,
int  bin_disp,
const polynomial order_locations,
cpl_image *  arcframe,
cpl_table *  linetable,
int *  ndetected,
int *  nrows 
) [static]

Find emission lines above a certain threshold.

Parameters:
spectrum The spectrum to search.
noise The noise (1 sigma) of the spectrum
flat_fielded was the spectrum flat-fielded?
spectrum_header Header of spectrum image. The (relative) order numbers are read from this header
RANGE Width (in pixels) of search window is 2 * RANGE + 1
THRESHOLD Detection threshold (relative to local background level)
CENTERING_METHOD The method used to calculate the line position.
bin_disp binning in dispersion direction
arcframe The (raw) arc frame image. Iff non-NULL, the results of the line search will be drawn on this image using the provided order_locations.
order_locations The order position bi-polynomial (defines the order locations as a function of x and relative order number)
linetable The line table where detected emission lines will be written
ndetected (output) Number of lines detected
nrows (output) Number of lines written to the line table
Returns:
CPL_ERROR_NONE iff okay

This function searches for emission lines at each row of the input image (an extracted spectrum).

The detected lines will be written to the columns "Y" (the relative order number), "X" (the line position), "Xwidth" (the detected line width (stddev in pixels)) and "Peak" (the detected emission line intensity).

A line is detected iff the flux at a position is more than THRESHOLD plus the local background level, which is defined as the median of a window of width (2 * RANGE + 1) centered on the current position. If flat_fielding was done THRESHOLD is the peak height in units of error bars.

After the line is initially detected, the line position is calculated depending on the CENTERING_METHOD. See also xcenter() .

The number of detected lines, ndetected, will be different from the number of rows written to the line table, nrows, only if the table is not large enough to contain all the detected lines.

Finally, doublets (i.e. lines with positions within 2.0 pixels) are removed from the line table. The maximum number of lines to search for is defined by the number of available rows in the line table.

If the provided arcframe is not NULL, the results of the line search will be marked on this image (used for debugging purposes).

Definition at line 393 of file uves_wavecal_search.c.

References check, passure, uves_max_int(), uves_min_int(), uves_msg, uves_msg_debug, uves_pfits_get_crval2(), uves_polynomial_evaluate_2d(), uves_round_double(), uves_tostring_cpl_type(), and xcenter().

Referenced by uves_wavecal_search().

cpl_table* uves_wavecal_search ( const cpl_image *  spectrum,
const cpl_image *  noise,
const uves_propertylist spectrum_header,
bool  flat_fielded,
const polynomial order_locations,
cpl_image *  arcframe,
int  RANGE,
int  MINLINES,
int  MAXLINES,
centering_method  CENTERING_METHOD,
int  bin_disp 
)

Search for a given number of emission lines.

Parameters:
spectrum The spectrum to search.
noise The noise (1 sigma) of the spectrum.
spectrum_header Header of spectrum image. The (relative) order numbers are read from this header
flat_fielded was the arclamp frame flat-fielded?
order_locations The order position bi-polynomial (defines the order locations as a function of x and relative order number)
arcframe The raw input image (before spectrum extraction). The results of the line search will be drawn on this image.
RANGE Width (in pixels) of search window is 2 * RANGE + 1
MINLINES The minimum number of lines to search for
MAXLINES The maximum number of lines to search for
CENTERING_METHOD The method used to calculate the line position.
Returns:
A newly allocated line table containing the detected lines, or NULL on error

This function tries to detect a number of emission lines in the already extracted spectrum. On success the number of detected lines is in the range from MINLINES to MAXLINES (both inclusive). This is achieved by calling detect_lines() . with varying detection threshold levels.

The result line table contains the columns "Y" (the relative order number), "X" (the line positions), "Xwidth" (the detected line widths (stddev in pixels)) and "Peak" (the detected emission line intensities).

Definition at line 184 of file uves_wavecal_search.c.

References check, detect_lines(), passure, uves_msg, and uves_tostring_cpl_type().

lt_type* uves_lt_new ( int  windows,
int  traces 
)

Allocate line table.

Parameters:
windows Number of extraction windows per trace
traces Number of traces
Returns:
newly allocated line table which must be deallocated using uves_lt_delete().

Definition at line 113 of file uves_wavecal_utils.c.

References lt_type::absolute_order, assure_mem, lt_type::dispersion_relation, lt_type::first_absolute_order, lt_type::last_absolute_order, lt_type::table, lt_type::traces, and lt_type::windows.

void uves_lt_delete ( lt_type **  lt  ) 

Deallocate line table.

Parameters:
lt the line table to deallocate

Definition at line 141 of file uves_wavecal_utils.c.

References uves_polynomial_delete().

cpl_table** uves_lt_get_table ( const lt_type lt,
int  window,
int  trace 
)

Get the table structure.

Parameters:
lt the line table
window window number
trace trace number
Returns:
the CPL table of this 'line table'

Definition at line 175 of file uves_wavecal_utils.c.

References lt_type::table, and lt_type::traces.

polynomial** uves_lt_get_disprel ( const lt_type lt,
int  window,
int  trace 
)

Get dispersion relation.

Parameters:
lt the line table
window window number
trace trace number
Returns:
the dispersion relation

Definition at line 190 of file uves_wavecal_utils.c.

References lt_type::dispersion_relation, and lt_type::traces.

polynomial** uves_lt_get_absord ( const lt_type lt,
int  window,
int  trace 
)

Get absolute order polynomial.

Parameters:
lt the line table
window window number
trace trace number
Returns:
the absolute order polynomial

Definition at line 204 of file uves_wavecal_utils.c.

References lt_type::absolute_order, and lt_type::traces.

int* uves_lt_get_firstabs ( const lt_type lt,
int  window,
int  trace 
)

Get first absolute order.

Parameters:
lt the line table
window window number
trace trace number
Returns:
the first absolute order

Definition at line 218 of file uves_wavecal_utils.c.

References lt_type::first_absolute_order, and lt_type::traces.

int* uves_lt_get_lastabs ( const lt_type lt,
int  window,
int  trace 
)

Get last absolute order.

Parameters:
lt the line table
window window number
trace trace number
Returns:
the last absolute order

Definition at line 232 of file uves_wavecal_utils.c.

References lt_type::last_absolute_order, and lt_type::traces.

int uves_wavecal_find_nearest ( const cpl_table *  line_refer,
double  lambda,
int  lo,
int  hi 
)

Find best matching catalogue wavelength.

Parameters:
line_refer Table to search
lambda Find catalogue wavelength nearest to this wavelength
lo First row (inclusive, counting from 0) of table to search
hi Last row (inclusive, counting from 0) of table to search
Returns:
The index of the row containing the wavelength closest to lambda

This function returns the index of the nearest wavelength in the range {lo, ..., hi} The input table must be sorted according to the column "Wave" which is the column that is searched.

Note:
The function implements a binary search (using recursion) which is why the input table must be sorted.

Definition at line 255 of file uves_wavecal_utils.c.

References uves_wavecal_find_nearest().

Referenced by identify_lines(), and uves_wavecal_find_nearest().

int uves_delete_bad_lines ( cpl_table *  table,
double  TOLERANCE,
double  kappa 
)

Delete bad lines from line table.

Parameters:
table The line table containing the bad lines to be removed
TOLERANCE If positive, this tolerance is considered in pixel units, otherwise in w.l.u.
kappa Used in kappa sigma clipping
Returns:
The number of deleted bad lines (i.e. deleted table rows), or undefined on error

Rows with invalid "Ident" (un-identified lines) are removed, as are rows with "Residual_pix" worse than the specified TOLERANCE. If the TOLERANCE is negative, rows with LINETAB_RESIDUAL worse than |TOLERANCE | are removed.

Definition at line 307 of file uves_wavecal_utils.c.

References check, check_nomsg, and uves_average_reject().

Referenced by calibrate_global().

cpl_error_code uves_draw_lines ( cpl_image *  image,
polynomial dispersion,
const polynomial order_locations,
const cpl_table *  t,
const char *  lambda_column,
const char *  abs_order,
const int *  relative_order,
int  minorder,
int  maxorder,
bool  vertical,
int  offset 
)

Draw lines in an echelle image.

Parameters:
image Image to draw on
dispersion The dispersion relation
order_locations The order locations (relative order numbering)
t Table containing wavelengths to plot
lambda_column Name of column containing the wavelengths to plot
abs_order Name of column containing the absolute order numbers
relative_order The map from absolute to relative orders
minorder The first absolute order used
maxorder The last absolute order used
vertical Flag indicating if plotted lines are vertical (true) or horizontal (false)
offset Horizontal offset (in pixels) of the plotted line
Returns:
CPL_ERROR_NONE iff OK

This function is used for debugging. It plots a line (length = 7 pixels) on the image corresponding to each wavelength listed in the specified table column.

If abs_order is non-NULL, each wavelength is drawn on the order specified by abs_order. In this case minorder and maxorder are ignored.

If abs_order is NULL, each wavelength is drawn on every order from minorder (inclusive) to maxorder (inclusive). This mode can be used to plot catalogue wavelengths that a priori are located in no specific order.

The y position of a wavelength is inferred from the polynomial order_locations (the order table), and the x position is inferred from the dispersion relation.

Definition at line 407 of file uves_wavecal_utils.c.

References check, passure, uves_error_reset, uves_max_int(), uves_min_int(), uves_msg_debug, uves_polynomial_evaluate_2d(), uves_polynomial_get_dimension(), and uves_polynomial_solve_2d().


Generated on 8 Mar 2011 for UVES Pipeline Reference Manual by  doxygen 1.6.1