SourceXtractorPlusPlus  0.8
Please provide a description of the project.
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Example_SimpleFit.cpp
Go to the documentation of this file.
1 
7 #include <iostream>
8 #include <vector>
12 #include <ModelFitting/utils.h>
13 #include <chrono>
15 #include "utils.h"
16 
17 using namespace ModelFitting;
18 
20 private:
23 
24 public:
27  : m_A{a}, m_lambda{lambda}, m_b{b}, m_y{y}, m_t{t} {
28  }
29 
30  virtual ~ExpResidualProvider() = default;
31 
32  std::size_t numberOfResiduals() const override {
33  return m_y.size();
34  }
35 
36  void populateResidualBlock(IterType iter) override {
37  std::cout << std::fixed << std::setprecision(4) << "A=" << m_A->getValue() << "\tlambda=" << m_lambda->getValue() << "\tb=" << m_b->getValue();
38  double e = 0;
39  for (size_t i = 0; i < m_y.size(); ++i) {
40  double Yi = m_b->getValue() + m_A->getValue() * exp(-m_lambda->getValue() * m_t[i]);
41  *iter = Yi - m_y[i];
42  e += *iter * *iter;
43  ++iter;
44  }
45  std::cout << std::scientific << "\tsum(e^2)=" << e << std::endl;
46  }
47 };
48 
49 /*
50  * Adapted from
51  * https://www.gnu.org/software/gsl/doc/html/nls.html#exponential-fitting-example
52  * "Simple" exponential function, so it is easier to understand and visualize what's going on
53  */
54 int main(int argc, char **argv) {
55  std::string engine_impl("levmar");
56  if (argc > 1) {
57  engine_impl = argv[1];
58  }
59 
60  // Generate data points, using the exponential function
61  int n = 100;
62  double tmax = 40.;
63  std::vector<double> y(n), t(n);
64 
65  for (int i = 0; i < n; ++i) {
66  double ti = i * tmax / (n - 1.0);
67  double yi = 1. + 5. * exp(-0.1 * ti);
68  t[i] = ti;
69  y[i] = yi;
70  }
71 
72  // Initialize and run engine
73  auto A = std::make_shared<EngineParameter>(2., make_unique<NeutralConverter>());
74  auto lambda = std::make_shared<EngineParameter>(1., make_unique<NeutralConverter>());
75  auto b = std::make_shared<EngineParameter>(0., make_unique<NeutralConverter>());
76 
77  EngineParameterManager manager {};
78  manager.registerParameter(A);
79  manager.registerParameter(lambda);
80  manager.registerParameter(b);
81 
82  ResidualEstimator res_estimator {};
83  res_estimator.registerBlockProvider(make_unique<ExpResidualProvider>(A, lambda, b, y, t));
84 
85  std::cout << "Using engine " << engine_impl << std::endl;
86  auto engine = LeastSquareEngineManager::create(engine_impl);
88  auto solution = engine->solveProblem(manager, res_estimator);
90 
92  std::cout << "Time of fitting: " << std::chrono::duration<double, std::milli>(t2 - t1).count() << " ms" << std::endl;
93  std::cout << "A (5.) : " << A->getValue() << std::endl;
94  std::cout << "lambda (0.1): " << lambda->getValue() << std::endl;
95  std::cout << "b (1.) : " << b->getValue() << std::endl;
96  printLevmarInfo(boost::any_cast<std::array<double, 10>>(solution.underlying_framework_info));
97 }
T exp(T...args)
int main()
Definition: Example1.cpp:47
constexpr double e
T endl(T...args)
ExpResidualProvider(std::shared_ptr< BasicParameter > a, std::shared_ptr< BasicParameter > lambda, std::shared_ptr< BasicParameter > b, const std::vector< double > &y, const std::vector< double > &t)
std::shared_ptr< DependentParameter< std::shared_ptr< EngineParameter > > > y
STL class.
Interface of a class which can provide a block of residuals for least square minimization solving...
void registerBlockProvider(std::unique_ptr< ResidualBlockProvider > provider)
Registers a ResidualBlockProvider to the ResidualEstimator.
std::size_t numberOfResiduals() const override
Returns the number of residuals provided by this provider.
void registerParameter(std::shared_ptr< EngineParameter > parameter)
Registers an EngineParameter to the EngineParameterManager.
std::vector< double > m_y
T count(T...args)
T fixed(T...args)
void populateResidualBlock(IterType iter) override
Provides the residual values.
STL class.
Class responsible for managing the parameters the least square engine minimizes.
static std::shared_ptr< LeastSquareEngine > create(const std::string &name, unsigned max_iterations=1000)
void printLevmarInfo(std::array< double, 10 > info)
Definition: utils.h:118
std::shared_ptr< BasicParameter > m_lambda
Provides to the LeastSquareEngine the residual values.
T setprecision(T...args)