17#include <gsl/gsl_multifit_nlin.h>
27 : m_dirtyVal(true), m_dirtyDeriv(true), m_dirtyHessian(true), m_includePenalty(true),
m_value(0), m_pushed(false),
28 m_pushedValue(false) {}
116 throw std::runtime_error(
"Fitting cost function isn't set");
127 const auto &j = J.
getJ();
128 assert(
static_cast<size_t>(j.cols()) == na);
137 tempJTr = j.transpose();
141 gsl_multifit_covar(&tempJTrGSL.matrix, epsrel, &covarTrGSL.matrix);
142 covar = covarTr.
tr();
157 bool isTransformationIdentity =
true;
158 for (
size_t i = 0; i < np; ++i) {
161 isTransformationIdentity =
165 if (isTransformationIdentity) {
172 covar = tm.
tr() * c * tm;
188 for (
size_t i = 0; i < np; ++i) {
193 for (
size_t j = 0; j < np; ++j) {
195 (*covarMatrix)[i][j] = covar.
get(ia, ja);
198 if (ja >= covar.
size2())
202 double err = sqrt(covar.
get(ia, ia));
206 if (ia >= covar.
size1())
211 m_function->setReducedChiSquared(chi2 /
static_cast<double>((
m_values->size() - np)));
219 const double epsilon = std::numeric_limits<double>::epsilon() * 100;
224 for (
size_t i = 0; i < np; ++i) {
228 for (
size_t j = 0; j < na; ++j) {
230 double step = ap == 0.0 ? epsilon : ap * epsilon;
269 if (np != params.
size()) {
272 for (
size_t i = 0; i < np; ++i) {
284 if (params.
size() != np) {
287 for (
size_t i = 0; i < np; ++i) {
304 auto seqDomain = std::dynamic_pointer_cast<SeqDomain>(
m_domain);
307 seqDomain->additiveCostFunctionVal(*
this);
310 throw std::runtime_error(
"CostFunction: undefined FunctionValues.");
317 for (
size_t i = 0; i <
m_function->nParams(); ++i) {
340 for (
size_t i = 0; i <
nParams(); ++i) {
355 for (
size_t i = 0; i <
nParams(); ++i) {
375 const size_t numParams =
nParams();
388 auto seqDomain = std::dynamic_pointer_cast<SeqDomain>(
m_domain);
391 seqDomain->additiveCostFunctionValDerivHessian(*
this, evalDeriv, evalHessian);
394 throw std::runtime_error(
"CostFunction: undefined FunctionValues.");
402 for (
size_t i = 0; i < np; ++i) {
414 for (
size_t ip = 0; ip < np; ++ip) {
429 for (
size_t ip = 0; ip < np; ++ip) {
478 throw std::runtime_error(
"Cost Function: double push.");
491 throw std::runtime_error(
"Cost Function: empty stack.");
506 throw std::runtime_error(
"Cost Function: empty stack.");
const std::string & m_value
double value
The value of the point.
An interface to a constraint.
virtual double checkDeriv()=0
Returns the derivative of the penalty for each active parameter.
virtual double checkDeriv2()=0
Returns the derivative of the penalty for each active parameter.
virtual double check()=0
Returns a penalty number which is bigger than or equal to zero If zero it means that the constraint i...
virtual void setParamToSatisfyConstraint()=0
Set the parameters of IFitFunction to satisfy constraint.
void setDirty()
Set all dirty flags.
bool m_includePenalty
Flag to include constraint in cost function value.
bool m_dirtyHessian
dirty hessian flag
API::FunctionValues_sptr m_values
Shared poinetr to the function values.
virtual void deriv(std::vector< double > &der) const override
Calculate the derivatives of the cost function.
void setParameter(size_t i, const double &value) override
Set i-th parameter.
virtual void calFittingErrors(const EigenMatrix &covar, double chi2)
Calculate fitting errors.
virtual double valDerivHessian(bool evalDeriv=true, bool evalHessian=true) const
Calculate the value, the first and the second derivatives of the cost function.
void pop()
Restore saved parameters, derivatives and hessian.
void reset() const
Reset the fitting function (neccessary if parameters get fixed/unfixed)
std::vector< size_t > m_indexMap
maps the cost function's parameters to the ones of the fitting function.
API::FunctionDomain_sptr m_domain
Shared pointer to the function domain.
void setParameters(const EigenVector ¶ms)
Set all parameters.
void push()
Save current parameters, derivatives and hessian.
virtual void setFittingFunction(API::IFunction_sptr function, API::FunctionDomain_sptr domain, API::FunctionValues_sptr values)
Set fitting function.
void getParameters(EigenVector ¶ms) const
Get all parameters into a GSLVector.
const EigenMatrix & getHessian() const
Return cached or calculate the Hessian.
virtual void addVal(API::FunctionDomain_sptr domain, API::FunctionValues_sptr values) const =0
Increment to the cost function by evaluating it on a domain.
size_t nParams() const override
Number of parameters.
void calTransformationMatrixNumerically(EigenMatrix &tm)
Calculate the transformation matrix T by numeric differentiation.
API::IFunction_sptr m_function
Shared pointer to the fitting function.
virtual double val() const override
Calculate value of cost function.
virtual double valAndDeriv(std::vector< double > &der) const override
Calculate the value and the derivatives of the cost function.
std::string parameterName(size_t i) const
Get name of i-th parameter.
bool isValid() const
Is the function set and valid?
CostFuncFitting()
Constructor.
virtual void calActiveCovarianceMatrix(EigenMatrix &covar, double epsrel=1e-8)
Calculates covariance matrix for fitting function's active parameters.
bool m_dirtyVal
dirty value flag
EigenVector m_pushedParams
double getParameter(size_t i) const override
Get i-th parameter.
virtual void addValDerivHessian(API::IFunction_sptr function, API::FunctionDomain_sptr domain, API::FunctionValues_sptr values, bool evalDeriv=true, bool evalHessian=true) const =0
Increments the cost function and its derivatives by evaluating them on a domain.
size_t m_numberFunParams
Number of all parameters in the fitting function.
void drop()
Discard saved parameters, derivatives and hessian.
void checkValidity() const
Throw a runtime_error if function is invalid.
const EigenVector & getDeriv() const
Return cached or calculate the drivatives.
bool m_dirtyDeriv
dirty derivatives flag
virtual void calCovarianceMatrix(EigenMatrix &covar, double epsrel=1e-8)
Calculates covariance matrix.
void applyTies()
Apply ties in the fitting function.
Two implementations of Jacobian.
map_type & getJ()
Get the map to the jacobian.
A wrapper around Eigen::Matrix.
bool isEmpty() const
Is matrix empty.
void set(size_t i, size_t j, double value)
Set an element.
double get(size_t i, size_t j) const
Get an element.
size_t size1() const
First size of the matrix.
EigenMatrix tr() const
Calculate the eigensystem of a symmetric matrix.
map_type & mutator()
Get the map to Eigen matrix.
void zero()
Set all elements to zero.
size_t size2() const
Second size of the matrix.
const map_type inspector() const
Get a const copy of the Eigen matrix.
void resize(const size_t nx, const size_t ny)
Resize the matrix.
A wrapper around Eigen::Vector.
void set(const size_t i, const double value)
Set an element.
double get(const size_t i) const
Get an element.
void resize(const size_t n)
Resize the vector.
size_t size() const
Size of the vector.
Exception thrown when a fitting function changes number of parameters during fit.
std::shared_ptr< FunctionValues > FunctionValues_sptr
typedef for a shared pointer
std::shared_ptr< IFunction > IFunction_sptr
shared pointer to the function base class
std::shared_ptr< FunctionDomain > FunctionDomain_sptr
typedef for a shared pointer
gsl_matrix_const_view const getGSLMatrixView_const(const map_type m)
take data from a constEigen Matrix and return a transposed gsl view.
gsl_matrix_view getGSLMatrixView(map_type &tr)
take data from an Eigen Matrix and return a transposed a gsl view.