26 #ifndef __Minpack_hh__
27 #define __Minpack_hh__
29 #include <bias_config.h>
31 #include <Base/Math/Matrix.hh>
32 #include <Base/Math/Vector.hh>
33 #include <MathAlgo/minpack/cminpack.h>
37 #define MINPACK_DEF_TOLERANCE 1e-7
38 #define LM_DEF_MAX_ITER -1
39 #define LM_DEF_EPSILON 1e-15
41 #define LM_IMPROPER_INPUT -1
43 #define LM_TOO_MANY_IT 5
44 #define LM_SMALL_FTOL 6.
45 #define LM_SMALL_XTOL 7
46 #define LM_SMALL_GTOL 8
96 long int Powell(minpack_func_nn fun,
void *p,
97 Vector<double>& InitialGuess, Vector<double>& Result,
98 double Tolerance = MINPACK_DEF_TOLERANCE);
113 Vector<double>& InitialGuess, Vector<double>& Result,
114 double EpsilonFun = LM_DEF_EPSILON,
115 double Tolerance = MINPACK_DEF_TOLERANCE);
169 Vector<double>& InitialGuess,
170 Vector<double>& Result,
171 double Tolerance = MINPACK_DEF_TOLERANCE);
222 Vector<double>& InitialGuess,
223 Vector<double>& Result,
224 double Tolerance = MINPACK_DEF_TOLERANCE,
225 long int MaxIter = LM_DEF_MAX_ITER);
265 Vector<double>& InitialGuess,
266 Vector<double>& Result,
267 double Tolerance = MINPACK_DEF_TOLERANCE,
268 long int MaxIter = LM_DEF_MAX_ITER);
325 Vector<double>& InitialGuess,
326 Vector<double>& Result,
327 double EpsilonFun = LM_DEF_EPSILON,
328 double Tolerance = MINPACK_DEF_TOLERANCE,
329 long int MaxIter = LM_DEF_MAX_ITER,
330 double* ResultingJacobian = NULL,
332 double *SumOfSquaredErrors = NULL);
351 int *ipvt,
double &SumOfSquaredErrors,
352 Matrix<double>& Cov);
373 const Matrix<double> &Jac,
374 const Vector<int>& Permutation,
375 const double &SumOfSquaredErrors,
376 Matrix<double>& Cov);
413 const long int NumErrors,
const Vector<double> &x,
414 Vector<double> &error);
428 const int NumErrors,
const Vector<double> &x,
430 const double gradientEpsilon = LM_DEF_EPSILON);
441 void QRFrac(
const Matrix<double> &A, Matrix<double>& Q, Matrix<double>& R);
445 #endif // __Minpack_hh__
long int LevenbergMarquardtExtended(minpack_funcder_mn fun, void *p, long int m, long int n, Vector< double > &InitialGuess, Vector< double > &Result, double Tolerance, long int MaxIter)
Solves an overdetermined system f(x) = 0 of m non-linear equations with n parameters using the extend...
long int ComputeCovariance(long int NumErrors, long int NumParams, double *Jac, int *Permutation, double &SumOfSquaredErrors, Matrix< double > &Cov)
Compute covariance matrix from Levenberg-Marquardt resulting Jacobian matrix J(x) and permutation mat...
void QRFrac(const Matrix< double > &A, Matrix< double > &Q, Matrix< double > &R)
Compute the QR decomposition of input matrix A using Householder transformations without column pivot...
void CheckJacobian(minpack_funcder_mn fun, void *p, const long int NumErrors, const Vector< double > &x, Vector< double > &error)
Convenience function to check correctness of user supplied Jacobians J(x) needed by the analytical ve...
long int ComputeJacobian(minpack_func_mn fun, void *p, const int NumErrors, const Vector< double > &x, Matrix< double > &fjac, const double gradientEpsilon)
Compute the Jacobian J(x) of function f at x.
long int LevenbergMarquardt(minpack_funcder_mn fun, void *p, long int m, long int n, Vector< double > &InitialGuess, Vector< double > &Result, double Tolerance)
Solves an overdetermined system f(x) = 0 of m non-linear equations with n parameters using the Levenb...
long int Powell(minpack_func_nn fun, void *p, Vector< double > &InitialGuess, Vector< double > &Result, double Tolerance)
Solve a system f(x) = 0 of n non-linear equations with n unknowns using the Powell hybrid method...
long int PowellExtended(minpack_func_nn fun, void *p, Vector< double > &InitialGuess, Vector< double > &Result, double EpsilonFun, double Tolerance)
Solve a system f(n) = 0 of n non-linear equations with n unknowns using the Powell hybrid method...