MUSE Pipeline Reference Manual  2.1.1
Data Structures | Typedefs | Enumerations | Functions | Variables
Utilities

Data Structures

struct  muse_cpl_optimize_control_t
 Optimization control parameters. More...
 

Typedefs

typedef cpl_error_code( muse_cpl_evaluate_func) (void *aData, cpl_array *aPar, cpl_array *aRetval)
 User provided function to evaluate in muse_cpl_optimize_lvmq(). More...
 

Enumerations

Functions

static void muse_cpl_optimize_print (const char *func, int n_par, double *par, int m_dat, double *fvec, void *data, int iflag, int iter, int nfev)
 Default printout method adopted to CPL. More...
 
cpl_error_code muse_cpl_optimize_lvmq (void *aData, cpl_array *aPar, int aSize, muse_cpl_evaluate_func *aFunction, muse_cpl_optimize_control_t *aCtrl)
 Minimize a function with the Levenberg-Marquardt algorithm. More...
 
const char * muse_get_license (void)
 Get the pipeline copyright and license. More...
 
unsigned char muse_utils_get_ifu (const cpl_propertylist *aHeaders)
 Find out the IFU/channel from which this header originated. More...
 
int muse_utils_get_extension_for_ifu (const char *aFilename, unsigned char aIFU)
 Return extension number that corresponds to this IFU/channel number. More...
 
cpl_frameset * muse_frameset_find (const cpl_frameset *aFrames, const char *aTag, unsigned char aIFU, cpl_boolean aInvert)
 return frameset containing data from an IFU/channel with a certain tag More...
 
cpl_frameset * muse_frameset_find_tags (const cpl_frameset *aFrames, const cpl_array *aTags, unsigned char aIFU, cpl_boolean aInvert)
 return frameset containing data from an IFU/channel with the given tag(s) More...
 
cpl_frameset * muse_frameset_check_raw (const cpl_frameset *aFrames, const cpl_array *aTags, unsigned char aIFU)
 return frameset containing good raw input data More...
 
cpl_frameset * muse_frameset_sort_raw_other (const cpl_frameset *aFrames, int aIndex, const char *aDateObs, cpl_boolean aSequence)
 Create a new frameset containing all relevant raw frames first then all other frames. More...
 
cpl_frame * muse_frameset_find_master (const cpl_frameset *aFrames, const char *aTag, unsigned char aIFU)
 find the master frame according to its CCD number and tag More...
 
cpl_error_code muse_utils_frameset_merge_frames (cpl_frameset *aFrames, cpl_boolean aDelete)
 Merge IFU-specific files from a frameset to create multi-IFU outputs. More...
 
muse_tablemuse_table_load_filter (muse_processing *aProcessing, const char *aFilterName)
 Load a table for a given filter name. More...
 
double muse_utils_filter_fraction (const muse_table *aFilter, double aLambda1, double aLambda2)
 Compute fraction of filter area covered by the data range. More...
 
cpl_error_code muse_utils_filter_copy_properties (cpl_propertylist *aHeader, const muse_table *aFilter, double aFraction)
 Add/propagate filter properties to header of collapsed image. More...
 
char * muse_utils_header_get_lamp_names (cpl_propertylist *aHeader, char aSep)
 Concatenate names of all active calibration lamps. More...
 
cpl_array * muse_utils_header_get_lamp_numbers (cpl_propertylist *aHeader)
 List numbers of all active calibration lamps. More...
 
cpl_matrix * muse_matrix_new_gaussian_2d (int aXHalfwidth, int aYHalfwidth, double aSigma)
 Create a matrix that contains a normalized 2D Gaussian. More...
 
cpl_image * muse_utils_image_fit_polynomial (const cpl_image *aImage, unsigned short aXOrder, unsigned short aYOrder)
 Create a smooth version of a 2D image by fitting it with a 2D polynomial. More...
 
cpl_error_code muse_utils_image_get_centroid_window (cpl_image *aImage, int aX1, int aY1, int aX2, int aY2, double *aX, double *aY, muse_utils_centroid_type aBgType)
 Compute centroid over an image window, optionally marginalizing over the background. More...
 
cpl_error_code muse_utils_get_centroid (const cpl_matrix *aPositions, const cpl_vector *aValues, const cpl_vector *aErrors, double *aX, double *aY, muse_utils_centroid_type aBgType)
 Compute centroid of a two-dimensional dataset. More...
 
cpl_error_code muse_utils_fit_multigauss_1d (const cpl_vector *aX, const cpl_bivector *aY, cpl_vector *aCenter, double *aSigma, cpl_vector *aFlux, cpl_vector *aPoly, double *aMSE, double *aRedChisq, cpl_matrix **aCovariance)
 Carry out a multi-Gaussian fit of data in a vector. More...
 
cpl_error_code muse_utils_fit_moffat_2d (const cpl_matrix *aPositions, const cpl_vector *aValues, const cpl_vector *aErrors, cpl_array *aParams, cpl_array *aPErrors, const cpl_array *aPFlags, double *aRMS, double *aRedChisq)
 Fit a 2D Moffat function to a given set of data. More...
 
cpl_polynomial * muse_utils_iterate_fit_polynomial (cpl_matrix *aPos, cpl_vector *aVal, cpl_vector *aErr, cpl_table *aExtra, const unsigned int aOrder, const double aRSigma, double *aMSE, double *aChiSq)
 Iterate a polynomial fit. More...
 
double muse_utils_pixtable_fit_line_gaussian (muse_pixtable *aPixtable, double aLambda, double aHalfWidth, double aBinSize, float aLo, float aHi, unsigned char aIter, cpl_array *aResults, cpl_array *aErrors)
 Fit a 1D Gaussian to a given wavelength range in a pixel table. More...
 
cpl_error_code muse_utils_copy_modified_header (cpl_propertylist *aH1, cpl_propertylist *aH2, const char *aKeyword, const char *aString)
 Copy a modified header keyword from one header to another. More...
 
cpl_error_code muse_utils_set_hduclass (cpl_propertylist *aHeader, const char *aClass2, const char *aExtData, const char *aExtDQ, const char *aExtStat)
 Set HDU headers for the ESO FITS data format. More...
 
void muse_utils_memory_dump (const char *aMarker)
 Display the current memory usage of the given program. More...
 
cpl_image * muse_convolve_image (const cpl_image *aImage, const cpl_matrix *aKernel)
 Compute the convolution of an image with a kernel. More...
 
muse_imagemuse_fov_load (const char *aFilename)
 Load a FOV image into a MUSE image. More...
 

Variables

const muse_cpltable_def muse_filtertable_def []
 MUSE filter table definition. More...
 

Detailed Description

This group handles several utility functions.

Typedef Documentation

typedef cpl_error_code( muse_cpl_evaluate_func) (void *aData, cpl_array *aPar, cpl_array *aRetval)

User provided function to evaluate in muse_cpl_optimize_lvmq().

Parameters
aDataData forwarded from the fit algorithm call.
aParCurrent fit parameter
aRetvalReturn value vector
Returns
CPL_ERROR_NONE if everything went OK. Any other value is considered as error and will stop the minimization muse_cpl_optimize_lvmq().
Sample prototype:
cpl_error_code
muse_evaluate_sample(void *aData, cpl_array *aPar, cpl_array *aRetval);

Definition at line 88 of file muse_optimize.h.

Enumeration Type Documentation

Background handling when computing centroids.

Enumerator
MUSE_UTILS_CENTROID_NORMAL 

no background treatment

MUSE_UTILS_CENTROID_MEAN 

subtract average as background

MUSE_UTILS_CENTROID_MEDIAN 

subtract median value as background

Definition at line 102 of file muse_utils.h.

Function Documentation

cpl_image* muse_convolve_image ( const cpl_image *  aImage,
const cpl_matrix *  aKernel 
)

Compute the convolution of an image with a kernel.

Parameters
aImageThe image to convolve.
aKernelThe convolution kernel.
Returns
On success the convoluted image is returned, and NULL if an error occurred.

The function convolves the input image aImage with the convolution kernel aKernel using an FFT. If aImage has an odd number of pixels along the y-axis, the last image row is discarded to make it a size even, however, the size of aImage along the x-axis must be even and an error is returned if it is not.

Exceptions
set CPL_ERROR_NULL_INPUT, return NULLa parameter is NULL
set CPL_ERROR_IVALID_TYPE, return NULLaImage is not if type double
set CPL_ERROR_INCOMPATIBLE_INPUT, return NULLinvalid image size

Definition at line 2748 of file muse_utils.c.

Referenced by muse_find_stars().

cpl_error_code muse_cpl_optimize_lvmq ( void *  aData,
cpl_array *  aPar,
int  aSize,
muse_cpl_evaluate_func aFunction,
muse_cpl_optimize_control_t aCtrl 
)

Minimize a function with the Levenberg-Marquardt algorithm.

Parameters
aDataAny data to be forwarded to the evaluation function.
aParArray containing the fits parameters.
aSizeNumber of data points.
aFunctionFunction that evaluates the fit.
aCtrlStructure that controls the algorithm, or NULL for defaults.
Return values
CPL_ERROR_NONEEverything went OK.
CPL_ERROR_NULL_INPUT(aPar != NULL) not fulfilled
CPL_ERROR_ILLEGAL_INPUT(cpl_array_get_size(aPar) > 0) not fulfilled
CPL_ERROR_ILLEGAL_INPUT(aSize > 0) not fulfilled
CPL_ERROR_NULL_INPUT(aFunction != NULL) not fulfilled
CPL_ERROR_SINGULAR_MATRIXReturn vector is orthogonal to the columns of the jacobian
CPL_ERROR_CONTINUEThe algorithm did not converge.
CPL_ERROR_ILLEGAL_OUTPUTMemory allocation failed
Remarks
If the evaluation function did not return CPL_ERROR_NONE, the minimization stops and returns with the return value of the evaluation function.

The interface is loosely modeled after cpl_fit_lvmq(). However it can be used with functions that calc the whole result vector at once, it does not need to provide a dervitative, and in may forward any data to the evaluation function.

Note
The implementation used here is copied from http://www.messen-und-deuten.de/lmfit/ .

Definition at line 251 of file muse_optimize.c.

References muse_cpl_optimize_control_t::debug, muse_cpl_optimize_control_t::ftol, muse_cpl_optimize_control_t::gtol, muse_cpl_optimize_control_t::maxcall, muse_cpl_optimize_print(), and muse_cpl_optimize_control_t::xtol.

Referenced by muse_lsf_params_fit(), muse_sky_lines_fit(), muse_sky_lines_fit_old(), and muse_utils_fit_moffat_2d().

static void muse_cpl_optimize_print ( const char *  func,
int  n_par,
double *  par,
int  m_dat,
double *  fvec,
void *  data,
int  iflag,
int  iter,
int  nfev 
)
static

Default printout method adopted to CPL.

Parameters
funcFunction name
n_parNumber of parameters
parVector of parameters
m_datNumber of datapoints
fvecVector of datapoints
dataFor soft control of printout behaviour, add control variables to the data struct.
iflag0 (init) 1 (outer loop) 2(inner loop) -1(terminated)
iterouter loop counter
nfevnumber of calls to *evaluate

Definition at line 177 of file muse_optimize.c.

Referenced by muse_cpl_optimize_lvmq().

muse_image* muse_fov_load ( const char *  aFilename)

Load a FOV image into a MUSE image.

Parameters
aFilenamename of the file to load
Returns
a muse_image or NULL on error

Load the DATA extension containing the FOV image data from the input file. Bad pixels may be encoded as NaN values in the image data directly, or, if a DQ extension is found in the input file, it is used to replace the bad pixels in the FOV image data with NaN values.

Alternatively, if the DATA extension is an extension with higher-dimensional data or no image content (NAXIS != 2), load the first 2D image extension. In that case, search for the related STAT and DQ extension by prefixing them with the filter name.

Definition at line 2864 of file muse_utils.c.

References muse_image::data, muse_image::dq, muse_image::header, muse_image_delete(), muse_image_dq_to_nan(), muse_image_new(), muse_pfits_get_bunit(), muse_pfits_get_extname(), muse_pfits_get_naxis(), MUSE_WCS_KEYS, and muse_image::stat.

cpl_frameset* muse_frameset_check_raw ( const cpl_frameset *  aFrames,
const cpl_array *  aTags,
unsigned char  aIFU 
)

return frameset containing good raw input data

Parameters
aFramesthe input frameset
aTagsthe required DFS tag(s) of the frame (optional)
aIFUthe IFU/channel number
Returns
the resulting frameset or NULL on failure

Returns the good raw input data in form of a frameset. If aTag is NULL, tags are not checked.

Definition at line 284 of file muse_utils.c.

References muse_frameset_find_tags(), muse_pfits_get_binx(), muse_pfits_get_biny(), muse_pfits_get_chip_id(), muse_pfits_get_chip_name(), muse_pfits_get_read_id(), muse_pfits_get_read_name(), and muse_utils_get_extension_for_ifu().

Referenced by muse_basicproc_params_delete().

cpl_frameset* muse_frameset_find ( const cpl_frameset *  aFrames,
const char *  aTag,
unsigned char  aIFU,
cpl_boolean  aInvert 
)

return frameset containing data from an IFU/channel with a certain tag

Parameters
aFramesthe input frameset
aTagthe required DFS tag of the frame (optional)
aIFUthe IFU/channel number to check for (optional)
aInvertif true, invert the selection
Returns
the resulting frameset or NULL on failure

Returns the data in form of a frameset. Note that this frameset may be empty (but non-NULL), if no frame matching the parameters was found! If aTag is NULL, tags are not checked. If aIFU is 0, the IFU number is not checked.

If aTag is MUSE_TAG_GEOMETRY_TABLE or MUSE_TAG_TWILIGHT_CUBE, then all frames with that tag are matched, irrespective of the IFU number found it the header of the corresponding file.

The IFU number is checked using muse_utils_get_ifu(). If this fails, then the respective file is expected to be a generic calibration file.

Exceptions
set CPL_ERROR_NULL_INPUT, return NULLaFrames is NULL
propagate error code, continue with next frameheader of a file cannot be loaded

Definition at line 160 of file muse_utils.c.

References muse_pfits_get_pipefile(), muse_utils_get_extension_for_ifu(), and muse_utils_get_ifu().

Referenced by muse_basicproc_params_delete(), muse_frameset_find_master(), muse_frameset_find_tags(), muse_lsf_cube_load_all(), muse_postproc_cube_load_output_wcs(), muse_processing_check_input(), muse_processing_load_mask(), muse_processing_lsf_params_load(), muse_sky_continuum_load(), and muse_sky_lines_load().

cpl_frame* muse_frameset_find_master ( const cpl_frameset *  aFrames,
const char *  aTag,
unsigned char  aIFU 
)

find the master frame according to its CCD number and tag

Parameters
aFramesthe frames list
aTagthe tag
aIFUthe IFU/channel number
Returns
the frame of the master image

The returned frame is duplicated from the input frameset and therefore has to be deallocated after use.

Definition at line 505 of file muse_utils.c.

References muse_frameset_find().

Referenced by muse_basicproc_params_delete(), muse_processing_load_header(), muse_processing_load_table(), and muse_table_load_filter().

cpl_frameset* muse_frameset_find_tags ( const cpl_frameset *  aFrames,
const cpl_array *  aTags,
unsigned char  aIFU,
cpl_boolean  aInvert 
)

return frameset containing data from an IFU/channel with the given tag(s)

Parameters
aFramesthe input frameset
aTagsthe required DFS tag(s) of the frame
aIFUthe IFU/channel number to check for (optional)
aInvertif true, invert the selection
Returns
the resulting frameset or NULL on failure

Returns the data in form of a frameset. Note that this frameset may be empty (but non-NULL), if no frame matching the parameters was found! If aIFU is 0, the IFU number is not checked.

This function just loops over muse_frameset_find() for all tags in aTags.

Exceptions
set CPL_ERROR_NULL_INPUT, return NULLaFrames is NULL

Definition at line 253 of file muse_utils.c.

References muse_frameset_find().

Referenced by muse_basicproc_combine_images_lampwise(), muse_basicproc_load_reduced(), and muse_frameset_check_raw().

cpl_frameset* muse_frameset_sort_raw_other ( const cpl_frameset *  aFrames,
int  aIndex,
const char *  aDateObs,
cpl_boolean  aSequence 
)

Create a new frameset containing all relevant raw frames first then all other frames.

Parameters
aFramesthe frames list
aIndexthe frame index (use a negative value to signify to take all raw frames or use aDateObs for comparison)
aDateObsthe DATE-OBS string to search for (only when aIndex is negative, use NULL to signify to take all frames)
aSequenceif target frame is part of a sequence
Returns
The new ordered frameset or NULL on error.

Normal raw frames are sorted first, if their location in aFrames match the aIndex, or aDateObs is given and the DATE-OBS in their header match it. If aSequence is true, then it is expected that they are part of an exposure sequence (e.g. calibration flats), and they are taken irrespective of aIndex.

Raw frames tagged "ILLUM" also get a special handling, in that the first of them – i.e. the one that got used in muse_basicproc_get_illum() – is sorted into the output frameset after all other raw frames, but before all other frames.

The returned frameset is duplicated from the input frameset and therefore has to be deallocated using cpl_frameset_delete() after use.

Exceptions
set CPL_ERROR_NULL_INPUT, return NULLinput frameset was NULL

Definition at line 418 of file muse_utils.c.

References muse_pfits_get_dateobs().

Referenced by muse_processing_append_used().

const char* muse_get_license ( void  )

Get the pipeline copyright and license.

Returns
The copyright and license string

The function wraps the cpl_get_license() macro and hence returns a pointer to the statically allocated license string. This string should not be modified using the returned pointer.

Definition at line 83 of file muse_utils.c.

cpl_matrix* muse_matrix_new_gaussian_2d ( int  aXHalfwidth,
int  aYHalfwidth,
double  aSigma 
)

Create a matrix that contains a normalized 2D Gaussian.

Parameters
aXHalfwidthhorizontal half width of the kernel matrix
aYHalfwidthvertical half width of the kernel matrix
aSigmasigma of Gaussian function
Returns
The new matrix or NULL on error

The 2*halfwidth+1 gives the size of one side of the matrix, i.e. halfwidth=3 for a 7x7 matrix. The created matrix has to be deallocated by cpl_matrix_delete().

Definition at line 1099 of file muse_utils.c.

Referenced by muse_geo_measure_spots().

muse_table* muse_table_load_filter ( muse_processing aProcessing,
const char *  aFilterName 
)

Load a table for a given filter name.

Parameters
aProcessingthe processing structure, includes the input frames list
aFilterNamethe filter name
Returns
the cpl_table * or NULL on error

This function tries to find a filter table and in it an extension for the filter of the given name. The table from that extension is loaded.

A filter called "white" is built into this function, it covers the full MUSE wavelength range with constant throughput.

Exceptions
set CPL_ERROR_NULL_INPUT, return NULLaFilterName is NULL
return NULL (and output debug message)aFilterName is "none"
set CPL_ERROR_FILE_NOT_FOUND, return NULLinput frame with tag MUSE_TAG_FILTER_LIST was not found
set CPL_ERROR_DATA_NOT_FOUND, return NULLvalid file with tag MUSE_TAG_FILTER_LIST does not contain filter named aFilterName
propagate error code, return NULLloading table for filter aFilterName from valid file with tag MUSE_TAG_FILTER_LIST failed

Definition at line 799 of file muse_utils.c.

References muse_table::header, muse_processing::inframes, muse_cpltable_new(), muse_frameset_find_master(), muse_processing_append_used(), muse_table_delete(), muse_table_new(), and muse_table::table.

Referenced by muse_postproc_cube_resample_and_collapse(), muse_postproc_process_exposure(), and muse_wcs_locate_sources().

cpl_error_code muse_utils_copy_modified_header ( cpl_propertylist *  aH1,
cpl_propertylist *  aH2,
const char *  aKeyword,
const char *  aString 
)

Copy a modified header keyword from one header to another.

Parameters
aH1input header
aH2output header
aKeywordheader keyword to copy
aStringstring to be used for the modification
Returns
CPL_ERROR_NONE on success or another CPL error on failure.

The keyword is modified by appending " (aString)" to it.

Exceptions
return CPL_ERROR_NULL_INPUTone of the input arguments is NULL
return CPL_ERROR_ILLEGAL_INPUTthe string in the input header cannot be read

Definition at line 2570 of file muse_utils.c.

Referenced by muse_datacube_save(), muse_datacube_save_recimages(), and muse_postproc_cube_resample_and_collapse().

cpl_error_code muse_utils_filter_copy_properties ( cpl_propertylist *  aHeader,
const muse_table aFilter,
double  aFraction 
)

Add/propagate filter properties to header of collapsed image.

Parameters
aHeaderthe header to add to
aFilterthe filter response curve
aFractionthe filter covering fraction
Returns
CPL_ERROR_NONE or a CPL error code on failure.
Exceptions
set and return CPL_ERROR_NULL_INPUTaHeader, aFilter, or the header component of aFilter are NULL

Definition at line 934 of file muse_utils.c.

References muse_table::header.

Referenced by muse_datacube_collapse(), muse_euro3dcube_collapse(), and muse_resampling_collapse_pixgrid().

double muse_utils_filter_fraction ( const muse_table aFilter,
double  aLambda1,
double  aLambda2 
)

Compute fraction of filter area covered by the data range.

Parameters
aFilterthe filter response curve
aLambda1the starting (data) wavelength
aLambda2the end (data) wavelength
Returns
the fraction or -1. on error.
Exceptions
set CPL_ERROR_NULL_INPUT, return -1.aFilter or its table component are NULL

Definition at line 893 of file muse_utils.c.

References MUSE_FLUX_RESP_FILTER, muse_flux_response_interpolate(), and muse_table::table.

Referenced by muse_resampling_collapse_pixgrid().

cpl_error_code muse_utils_fit_moffat_2d ( const cpl_matrix *  aPositions,
const cpl_vector *  aValues,
const cpl_vector *  aErrors,
cpl_array *  aParams,
cpl_array *  aPErrors,
const cpl_array *  aPFlags,
double *  aRMS,
double *  aRedChisq 
)

Fit a 2D Moffat function to a given set of data.

Parameters
aPositionsinput image to use for the fit
aValuesinput values to use for the fit
aErrorsinput errors to use for the fit (can be NULL)
aParamsparameter array (type double)
aPErrorsparameter error double array (can be NULL)
aPFlagsparameter flags integer array (can be NULL)
aRMSroot mean square of the fit (can be NULL)
aRedChisqreduced chi^2 of the fit (can be NULL)
Returns
CPL_ERROR_NONE on success any other CPL error code on failure.

This fits a Moffat (1969 A&A 3, 455, Eq. 7) function in two-dimensional cartesian form allowing for non-circular isophotes

  moffat(x,y) = B + A (beta - 1) / (pi alphax alphay sqrt(1 - rho^2))
                [ 1 + { ((x - xc) / alphax)^2
                       + 2 rho (x-xc)(y-yc) / (alphax alphay)
                       + ((y - yc) / alphay)^2} / {1 - rho^2} ]^(-beta)

to the input data. The parameters and their starting values are

  • aParams[0]: B, constant background level; median of the dataset
  • aParams[1]: A, volume of the function (flux); integrated value of the dataset in excess of the median
  • aParams[2]: xc, horizontal center; marginal centroid of the dataset
  • aParams[3]: yc, vertical center; see xc
  • aParams[4]: alphax, related to horizontal HWHM; computed from the HWHM of the dataset, using at least three data values near the half peak around the (marginal) centroid
  • aParams[5]: alphay, related to vertical HWHM; see alphax
  • aParams[6]: beta, seeing-related Moffat parameter; 2.5
  • aParams[7]: rho, the x-y correlation parameter; 0.0

Valid inputs in aParams override the first-guess parameters given above. If only one valid alpha parameter was given, the other one is used as starting value.

In case of uniform data, no error will be returned but a flat Moffat function with undefined center and widths will be created.

References:

Exceptions
return CPL_ERROR_NULL_INPUTaPositions, aValues, or aParams are NULL
return CPL_ERROR_INCOMPATIBLE_INPUTaPositions does not give the same number of positions as aValues
return CPL_ERROR_INCOMPATIBLE_INPUTaPositions does not contain two columns
return CPL_ERROR_INCOMPATIBLE_INPUTaErrors is given but does not have the same size as aValues
return CPL_ERROR_INVALID_TYPEaParams is not of type double
return CPL_ERROR_INVALID_TYPEaPErrors is given but not of type double
return CPL_ERROR_INVALID_TYPEaPFlags is given but not of type int
return CPL_ERROR_DATA_NOT_FOUNDaPErrors and/or aRedChisq are given but aErrors is NULL
return CPL_ERROR_SINGULAR_MATRIXinput data has less than 8 positions
return CPL_ERROR_ILLEGAL_INPUTaPFlags did not leave any parameter free
return CPL_ERROR_ILLEGAL_INPUTaPFlags was set to freeze a parameter but this parameter does not have an input value in aParams
propagate return code of cpl_fit_lvmq()the actual fit fails

Definition at line 1878 of file muse_utils.c.

References muse_cpl_optimize_lvmq(), MUSE_UTILS_CENTROID_MEDIAN, and muse_utils_get_centroid().

Referenced by muse_dar_check(), muse_flux_response_interpolate(), and muse_wcs_centroid_stars().

cpl_error_code muse_utils_fit_multigauss_1d ( const cpl_vector *  aX,
const cpl_bivector *  aY,
cpl_vector *  aCenter,
double *  aSigma,
cpl_vector *  aFlux,
cpl_vector *  aPoly,
double *  aMSE,
double *  aRedChisq,
cpl_matrix **  aCovariance 
)

Carry out a multi-Gaussian fit of data in a vector.

Parameters
aXPositions to fit
aYValues and associated errors to fit
aCentercenters of the Gaussian peaks
aSigmathe common Gaussian sigma of all peaks; use a negative value to have sigma be a fixed parameter in the fit
aFluxfluxes of the Gaussian peaks
aPolycoefficients of the background polynomial
aMSEcomputed mean squared error of the fit
aRedChisqreduced chi square of the fit
aCovarianceoutput covariance matrix of the fit
Returns
CPL_ERROR_NONE on success, another CPL error code on failure.

This function fits multiple Gaussians plus a background polynomial to the dataset given by aX and aY. The number of Gaussian peaks is determined by the length of the vector aCenter. It is mandatory to use it to give first-guess values for the centers of the peaks and an estimate of the expected Gaussian sigma in aSigma. aFlux has to be of the same length as aCenter, if given its values are also taken as first guesses of the areas of the Gaussian peaks (the fluxes). The background polynomial is defined using the aPoly vector, its contents are first-guess values for the polynomial coefficients, its length is used to define the order of the polynomial. If aPoly is NULL, the background is assumed to be a constant zero level.

If successful, fit results are returned in aCenter, aSigma, and – if given – aFlux. If aSigma was negative on input, the absolute value returned is same as the input, but positive. Then, the covariance matrix aCovariance is created and filled; it has to be deallocated using cpl_matrix_delete() after use. The order of parameters in the covariance matrix is the following: polynomial coefficients, sigma, centers, fluxes (where the number of elements, except for sigma, is determined by the lengths of the input objects).

This function is loosely modeled after cpl_vector_fit_gaussian(), it also uses cpl_fit_lvmq() to do that actual fit. But since computing first-guess parameters from the data is not possible in a meaningful way for multi-Gaussians, at least not without knowing the number of peaks, this function requires first guesses to be given, with the exception of the fluxes.

Exceptions
return CPL_ERROR_NULL_INPUTaX, aY, aCenter, or aSigma are NULL
return CPL_ERROR_INCOMPATIBLE_INPUTaY contains a different number of points as aX
return CPL_ERROR_INCOMPATIBLE_INPUTaFlux is given but contains a different number of parameters (peaks) than aCenter
return CPL_ERROR_ILLEGAL_INPUTaRedChisq is not NULL but there are less points than parameters, so that it cannot be computed
propagate errorcpl_fit_lvmq() fails

Definition at line 1558 of file muse_utils.c.

Referenced by muse_wave_line_fit_multiple().

cpl_error_code muse_utils_frameset_merge_frames ( cpl_frameset *  aFrames,
cpl_boolean  aDelete 
)

Merge IFU-specific files from a frameset to create multi-IFU outputs.

Parameters
aFramesthe frames list
aDeleteif CPL_TRUE, delete input files from framesete after merging
Returns
CPL_ERROR_NONE on success another CPL error code on failure

This function groups all frames in the input frameset regarding the base of the filenames (without -nn.fits), and then merges them, taking care to put IFU-specific information into the extension header (like QC and DRS parameters). The output filename for each group is the common basename of files (usually the PRO.CATG tag given to it by the pipeline during processing).

After all products were merged, the related files in the input frameset are physically removed from the filesystem and purged from aFrames.

The merged files in the frameset aFrames are set to the cpl_frame_group of CPL_FRAME_GROUP_PRODUCT, even if the input frames were not.

Note
Files tagged with "PIXTABLE_*" are not merged, but silently ignored!
Exceptions
set and return CPL_ERROR_NULL_INPUTaFrames is NULL
set and return CPL_ERROR_ILLEGAL_INPUTaFrames is empty

Definition at line 589 of file muse_utils.c.

cpl_error_code muse_utils_get_centroid ( const cpl_matrix *  aPositions,
const cpl_vector *  aValues,
const cpl_vector *  aErrors,
double *  aX,
double *  aY,
muse_utils_centroid_type  aBgType 
)

Compute centroid of a two-dimensional dataset.

Parameters
aPositionsmatrix containing the positions
aValuesvector containing the data values
aErrorsvector containing the data errors (1 sigma, can be NULL)
aXcomputed horizontal centroid position
aYcomputed vertical centroid position
aBgTypespecifies how to handle the background
Returns
CPL_ERROR_NONE on success, another cpl_error_code on failure.

This computes the center of gravity, using the data values and optionally the data errors as weights.

The input aPositions matrix must be for two-dimensional data, i.e. contain two columns, and the number of rows must be equal to the number of entries in both vectors.

Exceptions
return CPL_ERROR_NULL_INPUTaPositions or aValues are NULL
return CPL_ERROR_INCOMPATIBLE_INPUTaPositions does not contain two columns
return CPL_ERROR_INCOMPATIBLE_INPUTaPositions does not give the same number of positions as aValues
return CPL_ERROR_INCOMPATIBLE_INPUTaValues does not contain the same number of positions as aErrors
return CPL_ERROR_NULL_INPUTboth aX and aY are NULL so that the centroid cannot be returned
return CPL_ERROR_ILLEGAL_INPUTaBgType does not contain a valid centroid type

Definition at line 1328 of file muse_utils.c.

References MUSE_UTILS_CENTROID_MEAN, MUSE_UTILS_CENTROID_MEDIAN, and MUSE_UTILS_CENTROID_NORMAL.

Referenced by muse_dar_check(), and muse_utils_fit_moffat_2d().

int muse_utils_get_extension_for_ifu ( const char *  aFilename,
unsigned char  aIFU 
)

Return extension number that corresponds to this IFU/channel number.

Parameters
aFilenamename of the FITS file
aIFUthe IFU/channel number
Returns
Corresponding extension number, or -1 if not found.

Definition at line 118 of file muse_utils.c.

References muse_pfits_has_ifu().

Referenced by muse_basicproc_params_delete(), muse_frameset_check_raw(), muse_frameset_find(), and muse_table_load().

unsigned char muse_utils_get_ifu ( const cpl_propertylist *  aHeaders)
char* muse_utils_header_get_lamp_names ( cpl_propertylist *  aHeader,
char  aSep 
)

Concatenate names of all active calibration lamps.

Parameters
aHeaderthe FITS header to search for lamp entries
aSepthe separator to use for name concatenation
Returns
A newly allocated string with the filter name(s) or NULL on error or if no lamps are switched on.

This function iterates through all lamp shutters found in the header, and for those where the lamp is switched on and the shutter is open appends the name to the output string. Prefixes of lamp names (currently only "CU-LAMP-") are cut off, e.g. "CU-LAMP-HgCd" becomes "HgCd". XXX to aid AIT, lamp names like "CU-LAMP[3-6]" are converted to the respective element names (Ne, Xe, HgCd).

The returned string has to be deallocated with cpl_free().

Exceptions
set CPL_ERROR_NULL_INPUT, return NULLaHeader is NULL

Definition at line 989 of file muse_utils.c.

References muse_pfits_get_lamp_name(), muse_pfits_get_lamp_status(), muse_pfits_get_lampnum(), and muse_pfits_get_shut_status().

Referenced by muse_lsf_create_arcpixtable(), and muse_wave_calib_lampwise().

cpl_array* muse_utils_header_get_lamp_numbers ( cpl_propertylist *  aHeader)

List numbers of all active calibration lamps.

Parameters
aHeaderthe FITS header to search for lamp entries
Returns
A cpl_array of type CPL_TYPE_INT, or NULL on error or if no lamps are switched on.

This function iterates through all lamp shutters found in the header, and for those where the lamp is switched on and the shutter is open appends the lamp number to the output array.

The lamp numbers are the i in the ESO.INS.LAMPi keywords in the MUSE FITS headers.

The returned array has to be deallocated with cpl_array_delete().

Exceptions
set CPL_ERROR_NULL_INPUT, return NULLaHeader is NULL

Definition at line 1056 of file muse_utils.c.

References muse_pfits_get_lamp_status(), muse_pfits_get_lampnum(), and muse_pfits_get_shut_status().

Referenced by muse_wave_calib_lampwise().

cpl_image* muse_utils_image_fit_polynomial ( const cpl_image *  aImage,
unsigned short  aXOrder,
unsigned short  aYOrder 
)

Create a smooth version of a 2D image by fitting it with a 2D polynomial.

Parameters
aImageinput image to fit
aXOrderpolynomial order to use in x direction
aYOrderpolynomial order to use in y direction
Returns
a new cpl_image * constructed from the polynomial fit or NULL on error

This function transfers all pixels (their positions and values) into suitable structures and calls cpl_polynomial_fit(). If that polynomial fit succeeds, cpl_image_fill_polynomial() is used to create the output image. Rejected (bad) pixels in the input image are marked as rejected in the output image as well.

Exceptions
set CPL_ERROR_NULL_INPUT, return NULLaImage is NULL
set CPL_ERROR_DATA_NOT_FOUND, return NULLaImage does not contain good pixels
propagate error code, return NULLcpl_polynomial_fit() fails

Definition at line 1147 of file muse_utils.c.

cpl_error_code muse_utils_image_get_centroid_window ( cpl_image *  aImage,
int  aX1,
int  aY1,
int  aX2,
int  aY2,
double *  aX,
double *  aY,
muse_utils_centroid_type  aBgType 
)

Compute centroid over an image window, optionally marginalizing over the background.

Parameters
aImageinput image to use
aX1lower left x position of window
aY1lower left y position of window
aX2upper right x position of window
aY2upper right y position of window
aXcomputed centroid in x-direction
aYcomputed centroid in x-direction
aBgTypespecifies how to handle the background
Returns
CPL_ERROR_NONE on success, another CPL error code on failure

This computes the ("marginal") centroid over an image window, similar to IRAF's imcentroid algorithm but without iterations. Before the normal centroid is actually computed, the average level of the image is subtracted. Only the positions with positive values are then used in the computation of the output centroid.

Exceptions
return CPL_ERROR_NULL_INPUTaImage or both aXCen and aYCen are NULL
propagate error codecpl_image_extract() fails for the given coordinates
return CPL_ERROR_ILLEGAL_INPUTaBgType does not contain a valid centroid type

Definition at line 1234 of file muse_utils.c.

References MUSE_UTILS_CENTROID_MEAN, MUSE_UTILS_CENTROID_MEDIAN, and MUSE_UTILS_CENTROID_NORMAL.

Referenced by muse_dar_check(), and muse_wcs_centroid_stars().

cpl_polynomial* muse_utils_iterate_fit_polynomial ( cpl_matrix *  aPos,
cpl_vector *  aVal,
cpl_vector *  aErr,
cpl_table *  aExtra,
const unsigned int  aOrder,
const double  aRSigma,
double *  aMSE,
double *  aChiSq 
)

Iterate a polynomial fit.

Parameters
aPosmatrix with input positions
aValvector with input values
aErrvector with input errors (optional, currently unused)
aExtratable with extra information (optional, same number of rows as aVal)
aOrderpolynomial order to use for the fit
aRSigmarejection sigma
aMSEoutput mean squared error of final fit (optional)
aChiSqoutput chi square of final fit (optional)
Returns
CPL_ERROR_NONE on success any other CPL error code on failure.

The position matrix aPos contains the same number of columns as the values vector aVal. The number of rows in aPos defines the number of dimensions of the output polynomial.

The table aExtra contains extra information, with each matrix row related to each entry in the aVal vector.

Infinite entries (NANs and INFs) in aVal are rejected at the start, so can be used as "bad pixel" markers. The procedure will fail, if only infinite values are present.

Note
All input structures (i.e. aPos, aVal, aErr, and aExtra) are changed by this routine, so that entries rejected during the iterations are erased. This means that pointers to the structures may be changed afterwards.

If an error occurs, both aChiSq and aMSE are set to high values (DBL_MAX).

Exceptions
set CPL_ERROR_NULL_INPUT, return NULLaPos and/or aVal are NULL
set CPL_ERROR_INCOMPATIBLE_INPUT, return NULLthe number of aPos does not match the size of aVal
set CPL_ERROR_INCOMPATIBLE_INPUT, return NULLaErr was specified but does not match the size of aVal
set CPL_ERROR_INCOMPATIBLE_INPUT, return NULLaExtra was specified but the number of its rows does not match the size of aVal
set CPL_ERROR_ILLEGAL_INPUT, return NULLaVal only contains infinite values
set CPL_ERROR_ILLEGAL_OUTPUT, return NULLthe iterative procedure tries to remove the last vector/matrix element

Definition at line 2234 of file muse_utils.c.

References muse_cplvector_erase_element().

Referenced by muse_dar_check(), muse_flux_response_interpolate(), muse_geo_finalize(), muse_quadrants_overscan_polyfit_vertical(), muse_quadrants_overscan_stats(), muse_trace_iterate_fit(), muse_wave_calib_lampwise(), and muse_wave_line_fit_iterate().

void muse_utils_memory_dump ( const char *  aMarker)

Display the current memory usage of the given program.

Parameters
aMarkermarker string to use in the debug output
Note
This function prints CPL memory info to stderr and other data to stdout!
This function only outputs anything if an executable name was given using the environment variable MUSE_DEBUG_MEMORY_PROGRAM.

This function uses the cpl_memory_dump() function to list pointers allocated through CPL and the Unix ps command to display system information about the process.

Definition at line 2702 of file muse_utils.c.

Referenced by muse_resampling_cube(), and muse_xcombine_tables().

double muse_utils_pixtable_fit_line_gaussian ( muse_pixtable aPixtable,
double  aLambda,
double  aHalfWidth,
double  aBinSize,
float  aLo,
float  aHi,
unsigned char  aIter,
cpl_array *  aResults,
cpl_array *  aErrors 
)

Fit a 1D Gaussian to a given wavelength range in a pixel table.

Parameters
aPixtableinput pixel table
aLambdathe wavelength around which to extract the range
aHalfWidththe half width of the range in wavelength to use
aBinSizethe size of the bins in wavelength for the spectrum
aLolow sigma-clipping limit for the spectrum
aHihigh sigma-clipping limit for the spectrum
aIternumber of iterations for sigma-clipping the spectrum
aResultsparameter error double array (can be NULL)
aErrorsparameter flags integer array (can be NULL)
Returns
the center of the Gaussian on success or 0. on error.

This fits a 1D Gaussian of the form

  gauss(x) = A / sqrt(2 pi sigma^2) * exp( -(x - xc)^2/(2 sigma^2)) + B

to the input pixel table. It does this by resampling the pixel table in the specified wavelength range into a temporary spectrum, sampled at aBinSize, using muse_resampling_spectrum(), which is then given to cpl_vector_fit_gaussian().

This function directly returns the center xc of the Gaussian fit. If the other properties are required, the optional inputs aResults (for the values) and aErrors (for the sigma values) need to be passed as CPL arrays of type CPL_TYPE_DOUBLE. These arrays are resized to 4 elements and overwritten with output. The Gaussian parameters are always written to aResults (if given), the errors returned in aErrors may be invalid, except for the 0th element (the centroid error) which is always written. The order of the elements is:

  • aParams[0]: xc, the Gaussian center
  • aParams[1]: sigma
  • aParams[2]: A, the area under the Gaussian curve
  • aParams[3]: B, the background level

First-guess parameters cannot be given, but the wavelength range should be given to restrict the data to a set of data that contains one Gaussian-like peak.

Note
The input pixel table is preserved, except for the table selection flags, which are reset to none selected when the function returns.
Exceptions
set CPL_ERROR_NULL_INPUT, return 0.aPixtable or one of its components is NULL
set CPL_ERROR_DATA_NOT_FOUND, return 0.aPixtable does not have data in the range |aLambda| +/- aHalfWidth
propagate error code, return 0.muse_resampling_spectrum() fails

Definition at line 2442 of file muse_utils.c.

References muse_pixtable::header, MUSE_PIXTABLE_LAMBDA, muse_resampling_spectrum_iterate(), muse_utils_get_ifu(), and muse_pixtable::table.

Referenced by muse_basicproc_shift_pixtable().

cpl_error_code muse_utils_set_hduclass ( cpl_propertylist *  aHeader,
const char *  aClass2,
const char *  aExtData,
const char *  aExtDQ,
const char *  aExtStat 
)

Set HDU headers for the ESO FITS data format.

Parameters
aHeaderthe FITS headers to modify
aClass2type of class to set, i.e. contents of HDUCLAS2
aExtDataextension name for the data extension
aExtDQextension name for the bad pixel extension (or NULL)
aExtStatextension name for the variance extension (or NULL)
Returns
CPL_ERROR_NONE on success, another CPL error code on failure.

This function adds special FITS headers to support the ESO format. The information was taken from v2.5.1 (dated Feb. 15th, 2012) of the document "FITS format description for pipeline products with data, error and data quality information" and further info by Martin Kuemmel and Harald Kuntschner.

This function tries to set these header keywords directly after EXTNAME, if available, otherwise at the end of the input header.

Quality Assessment:
This function is used by muse_image_save() and muse_datacube_save(). Tests of their functionality include tests of the headers that are created using this function.
Exceptions
return CPL_ERROR_NULL_INPUTaHeader, aClass2, or aExtData are NULL
return CPL_ERROR_ILLEGAL_INPUTaClass2 is neither "DATA", "ERROR", nor "QUALITY"

Definition at line 2611 of file muse_utils.c.

Referenced by muse_datacube_save(), muse_datacube_save_recimages(), and muse_image_save().

Variable Documentation

const muse_cpltable_def muse_filtertable_def[]
Initial value:
= {
{ "lambda", CPL_TYPE_DOUBLE, "Angstrom", "%7.2f", "wavelength", CPL_TRUE},
{ "throughput", CPL_TYPE_DOUBLE, "", "%.4e",
"filter response (in fractions of 1)", CPL_TRUE},
{ NULL, 0, NULL, NULL, NULL, CPL_FALSE }
}

MUSE filter table definition.

A MUSE filter table has the following columns:

  • 'lambda': the wavelength in Angstrom
  • 'throughput': the filter response in fractions of 1

Definition at line 768 of file muse_utils.c.