Marquardt-Levenberg fit to a generic function. More...

Public Member Functions
	LMFit (const TFUNC &func, const int itermax=1000, const TPAR tol=1e-7, const FVector< bool, NPAR > ignin=false, const TPAR astart=1e-3, const TPAR astep=10.0)
	Construct class with fit constraints. More...

void	setIgnore (const FVector< bool, NPAR > &ignin)
	specify which parameters to leave constant during fit [ignin(i)=true]. More...

FVector< TPAR, NPAR >	getResult ()
	Return result vector. More...

double	getChiSquare () const
	Return final $\chi^2$ . More...

int	getNIteration () const
	Return No needed iterations. More...

FVector< TPAR, NPAR >	getVariance ()
	Return diagonal of covariance matrix. More...

FMatrix< TPAR, NPAR, NPAR >	getCovarianceMatrix ()
	Return diagonal of covariance matrix. More...

FVector< TPAR, NPAR >	getErrors ()
	Return the formal errors of the fit parameters. More...

template<typename TDATAX , typename TDATAY >
void	eval (const TDATAX &x, const TDATAY &y, const TDATAY &dy2, const typename TDATAY::value_type nan_y, const FVector< TPAR, NPAR > &inpar)
	Fit to data and ignoring nan, start with inpar. More...

Detailed Description

template<typename TFUNC, typename TPAR, int NPAR, typename Solver = GaussJ<TPAR,NPAR>>
class ltl::LMFit< TFUNC, TPAR, NPAR, Solver >

Marquardt-Levenberg fit to a generic function.

Class to fit an arbitrary function to data and errors by non-linear least squares using the Marquardt-Levenberg algorithm. Either Gauss-Jordan or LU decomposition can be used to solve linear systems during fitting, controlled by a template parameter.

The function object should provide a typedef of value_type and

operator()

. Using TPAR to denote the type of the parameters, NPAR their number, and T the type of the argument:

class function
{
   public:
      typedef T value_type;
      value_type operator()( const T x, const FVector<TPAR,NPAR>& p, FVector<TPAR,NPAR>& df_dpi ) const
      {
         value_type value = ...;
         for( int i=1; i<=NPAR; ++i )
            df_dpi(i) = ...;
         return value;
      }
}

Here's an example class that fits a quadratic function:

class function
{
   public:
      typedef float value_type;
      float operator()( const float x, const FVector<float,3>& p, FVector<float,3>& df_dpi ) const
      {
         float value = p(1)*x*x + p(2)*x + p(3);
         df_dpi(1) = x*x;
         df_dpi(2) = x;
         df_dpi(3) = 1.0f;
         return value;
      }
};

This function is used in the following way:

// just to illustrate the template parameters:
// TFUNC is the function object,
// TPAR the type of the parameters (needs not be the same as the input or return type of the function),
// NPAR the number of parameters.
// LinSolver is the solver for the linear systems during fitting. Compatible objects are
// GaussJ (Gauss-Jordan elimination) and LUDecomposition (LU decomposition/SVD).
// The former is faster, the latter stable against numerically singular matrices.
template<typename TFUNC, typename TPAR, int NPAR, typename LinSolver=GaussJ<NPAR,NPAR> >
class LMFit;
function F;
LMFit< function, float, 3> M( F );
MArray<float,1> X(10);
X = 1.27355, -0.654883, 3.7178, 2.31818, 2.6652, -2.02182, 4.82368, -4.36208, 4.84084, 2.44391;
MArray<float,1> Y(10);
Y = 4.30655,0.420333,18.6992,8.27967,10.8136,3.3598,29.3496,15.9234,29.4432,9.6158;
MArray<float,1> dY2(10);
dY2 = 1.0f;
FVector<float,3> P0;   // initial guess.
P0 = 0.3, -1.0, 0.0;
M.eval( X, Y, dY2, -9999.0f, P0 );   // perform the fitting.
cout << M.getResult() << endl;
// the same, but using LU decomposition to invert and solve linear systems:
LMFit< function, float, 3, LUDecomposition<float,3> > M( F );
...

About any container can be used for passing X, Y, and dY data, as long as the container provides

typedef Container::const_itertor ...

,

typedef Container::value_type ...

,

Container::value_type

is the same as the type of the first argument to operator() of the function to be fit. For example, this might well be

FVector<float,2>

, and the container an

MArray<FVector<float,2> >

to fit a function of 2 variables (and N parameters). This might be a polynomial with crossterms.

Constructor & Destructor Documentation

template<typename TFUNC , typename TPAR , int NPAR, typename Solver = GaussJ<TPAR,NPAR>>

ltl::LMFit< TFUNC, TPAR, NPAR, Solver >::LMFit	(	const TFUNC &	func,
		const int	itermax = `1000`,
		const TPAR	tol = `1e-7`,
		const FVector< bool, NPAR >	ignin = `false`,
		const TPAR	astart = `1e-3`,
		const TPAR	astep = `10.0`
	)

inline

Construct class with fit constraints.

Member Function Documentation

template<typename TFUNC , typename TPAR , int NPAR, typename Solver = GaussJ<TPAR,NPAR>>

void ltl::LMFit< TFUNC, TPAR, NPAR, Solver >::setIgnore ( const FVector< bool, NPAR > & ignin )

inline

specify which parameters to leave constant during fit [ignin(i)=true].

References ltl::anyof().

template<typename TFUNC , typename TPAR , int NPAR, typename Solver >

template<typename TDATAX , typename TDATAY >

void ltl::LMFit< TFUNC, TPAR, NPAR, Solver >::eval	(	const TDATAX &	x,
		const TDATAY &	y,
		const TDATAY &	dy2,
		const typename TDATAY::value_type	nan_y,
		const FVector< TPAR, NPAR > &	inpar
	)

Fit to data and $error^2$ ignoring nan, start with inpar.

template<typename TFUNC , typename TPAR , int NPAR, typename Solver = GaussJ<TPAR,NPAR>>

FVector<TPAR,NPAR> ltl::LMFit< TFUNC, TPAR, NPAR, Solver >::getResult ( )

inline

Return result vector.

template<typename TFUNC , typename TPAR , int NPAR, typename Solver = GaussJ<TPAR,NPAR>>

double ltl::LMFit< TFUNC, TPAR, NPAR, Solver >::getChiSquare ( ) const

inline

Return final $\chi^2$ .

template<typename TFUNC , typename TPAR , int NPAR, typename Solver = GaussJ<TPAR,NPAR>>

int ltl::LMFit< TFUNC, TPAR, NPAR, Solver >::getNIteration ( ) const

inline

Return No needed iterations.

template<typename TFUNC , typename TPAR , int NPAR, typename Solver = GaussJ<TPAR,NPAR>>

FVector<TPAR,NPAR> ltl::LMFit< TFUNC, TPAR, NPAR, Solver >::getVariance ( )

inline

Return diagonal of covariance matrix.

References ltl::FMatrix< T, M, N >::traceVector().

template<typename TFUNC , typename TPAR , int NPAR, typename Solver = GaussJ<TPAR,NPAR>>

FMatrix<TPAR, NPAR, NPAR> ltl::LMFit< TFUNC, TPAR, NPAR, Solver >::getCovarianceMatrix ( )

inline

Return diagonal of covariance matrix.

template<typename TFUNC , typename TPAR , int NPAR, typename Solver = GaussJ<TPAR,NPAR>>

FVector<TPAR,NPAR> ltl::LMFit< TFUNC, TPAR, NPAR, Solver >::getErrors ( )

inline

Return the formal errors of the fit parameters.

References ltl::FMatrix< T, M, N >::traceVector().

Public Member Functions

Detailed Description

template<typename TFUNC, typename TPAR, int NPAR, typename Solver = GaussJ<TPAR,NPAR>> class ltl::LMFit< TFUNC, TPAR, NPAR, Solver >

Constructor & Destructor Documentation

Member Function Documentation

template<typename TFUNC, typename TPAR, int NPAR, typename Solver = GaussJ<TPAR,NPAR>>
class ltl::LMFit< TFUNC, TPAR, NPAR, Solver >