23 #ifndef O2SCL_FIT_BASE_H
24 #define O2SCL_FIT_BASE_H
32 #include <o2scl/jacobian.h>
33 #include <o2scl/fparser.h>
34 #include <o2scl/mm_funct.h>
36 #ifndef DOXYGEN_NO_O2NS
41 typedef std::function<
50 template<
class vec_t=boost::numeric::ublas::vector<
double> >
58 std::string var,
int nauxp=0,
59 std::string auxp=
"") {
61 std::string both=parms+
","+var;
62 fpw.Parse(formula,both);
66 std::string all=parms+
","+var+
","+auxp;
67 fpw.Parse(formula,all);
87 for(
int i=0;i<naux;i++) {
95 virtual double operator()(
size_t np,
const vec_t &p,
double x) {
100 double *both=
new double[np+1];
101 for(i=0;i<np;i++) both[i]=p[i];
106 double *all=
new double[np+1+naux];
107 for(i=0;i<np;i++) all[i]=p[i];
109 for(i=np+1;i<np+1+naux;i++) all[i]=aux[i-np-1];
116 #ifndef DOXYGEN_INTERNAL
140 int set_function(std::string formula, std::string parms,
141 std::string var,
int nauxp=0, std::string auxp=
"");
158 template<
class vec_t=boost::numeric::ublas::vector<
double>,
159 class mat_t=boost::numeric::ublas::matrix<
double> >
171 virtual void operator()(
size_t np,
const vec_t &p,
size_t nd,
177 virtual void jac(
size_t np, vec_t &p,
size_t nd, vec_t &f,
183 #ifndef DOXYGEN_INTERNAL
222 template<
class vec_t=boost::numeric::ublas::vector<
double>,
223 class mat_t=boost::numeric::ublas::matrix<
double>,
233 const vec_t &yerr, fit_func_t &fun) {
240 mfm=std::bind(std::mem_fn<
int(
size_t,
const vec_t &,vec_t &)>
242 std::placeholders::_1,std::placeholders::_2,
243 std::placeholders::_3);
247 double sqrt_dbl_eps=sqrt(std::numeric_limits<double>::epsilon());
253 void set_data(
size_t ndat,
const vec_t &xdat,
const vec_t &ydat,
271 virtual double chi2(
size_t np,
const vec_t &p) {
273 for(
size_t i=0;i<ndat_;i++) {
274 double yi=((*fun_)(np,p,(*xdat_)[i])-(*ydat_)[i])/((*yerr_)[i]);
283 virtual void operator()(
size_t np,
const vec_t &p,
size_t nd, vec_t &f) {
285 for(
size_t i=0;i<nd;i++) {
286 double yi=(*fun_)(np,p,(*xdat_)[i]);
287 f[i]=(yi-(*ydat_)[i])/((*yerr_)[i]);
295 virtual void jac(
size_t np, vec_t &p,
size_t nd, vec_t &f,
312 #ifndef DOXYGEN_INTERNAL
323 std::function<int(size_t,const vec_t &,vec_t &)>
mfm;
384 virtual int print_iter(
size_t nv, vec_t &x,
double y,
int iter,
385 double value=0.0,
double limit=0.0) {
391 std::cout <<
"Iteration: " << iter << std::endl;
393 for(i=0;i<nv;i++) std::cout << x[i] <<
" ";
394 std::cout << std::endl;
395 std::cout <<
"y: " << y <<
" Val: " << value <<
" Lim: " << limit
398 std::cout <<
"Press a key and type enter to continue. ";
412 virtual int fit(
size_t npar, vec_t &parms, mat_t &covar,
413 double &chi2, func_t &fitfun)=0;
420 virtual const char *
type() {
return "fit_base"; }
431 #ifndef DOXYGEN_NO_O2NS
Parse a mathematical function specified in a string.
Generalized fitting function [abstract base].
virtual int set_function(func_t &f)
Set the function to compute the Jacobian of.
virtual void jac(size_t np, vec_t &p, size_t nd, vec_t &f, mat_t &J)
Using parameters in p, compute the Jacobian in J.
void set_func(fit_func_t &fun)
Set the fitting function.
Standard fitting function based on one-dimensional data with a numerical Jacobian.
int set_aux_parms(vec_t &ap)
Set the values of the auxilliary parameters that were specified in auxp in the constructor.
virtual void jac(size_t np, vec_t &p, size_t nd, vec_t &f, mat_t &J)=0
Using parameters in p, compute the Jacobian in J.
virtual const char * type()
Return string denoting type ("fit_base")
size_t n_dat
The number of data points.
void set_data(size_t ndat, const vec_t &xdat, const vec_t &ydat, const vec_t &yerr)
Set the data to be fit.
std::function< double(size_t, const boost::numeric::ublas::vector< double > &, double)> fit_funct11
Array of multi-dimensional functions typedef (C++11 version)
double epsrel
The relative stepsize for finite-differencing (default )
size_t ntrial
Maximum number of iterations (default 500)
virtual void operator()(size_t np, const vec_t &p, size_t nd, vec_t &f)=0
Using parameters in p, compute the relative deviations in f.
int set_function(std::string formula, std::string parms, std::string var, int nauxp=0, std::string auxp="")
Specify the strings which define the fitting function.
fit_func_t * fun_
Fitting function.
std::function< int(size_t, const vec_t &, vec_t &)> mfm
Function object for Jacobian object.
double tol_rel
(default 1.0e-4)
virtual int print_iter(size_t nv, vec_t &x, double y, int iter, double value=0.0, double limit=0.0)
Print out iteration information.
int verbose
An integer describing the verbosity of the output.
std::string st_form
The formula.
std::string st_parms
The parameters.
size_t n_par
The number of parameters.
std::string st_var
The variables.
FunctionParser fpw
The function evaluation object.
virtual size_t get_ndata()
Return the number of data points.
virtual size_t get_ndata()=0
Return the number of data points.
virtual void operator()(size_t np, const vec_t &p, size_t nd, vec_t &f)
Using parameters in p, compute the relative deviations in f.
Non-linear least-squares fitting [abstract base].
virtual double chi2(size_t np, const vec_t &p)
Return .
chi_fit_funct(size_t ndat, const vec_t &xdat, const vec_t &ydat, const vec_t &yerr, fit_func_t &fun)
Create an object with specified data and specified fitting function.
jacobian_gsl< std::function< int(size_t, const vec_t &, vec_t &)>, vec_t, mat_t > auto_jac
Automatic Jacobian object.
double tol_abs
Absolute tolerance (default 1.0e-4)
int jac_mm_funct(size_t np, const vec_t &p, vec_t &f)
Reformulate operator() into a mm_funct11 object.
fit_funct11_strings(std::string formula, std::string parms, std::string var, int nauxp=0, std::string auxp="")
Specify a fitting function through a string.
Simple automatic Jacobian.
virtual int fit(size_t npar, vec_t &parms, mat_t &covar, double &chi2, func_t &fitfun)=0
Fit function fitfun using parameters in parms as initial guesses.
virtual double operator()(size_t np, const vec_t &p, double x)
Using parameters in p, predict y given x.