gsl_fit Class Template Reference

Non-linear least-squares fitting class (GSL). More...

#include <gsl_fit.h>

Inheritance diagram for gsl_fit:

fit_base

Detailed Description

template<class param_t, class func_t = fit_funct<param_t>, class vec_t = ovector_base, class mat_t = omatrix_base, class bool_vec_t = bool *>
class gsl_fit< param_t, func_t, vec_t, mat_t, bool_vec_t >

The GSL-based fitting class using a Levenberg-Marquardt type algorithm. The algorithm stops when

\[ |dx_i| < \mathrm{epsabs}+\mathrm{epsrel}\times|x_i| \]

where $dx$ is the last step and $x$ is the current position. If test_gradient is true, then additionally fit() requires that

\[ \sum_i |g_i| < \mathrm{epsabs} \]

where $g_i$ is the $i$-th component of the gradient of the function $\Phi(x)$ where

\[ \Phi(x) = || F(x) ||^2 \]

Todo:
Properly generalize other vector types than ovector_base
Todo:
Allow the user to specify the derivatives
Todo:
Fix so that the user can specify automatic scaling of the fitting parameters, where the initial guess are used for scaling so that the fitting parameters are near unity.

Definition at line 87 of file gsl_fit.h.


Data Structures

struct  func_par
 A structure for passing to the functions func(), dfunc(), and fdfunc(). More...

Public Member Functions

virtual int fit (size_t ndat, vec_t &xdat, vec_t &ydat, vec_t &yerr, size_t npar, vec_t &par, mat_t &covar, double &chi2, param_t &pa, func_t &fitfun)
 Fit the data specified in (xdat,ydat) to the function fitfun with the parameters in par.
virtual const char * type ()
 Return string denoting type ("gsl_fit").

Data Fields

int max_iter
 (default 500)
double epsabs
 (default 1.0e-4)
double epsrel
 (default 1.0e-4)
bool test_gradient
 If true, test the gradient also (default false).
bool use_scaled
 Use the scaled routine if true (default true).

Protected Member Functions

virtual int print_iter (int nv, gsl_vector *x, gsl_vector *dx, int iter, double l_epsabs, double l_epsrel)
 Print the progress in the current iteration.

Static Protected Member Functions

static int func (const gsl_vector *x, void *pa, gsl_vector *f)
 Evaluate the function.
static int dfunc (const gsl_vector *x, void *pa, gsl_matrix *jac)
 Evaluate the jacobian.
static int fdfunc (const gsl_vector *x, void *pa, gsl_vector *f, gsl_matrix *jac)
 Evaluate the function and the jacobian.

Member Function Documentation

virtual int fit ( size_t  ndat,
vec_t &  xdat,
vec_t &  ydat,
vec_t &  yerr,
size_t  npar,
vec_t &  par,
mat_t &  covar,
double &  chi2,
param_t &  pa,
func_t &  fitfun 
) [inline, virtual]

The covariance matrix for the parameters is returned in covar and the value of $ \chi^2 $ is returned in chi2.

Implements fit_base.

Definition at line 257 of file gsl_fit.h.


The documentation for this class was generated from the following file:

Documentation generated with Doxygen and provided under the GNU Free Documentation License. See License Information for details.

Project hosting provided by SourceForge.net Logo, O2scl Sourceforge Project Page