Vector


Functions

cpl_error_code cpl_vector_add (cpl_vector *v1, const cpl_vector *v2)
 Add a cpl_vector to another.
cpl_error_code cpl_vector_add_scalar (cpl_vector *v, double addend)
 Elementwise addition of a scalar to a vector.
cpl_error_code cpl_vector_copy (cpl_vector *destination, const cpl_vector *source)
 This function copies contents of a vector into another vector.
int cpl_vector_correlate (cpl_vector *vxc, const cpl_vector *v1, const cpl_vector *v2)
 Cross-correlation of two vectors.
void cpl_vector_delete (cpl_vector *v)
 Delete a cpl_vector.
cpl_error_code cpl_vector_divide (cpl_vector *v1, const cpl_vector *v2)
 Divide two vectors element-wise.
cpl_error_code cpl_vector_divide_scalar (cpl_vector *v, double divisor)
 Elementwise division of a vector with a scalar.
void cpl_vector_dump (const cpl_vector *v, FILE *stream)
 Dump a cpl_vector as ASCII to a stream.
cpl_vector * cpl_vector_duplicate (const cpl_vector *v)
 This function duplicates an existing vector and allocates memory.
cpl_error_code cpl_vector_exponential (cpl_vector *v, double base)
 Compute the exponential of all vector elements.
cpl_vector * cpl_vector_extract (const cpl_vector *v, int istart, int istop, int istep)
 Extract a sub_vector from a vector.
cpl_error_code cpl_vector_fill (cpl_vector *v, double val)
 Fill a cpl_vector.
cpl_error_code cpl_vector_fill_kernel_profile (cpl_vector *profile, cpl_kernel type, double radius)
 Fill a vector with a kernel profile.
cpl_vector * cpl_vector_filter_lowpass_create (const cpl_vector *v, cpl_lowpass filter_type, int hw)
 Apply a low-pass filter to a cpl_vector.
cpl_vector * cpl_vector_filter_median_create (const cpl_vector *v, int hw)
 Apply a 1d median filter of given half-width to a cpl_vector.
int cpl_vector_find (const cpl_vector *sorted, double key)
 In a sorted vector find the element closest to the given value.
cpl_error_code cpl_vector_fit_gaussian (const cpl_vector *x, const cpl_vector *sigma_x, const cpl_vector *y, const cpl_vector *sigma_y, cpl_fit_mode fit_pars, double *x0, double *sigma, double *area, double *offset, double *mse, double *red_chisq, cpl_matrix **covariance)
 Apply a 1d gaussian fit.
double cpl_vector_get (const cpl_vector *in, int idx)
 Get an element of the vector.
double * cpl_vector_get_data (cpl_vector *in)
 Get a pointer to the data part of the vector.
const double * cpl_vector_get_data_const (const cpl_vector *in)
 Get a pointer to the data part of the vector.
double cpl_vector_get_max (const cpl_vector *v)
 Get the maximum of the cpl_vector.
double cpl_vector_get_mean (const cpl_vector *v)
 Compute the mean value of vector elements.
double cpl_vector_get_median (cpl_vector *v)
 Compute the median of the elements of a vector.
double cpl_vector_get_median_const (const cpl_vector *v)
 Compute the median of the elements of a vector.
double cpl_vector_get_min (const cpl_vector *v)
 Get the minimum of the cpl_vector.
int cpl_vector_get_size (const cpl_vector *in)
 Get the size of the vector.
double cpl_vector_get_stdev (const cpl_vector *v)
 Compute the bias-corrected standard deviation of a vectors elements.
cpl_vector * cpl_vector_load (const char *filename, int xtnum)
 Load a list of values from a FITS file.
cpl_error_code cpl_vector_logarithm (cpl_vector *v, double base)
 Compute the element-wise logarithm.
cpl_error_code cpl_vector_multiply (cpl_vector *v1, const cpl_vector *v2)
 Multiply two vectors component-wise.
cpl_error_code cpl_vector_multiply_scalar (cpl_vector *v, double factor)
 Elementwise multiplication of a vector with a scalar.
cpl_vector * cpl_vector_new (int n)
 Create a new cpl_vector.
cpl_error_code cpl_vector_power (cpl_vector *v, double exponent)
 Compute the power of all vector elements.
double cpl_vector_product (const cpl_vector *v1, const cpl_vector *v2)
 Compute the vector dot product.
cpl_vector * cpl_vector_read (const char *filename)
 Read a list of values from an ASCII file and create a cpl_vector.
cpl_error_code cpl_vector_save (const cpl_vector *to_save, const char *filename, cpl_type_bpp bpp, const cpl_propertylist *pl, unsigned mode)
 Save a vector to a FITS file.
cpl_error_code cpl_vector_set (cpl_vector *in, int idx, double value)
 Set an element of the vector.
cpl_error_code cpl_vector_set_size (cpl_vector *in, int newsize)
 Resize the vector.
cpl_error_code cpl_vector_sort (cpl_vector *v, int c)
 Sort a cpl_vector by increasing/decreasing data.
cpl_error_code cpl_vector_sqrt (cpl_vector *v)
 Compute the sqrt of a cpl_vector.
cpl_error_code cpl_vector_subtract (cpl_vector *v1, const cpl_vector *v2)
 Subtract a cpl_vector from another.
cpl_error_code cpl_vector_subtract_scalar (cpl_vector *v, double subtrahend)
 Elementwise subtraction of a scalar from a vector.
void * cpl_vector_unwrap (cpl_vector *v)
 Delete a cpl_vector except the data array.
cpl_vector * cpl_vector_wrap (int n, double *data)
 Create a cpl_vector from existing data.

Detailed Description

This module provides functions to handle cpl_vector.

A cpl_vector is an object containing a list of values (as doubles) and the (always positive) number of values. The functionalities provided here are simple ones like sorting, statistics, or simple operations. The cpl_bivector object is composed of two of these vectors.

Synopsis:
   #include "cpl_vector.h"

Function Documentation

cpl_error_code cpl_vector_add ( cpl_vector *  v1,
const cpl_vector *  v2 
)

Add a cpl_vector to another.

Parameters:
v1 First cpl_vector (modified)
v2 Second cpl_vector
Returns:
CPL_ERROR_NONE or the relevant _cpl_error_code_ on error
The second vector is added to the first one. The input first vector is modified.

The input vectors must have the same size.

Possible _cpl_error_code_ set in this function:

cpl_error_code cpl_vector_add_scalar ( cpl_vector *  v,
double  addend 
)

Elementwise addition of a scalar to a vector.

Parameters:
v cpl_vector to modify
addend Number to add
Returns:
CPL_ERROR_NONE or the relevant _cpl_error_code_ on error
Add a number to each element of the cpl_vector.

Possible _cpl_error_code_ set in this function:

cpl_error_code cpl_vector_copy ( cpl_vector *  destination,
const cpl_vector *  source 
)

This function copies contents of a vector into another vector.

Parameters:
destination destination cpl_vector
source source cpl_vector
Returns:
CPL_ERROR_NONE or the relevant _cpl_error_code_ on error
See also:
cpl_vector_set_size() if source and destination have different sizes.
Possible _cpl_error_code_ set in this function:

int cpl_vector_correlate ( cpl_vector *  vxc,
const cpl_vector *  v1,
const cpl_vector *  v2 
)

Cross-correlation of two vectors.

Parameters:
vxc Odd-sized vector with the computed cross-correlations
v1 1st vector to correlate
v2 2nd vector to correlate
Returns:
Index of maximum cross-correlation, or negative on error
vxc must have an odd number of elements, 2*half_search+1, where half_search is the half-size of the search domain.

The length of v2 may not exceed that of v1. If the difference in length between v1 and v2 is less than half_search then this difference must be even (if the difference is odd resampling of v2 may be useful).

The cross-correlation is computed with shifts ranging from -half_search to half_search.

On succesful return element i (starting with 0) of vxc contains the cross- correlation at offset i-half_search. On error vxc is unmodified.

The cross-correlation is in fact the dot-product of two unit-vectors and ranges therefore from -1 to 1.

The cross-correlation is, in absence of rounding errors, commutative only for equal-sized vectors, i.e. changing the order of v1 and v2 will move element j in vxc to 2*half_search - j and thus change the return value from i to 2*half_search - i.

If, in absence of rounding errors, more than one shift would give the maximum cross-correlation, rounding errors may cause any one of those shifts to be returned. If rounding errors have no effect the index corresponding to the shift with the smallest absolute value is returned (with preference given to the smaller of two indices that correspond to the same absolute shift).

If v1 is longer than v2, the first element in v1 used for the resulting cross-correlation is max(0,shift + (cpl_vector_get_size(v1)-cpl_vector_get_size(v2))/2).

Cross-correlation with half_search == 0 requires about 8n FLOPs, where n = cpl_vector_get_size(v2). Each increase of half_search by 1 requires about 4n FLOPs more, when all of v2's elements can be cross-correlated, otherwise the extra cost is about 4m, where m is the number of elements in v2 that can be cross-correlated, n - half_search <= m < n.

In case of error, the _cpl_error_code_ code is set, and the returned delta and cross-correlation is undefined.

Example of 1D-wavelength calibration (error handling omitted for brevity):

  cpl_vector   * model = my_model(dispersion);
  cpl_vector   * vxc   = cpl_vector_new(1+2*maxshift);
  const int      shift = cpl_vector_correlate(vxc, model, observed) - maxshift;
  cpl_error_code error = cpl_polynomial_shift_1d(dispersion, 0, (double)shift);

  cpl_msg_info(cpl_func, "Shifted dispersion relation by %d pixels, shift);

Possible _cpl_error_code_ set in this function:

void cpl_vector_delete ( cpl_vector *  v  ) 

Delete a cpl_vector.

Parameters:
v cpl_vector to delete
Returns:
void
If the vector v is NULL, nothing is done and no error is set.

cpl_error_code cpl_vector_divide ( cpl_vector *  v1,
const cpl_vector *  v2 
)

Divide two vectors element-wise.

Parameters:
v1 First cpl_vector
v2 Second cpl_vector
Returns:
CPL_ERROR_NONE or the relevant _cpl_error_code_ on error
See also:
cpl_vector_add()
If an element in v2 is zero v1 is not modified and CPL_ERROR_DIVISION_BY_ZERO is returned.

Possible _cpl_error_code_ set in this function:

cpl_error_code cpl_vector_divide_scalar ( cpl_vector *  v,
double  divisor 
)

Elementwise division of a vector with a scalar.

Parameters:
v cpl_vector to modify
divisor Non-zero number to divide with
Returns:
CPL_ERROR_NONE or the relevant _cpl_error_code_ on error
Divide each element of the cpl_vector with a number.

Possible _cpl_error_code_ set in this function:

void cpl_vector_dump ( const cpl_vector *  v,
FILE *  stream 
)

Dump a cpl_vector as ASCII to a stream.

Parameters:
v Input cpl_vector to dump
stream Output stream, accepts stdout or stderr
Returns:
void
Each element is preceded by its index number (starting with 1!) and written on a single line.

Comment lines start with the hash character.

stream may be NULL in which case stdout is used.

Note:
In principle a cpl_vector can be saved using cpl_vector_dump() and re-read using cpl_vector_read(). This will however introduce significant precision loss due to the limited accuracy of the ASCII representation.

cpl_vector* cpl_vector_duplicate ( const cpl_vector *  v  ) 

This function duplicates an existing vector and allocates memory.

Parameters:
v the input cpl_vector
Returns:
a newly allocated cpl_vector or NULL in case of an error
The returned object must be deallocated using cpl_vector_delete()

Possible _cpl_error_code_ set in this function:

cpl_error_code cpl_vector_exponential ( cpl_vector *  v,
double  base 
)

Compute the exponential of all vector elements.

Parameters:
v Target cpl_vector.
base Exponential base.
Returns:
CPL_ERROR_NONE or the relevant _cpl_error_code_ on error
If the base is zero all vector elements must be positive and if the base is negative all vector elements must be integer, otherwise a cpl_error_code is returned and the vector is unmodified.

Possible _cpl_error_code_ set in this function:

cpl_vector* cpl_vector_extract ( const cpl_vector *  v,
int  istart,
int  istop,
int  istep 
)

Extract a sub_vector from a vector.

Parameters:
v Input cpl_vector
istart Start index (from 0 to number of elements - 1)
istop Stop index (from 0 to number of elements - 1)
istep Extract every step element
Returns:
A newly allocated cpl_vector or NULL in case of an error
The returned object must be deallocated using cpl_vector_delete()

FIXME: Currently istop must be greater than istart. FIXME: Currently istep must equal 1.

Possible _cpl_error_code_ set in this function:

cpl_error_code cpl_vector_fill ( cpl_vector *  v,
double  val 
)

Fill a cpl_vector.

Parameters:
v cpl_vector to be filled with the value val
val Value used to fill the cpl_vector
Returns:
CPL_ERROR_NONE or the relevant _cpl_error_code_ on error
Input vector is modified

Possible _cpl_error_code_ set in this function:

cpl_error_code cpl_vector_fill_kernel_profile ( cpl_vector *  profile,
cpl_kernel  type,
double  radius 
)

Fill a vector with a kernel profile.

Parameters:
profile cpl_vector to be filled
type Type of kernel profile
radius Radius of the profile in pixels
Returns:
CPL_ERROR_NONE or the relevant _cpl_error_code_ on error
See also:
cpl_image_get_interpolated
A number of predefined kernel profiles are available:

Possible _cpl_error_code_ set in this function:

cpl_vector* cpl_vector_filter_lowpass_create ( const cpl_vector *  v,
cpl_lowpass  filter_type,
int  hw 
)

Apply a low-pass filter to a cpl_vector.

Parameters:
v cpl_vector
filter_type Type of filter to use
hw Filter half-width
Returns:
Pointer to newly allocated cpl_vector or NULL in case of an error
This type of low-pass filtering consists in a convolution with a given kernel. The chosen filter type determines the kind of kernel to apply for convolution. Supported kernels are CPL_LOWPASS_LINEAR and CPL_LOWPASS_GAUSSIAN.

In the case of CPL_LOWPASS_GAUSSIAN, the gaussian sigma used is 1/sqrt(2). As this function is not meant to be general and cover all possible cases, this sigma is hardcoded and cannot be changed.

The returned smooth cpl_vector must be deallocated using cpl_vector_delete(). The returned signal has exactly as many samples as the input signal.

Possible _cpl_error_code_ set in this function:

cpl_vector* cpl_vector_filter_median_create ( const cpl_vector *  v,
int  hw 
)

Apply a 1d median filter of given half-width to a cpl_vector.

Parameters:
v cpl_vector
hw Filter half-width
Returns:
Pointer to newly allocated cpl_vector or NULL in case of an error
This function applies a median smoothing to a cpl_vector and returns a newly allocated cpl_vector containing a median-smoothed version of the input. The returned cpl_vector has exactly as many samples as the input one. It must be deallocated using cpl_vector_delete(). For half-widths of 1,2,3,4, the filtering is optimised for speed.

Possible _cpl_error_code_ set in this function:

int cpl_vector_find ( const cpl_vector *  sorted,
double  key 
)

In a sorted vector find the element closest to the given value.

Parameters:
sorted Sorted cpl_vector
key Value to find
Returns:
The index that minimizes fabs(sorted[index] - key) or -1 on error
Bisection is used to find the element.

If two (neighboring) elements with different values both minimize fabs(sorted[index] - key) the index of the larger element is returned.

If the vector contains identical elements that minimize fabs(sorted[index] - key) then it is undefined which element has its index returned.

Possible _cpl_error_code_ set in this function:

cpl_error_code cpl_vector_fit_gaussian ( const cpl_vector *  x,
const cpl_vector *  sigma_x,
const cpl_vector *  y,
const cpl_vector *  sigma_y,
cpl_fit_mode  fit_pars,
double *  x0,
double *  sigma,
double *  area,
double *  offset,
double *  mse,
double *  red_chisq,
cpl_matrix **  covariance 
)

Apply a 1d gaussian fit.

Parameters:
x Positions to fit
sigma_x Uncertainty (one sigma, gaussian errors assumed) assosiated with x. Taking into account the uncertainty of the independent variable is currently unsupported, and this parameter must therefore be set to NULL.
y Values to fit
sigma_y Uncertainty (one sigma, gaussian errors assumed) associated with y. If NULL, constant uncertainties are assumed.
fit_pars Specifies which parameters participate in the fit (any other parameters will be held constant). Possible values are CPL_FIT_CENTROID, CPL_FIT_STDEV, CPL_FIT_AREA, CPL_FIT_OFFSET and any bitwise combination of these. As a shorthand for including all four parameters in the fit, use CPL_FIT_ALL.
x0 (output) Center of best fit gaussian. If CPL_FIT_CENTROID is not set, this is also an input parameter.
sigma (output) Width of best fit gaussian. A positive number on success. If CPL_FIT_STDEV is not set, this is also an input parameter.
area (output) Area of gaussian. A positive number on succes. If CPL_FIT_AREA is not set, this is also an input parameter.
offset (output) Fitted background level. If CPL_FIT_OFFSET is not set, this is also an input parameter.
mse (output) If non-NULL, the mean squared error of the best fit is returned.
red_chisq (output) If non-NULL, the reduced chi square of the best fit is returned. This requires the noise vector to be specified.
covariance (output) If non-NULL, the formal covariance matrix of the best fit is returned. This requires sigma_y to be specified. The order of fit parameters in the covariance matrix is defined as (x0, sigma, area, offset), for example the (3,3) element of the matrix (counting from zero) is the variance of the fitted offset. The matrix must be deallocated by calling cpl_matrix_delete() . On error, NULL is returned.
Returns:
CPL_ERROR_NONE iff okay
This function fits to the input vectors a 1d gaussian function of the form

f(x) = area / sqrt(2 pi sigma^2) * exp( -(x - x0)^2/(2 sigma^2)) + offset

(area > 0) by minimizing chi^2 using a Levenberg-Marquardt algorithm.

The values to fit are read from the input vector x. Optionally, a vector sigma_x (of same size as x) may be specified.

Optionally, the mean squared error, the reduced chi square and the covariance matrix of the best fit are computed . Set corresponding parameter to NULL to ignore.

If the covariance matrix is requested and successfully computed, the diagonal elements (the variances) are guaranteed to be positive.

Occasionally, the Levenberg-Marquardt algorithm fails to converge to a set of sensible parameters. In this case (and only in this case), a CPL_ERROR_CONTINUE is set. To allow the caller to recover from this particular error, the parameters x0, sigma, area and offset will on output contain estimates of the best fit parameters, specifically estimated as the median position, the median of the absolute residuals multiplied by 1.4828, the minimum flux value and the maximum flux difference multiplied by sqrt(2 pi sigma^2), respectively. In this case, mse, red_chisq and covariance are not computed. Note that the variance of x0 (the (0,0) element of the covariance matrix) in this case can be estimated by sigma^2 / area .

A CPL_ERROR_SINGULAR_MATRIX occurs if the covariance matrix cannot be computed. In that case all other output parameters are valid.

Current limitations

Possible _cpl_error_code_ set in this function:

double cpl_vector_get ( const cpl_vector *  in,
int  idx 
)

Get an element of the vector.

Parameters:
in the input vector
idx the index of the element (0 to nelem-1)
Returns:
The element value
In case of error, the _cpl_error_code_ code is set, and the returned double is undefined.

Possible _cpl_error_code_ set in this function:

double* cpl_vector_get_data ( cpl_vector *  in  ) 

Get a pointer to the data part of the vector.

Parameters:
in the input vector
Returns:
Pointer to the data or NULL in case of an error
The returned pointer refers to already allocated data.

Note:
Use at your own risk: direct manipulation of vector data rules out any check performed by the vector object interface, and may introduce inconsistencies between the information maintained internally, and the actual vector data and structure.
Possible _cpl_error_code_ set in this function:

const double* cpl_vector_get_data_const ( const cpl_vector *  in  ) 

Get a pointer to the data part of the vector.

Parameters:
in the input vector
Returns:
Pointer to the data or NULL in case of an error
See also:
cpl_vector_get_data

double cpl_vector_get_max ( const cpl_vector *  v  ) 

Get the maximum of the cpl_vector.

Parameters:
v const cpl_vector
Returns:
the maximum value of the vector or undefined on error
See also:
cpl_vector_get_min()

double cpl_vector_get_mean ( const cpl_vector *  v  ) 

Compute the mean value of vector elements.

Parameters:
v Input const cpl_vector
Returns:
Mean value of vector elements or undefined on error.
See also:
cpl_vector_get_min()

double cpl_vector_get_median ( cpl_vector *  v  ) 

Compute the median of the elements of a vector.

Parameters:
v Input cpl_vector
Returns:
Median value of the vector elements or undefined on error.
See also:
cpl_vector_get_median_const()
Note:
For efficiency reasons, this function modifies the order of the elements of the input vector.

double cpl_vector_get_median_const ( const cpl_vector *  v  ) 

Compute the median of the elements of a vector.

Parameters:
v Input const cpl_vector
Returns:
Median value of the vector elements or undefined on error.
See also:
cpl_vector_get_min()
For a finite population or sample, the median is the middle value of an odd number of values (arranged in ascending order) or any value between the two middle values of an even number of values. The computed median should actually be a value of the input array, so, in case of an even number of values in the input array, the median would be the lowest of the two central values.

double cpl_vector_get_min ( const cpl_vector *  v  ) 

Get the minimum of the cpl_vector.

Parameters:
v const cpl_vector
Returns:
The minimum value of the vector or undefined on error
In case of error, the _cpl_error_code_ code is set, and the returned double is undefined.

Possible _cpl_error_code_ set in this function:

int cpl_vector_get_size ( const cpl_vector *  in  ) 

Get the size of the vector.

Parameters:
in the input vector
Returns:
The size or -1 in case of an error
Possible _cpl_error_code_ set in this function:

double cpl_vector_get_stdev ( const cpl_vector *  v  ) 

Compute the bias-corrected standard deviation of a vectors elements.

Parameters:
v Input const cpl_vector
Returns:
standard deviation of the elements or a negative number on error.
See also:
cpl_vector_get_min()
S(n-1) = sqrt((1/n-1) sum((xi-mean)^2) (i=1 -> n)

The length of v must be at least 2.

cpl_vector* cpl_vector_load ( const char *  filename,
int  xtnum 
)

Load a list of values from a FITS file.

Parameters:
filename Name of the input file
xtnum Extension number in the file.
Returns:
1 newly allocated cpl_vector or NULL in case of an error
See also:
cpl_vector_save
This function loads a vector from a FITS file (NAXIS=1), using cfitsio. The returned image has to be deallocated with cpl_vector_delete().

'xtnum' specifies from which extension the image should be loaded. This could be 0 for the main data section (files without extension), or any number between 1 and N, where N is the number of extensions present in the file.

Possible _cpl_error_code_ set in this function:

cpl_error_code cpl_vector_logarithm ( cpl_vector *  v,
double  base 
)

Compute the element-wise logarithm.

Parameters:
v cpl_vector to modify.
base Logarithm base.
Returns:
CPL_ERROR_NONE or the relevant _cpl_error_code_ on error
The base and all the vector elements must be positive and the base must be different from 1, or a cpl_error_code will be returned and the vector will be left unmodified.

Possible _cpl_error_code_ set in this function:

cpl_error_code cpl_vector_multiply ( cpl_vector *  v1,
const cpl_vector *  v2 
)

Multiply two vectors component-wise.

Parameters:
v1 First cpl_vector
v2 Second cpl_vector
Returns:
CPL_ERROR_NONE or the relevant _cpl_error_code_ on error
See also:
cpl_vector_add()

cpl_error_code cpl_vector_multiply_scalar ( cpl_vector *  v,
double  factor 
)

Elementwise multiplication of a vector with a scalar.

Parameters:
v cpl_vector to modify
factor Number to multiply with
Returns:
CPL_ERROR_NONE or the relevant _cpl_error_code_ on error
Multiply each element of the cpl_vector with a number.

Possible _cpl_error_code_ set in this function:

cpl_vector* cpl_vector_new ( int  n  ) 

Create a new cpl_vector.

Parameters:
n Number of element of the cpl_vector
Returns:
1 newly allocated cpl_vector or NULL in case of an error
The returned object must be deallocated using cpl_vector_delete(). There is no default values assigned to the created object, they are undefined until they are set.

Possible _cpl_error_code_ set in this function:

cpl_error_code cpl_vector_power ( cpl_vector *  v,
double  exponent 
)

Compute the power of all vector elements.

Parameters:
v Target cpl_vector.
exponent Constant exponent.
Returns:
CPL_ERROR_NONE or the relevant _cpl_error_code_ on error
If the exponent is negative all vector elements must be non-zero and if the exponent is non-integer all vector elements must be non-negative, otherwise a cpl_error_code is returned and the vector is unmodified.

Following the behaviour of C99 pow() function, this function sets 0^0 = 1.

Possible _cpl_error_code_ set in this function:

double cpl_vector_product ( const cpl_vector *  v1,
const cpl_vector *  v2 
)

Compute the vector dot product.

Parameters:
v1 One vector
v2 Another vector of the same size
Returns:
The (non-negative) product or negative on error.
The same vector may be passed twice, in which case the square of its 2-norm is computed.

Possible _cpl_error_code_ set in this function:

cpl_vector* cpl_vector_read ( const char *  filename  ) 

Read a list of values from an ASCII file and create a cpl_vector.

Parameters:
filename Name of the input ASCII file
Returns:
1 newly allocated cpl_vector or NULL in case of an error
See also:
cpl_vector_dump
Parse an input ASCII file values and create a cpl_vector from it Lines beginning with a hash are ignored, blank lines also. In valid lines the value is preceeded by an integer, which is ignored.

The returned object must be deallocated using cpl_vector_delete()

In addition to normal files, FIFO (see man mknod) are also supported.

Possible _cpl_error_code_ set in this function:

cpl_error_code cpl_vector_save ( const cpl_vector *  to_save,
const char *  filename,
cpl_type_bpp  bpp,
const cpl_propertylist pl,
unsigned  mode 
)

Save a vector to a FITS file.

Parameters:
to_save Vector to write to disk or NULL
filename Name of the file to write
bpp Bits per pixel
pl Property list for the output header or NULL
mode The desired output options (combined with bitwise or)
Returns:
CPL_ERROR_NONE or the relevant _cpl_error_code_ on error
This function saves a vector to a FITS file (NAXIS=1), using cfitsio. If a property list is provided, it is written to the named file before the pixels are written. If the image is not provided, the created file will only contain the primary header. This can be useful to create multiple extension files.

The requested pixel depth (bpp) follows the FITS convention. Possible values are CPL_BPP_8_UNSIGNED (8), CPL_BPP_16_SIGNED (16), CPL_BPP_32_SIGNED (32), CPL_BPP_IEEE_FLOAT (-32), CPL_BPP_IEEE_DOUBLE (-64)

Supported output modes are CPL_IO_CREATE (create a new file) and CPL_IO_EXTEND (append to an existing file)

If you are in append mode, make sure that the file has writing permissions. You may have problems if you create a file in your application and append something to it with the umask set to 222. In this case, the file created by your application would not be writable, and the append would fail.

Possible _cpl_error_code_ set in this function:

cpl_error_code cpl_vector_set ( cpl_vector *  in,
int  idx,
double  value 
)

Set an element of the vector.

Parameters:
in the input vector
idx the index of the element (0 to nelem-1)
value the value to set in the vector
Returns:
CPL_ERROR_NONE or the relevant _cpl_error_code_ on error
Possible _cpl_error_code_ set in this function:

cpl_error_code cpl_vector_set_size ( cpl_vector *  in,
int  newsize 
)

Resize the vector.

Parameters:
in The vector to be resized
newsize The new (positive) number of elements in the vector
Returns:
CPL_ERROR_NONE or the relevant _cpl_error_code_ on error
Note:
On succesful return the value of the elements of the vector is unchanged to the minimum of the old and new sizes; any newly allocated elements are undefined. The pointer to the vector data buffer may change, therefore pointers previously retrieved by calling cpl_vector_get_data() should be discarded. If the vector was created with cpl_vector_wrap() the argument pointer to that call should be discarded as well.
Possible _cpl_error_code_ set in this function:

cpl_error_code cpl_vector_sort ( cpl_vector *  v,
int  c 
)

Sort a cpl_vector by increasing/decreasing data.

Parameters:
v cpl_vector to sort
c Sorting criterion to use
Returns:
CPL_ERROR_NONE or the relevant _cpl_error_code_ on error
The input cpl_vector is modified to sort out all values following the sorting criterion. Possible sorting criteria are:

Possible _cpl_error_code_ set in this function:

cpl_error_code cpl_vector_sqrt ( cpl_vector *  v  ) 

Compute the sqrt of a cpl_vector.

Parameters:
v cpl_vector
Returns:
CPL_ERROR_NONE or the relevant _cpl_error_code_ on error
The sqrt of the data is computed. The input cpl_vector is modified

If an element in v is negative v is not modified and CPL_ERROR_ILLEGAL_INPUT is returned.

Possible _cpl_error_code_ set in this function:

cpl_error_code cpl_vector_subtract ( cpl_vector *  v1,
const cpl_vector *  v2 
)

Subtract a cpl_vector from another.

Parameters:
v1 First cpl_vector (modified)
v2 Second cpl_vector
Returns:
CPL_ERROR_NONE or the relevant _cpl_error_code_ on error
See also:
cpl_vector_add()

cpl_error_code cpl_vector_subtract_scalar ( cpl_vector *  v,
double  subtrahend 
)

Elementwise subtraction of a scalar from a vector.

Parameters:
v cpl_vector to modify
subtrahend Number to subtract
Returns:
CPL_ERROR_NONE or the relevant _cpl_error_code_ on error
Subtract a number from each element of the cpl_vector.

Possible _cpl_error_code_ set in this function:

void* cpl_vector_unwrap ( cpl_vector *  v  ) 

Delete a cpl_vector except the data array.

Parameters:
v cpl_vector to delete
Returns:
A pointer to the data array or NULL if the input is NULL.
Note:
The data array must subsequently be deallocated. Failure to do so will result in a memory leak.

cpl_vector* cpl_vector_wrap ( int  n,
double *  data 
)

Create a cpl_vector from existing data.

Parameters:
n Number of elements in the vector
data Pointer to array of n doubles
Returns:
1 newly allocated cpl_vector or NULL in case of an error
The returned object must be deallocated using cpl_vector_unwrap().

Possible _cpl_error_code_ set in this function:


Generated on Wed Mar 18 09:40:13 2009 for Common Pipeline Library Reference Manual by  doxygen 1.4.7