int lsqline( // Data in: const Vector& x, // ordinates const Vector& y, // y(x) data values const int fixb, // 1: define bhat / 0: fit it const double eps, // numerical bandwidth to prevent overflow on division // Data in/out: double& bhat, // best-fitting (or supplied) intercept // Data out: double& ahat, // best-fitting gradient double& sd // standard deviation )Returns 1 on success; 0 if failed because gradient or intercept calculation would involve division by small number.
ArrayReturns polynomial coefficients (order of increasing powers of x), or an Array of length zero if it fails.polyfitg( const Array & x, // Array of x-ordinates const Array & y, // Array of y-values const Array & w, // Array of weights const int n, // number of polynomial coefficients, i.e. order + 1 double& res, // residual of l2 fit (returned value) const double eps = 1.0e-8 // numerical bandwidth for spotting singular matrix )
polyfitg
is the recommended function.
ArrayReturns polynomial coefficients (order of increasing powers of x), or an Array of length zero if it fails.polyfitn( const Array & x, // Array of x-ordinates const Array & y, // Array of y-values const Array & w, // Array of weights const int n, // number of polynomial coefficients, i.e. order + 1 const double eps = 1.0e-8 // numerical bandwidth for spotting singular matrix )
Arraylinfit( const Matrix& xf, // xf[i][j]: function j's value at data point i const Vector& y, // data values const Vector& w, // weight of each data value double& res, // residual of fit const double eps = 1.0e-8 // numerical bandwidth: singular matrix )
ArrayReturns coefficients for the functions {1,x,y,xy} respectively, or an Array of no elements if the fitting fails.frectfit( const Array & x, // x ordinates of corners (2) const Array & y, // y ordinates of corners (2) const matrix& f, // function values (2x2) const eps = 1.0e-8 // numerical bandwidth )
int fitstart( const Array(returns index of first point)& r, // Array of ordinates const double rr, // value around which fit is to be performed const int nr // number of points required for the fit )
datafit
acts as a wrapper for the fitting of general
functions (especially ones nonlinear in the adjustable parameters) to data.
This class provides alternative parameter-fitting methods.
To use, a program defines a subclass inheriting from datafit
,
and providing instantiations of essential methods describing the fitting
functions and their derivatives.
The exact functions required depend on the method used.
Specifications of the functions are:
// evaluate function given coefficient array and index of data point double func(const Array& co, const int ipt) const // evaluate row of Jacobian matrix: derivative of function with respect // to coefficients, at specified data point Vector jacob(const Array & co, const int ipt) const // evaluate Hessian matrix: second derivatives of function with respect // to coefficients, at specified data point Matrix jacob(const Array & co, const int ipt) const
Fitting methods are described below:
int gaussnewton( const Array& y, // data values const Array & w, // weight of each data value Array & co, // coefficients: guess in, optimised out int& hitmax, // reached maximum iterations before finding solution int& sing, // Jacobian was close to singular: could not proceed int& stuck, // could not find valid parameter change int& nits, // iterations taken double& el2norm, // least squares norm in final answer const double acc = 1.0e-6, // accuracy const double eps = 1.0e-8, // bandwidth: singular matrix or tiny steps const double dampf = 0.1, // damping factor const double dampfac = 0.5, // change factor for damping factor const int maxits = 0, // maximum iterations (0 => no limit) const int verbose = 0 // diagnostic output )
int levmar( const Array& y, // data values const Array & w, // weight of each data value Array & co, // coefficients: guess in, optimised out int& hitmax, // reached maximum iterations before finding solution int& sing, // Jacobian was close to singular: could not proceed int& stuck, // could not find valid parameter change int& nits, // iterations taken double& el2norm, // least squares norm in final answer const double tregin, // initial trust region const double gamstep, // initial jump in growth parameter const double dfac = 2, // change factor for growth and trust region const double acc = 1.0e-6, // accuracy const double eps = 1.0e-8, // bandwidth: singular matrix or tiny steps const int bothways = 1, // to alternate direction of search for robustness const int maxits = 0, // maximum iterations (0 => no limit) const int verbose = 0 // diagnostic output )
int gaussnewtonhess( const Array& y, // data values const Array & w, // weight of each data value Array & co, // coefficients: guess in, optimised out int& hitmax, // reached maximum iterations before finding solution int& sing, // Jacobian was close to singular: could not proceed int& stuck, // could not find valid parameter change int& nits, // iterations taken double& el2norm, // least squares norm in final answer const double cutoff = 0.1, // error below which Hessian is not calculated const double acc = 1.0e-6, // accuracy const double eps = 1.0e-8, // bandwidth: singular matrix or tiny steps const double dampf = 0.1, // damping factor const double dampfac = 0.5, // change factor for damping factor const int maxits = 0, // maximum iterations (0 => no limit) const int verbose = 0 // diagnostic output )