Logo ROOT  
Reference Guide
FumiliFCNAdapter.h
Go to the documentation of this file.
1 // @(#)root/minuit2:$Id$
2 // Author: L. Moneta 10/2006
3 
4 /**********************************************************************
5  * *
6  * Copyright (c) 2006 ROOT Foundation, CERN/PH-SFT *
7  * *
8  **********************************************************************/
9 
10 #ifndef ROOT_Minuit2_FumiliFCNAdapter
11 #define ROOT_Minuit2_FumiliFCNAdapter
12 
13 #include "Minuit2/FumiliFCNBase.h"
14 
15 #include "Math/FitMethodFunction.h"
16 
17 #include "Minuit2/MnPrint.h"
18 
19 // #ifndef ROOT_Math_Util
20 // #include "Math/Util.h"
21 // #endif
22 
23 #include <cmath>
24 #include <cassert>
25 #include <vector>
26 
27 namespace ROOT {
28 
29 namespace Minuit2 {
30 
31 /**
32 
33 
34 template wrapped class for adapting to FumiliFCNBase signature
35 
36 @author Lorenzo Moneta
37 
38 @ingroup Minuit
39 
40 */
41 
42 template <class Function>
44 
45 public:
46  // typedef ROOT::Math::FitMethodFunction Function;
47  typedef typename Function::Type_t Type_t;
48 
49  FumiliFCNAdapter(const Function &f, unsigned int ndim, double up = 1.) : FumiliFCNBase(ndim), fFunc(f), fUp(up) {}
50 
52 
53  double operator()(const std::vector<double> &v) const { return fFunc.operator()(&v[0]); }
54  double operator()(const double *v) const { return fFunc.operator()(v); }
55  double Up() const { return fUp; }
56 
57  void SetErrorDef(double up) { fUp = up; }
58 
59  // virtual std::vector<double> Gradient(const std::vector<double>&) const;
60 
61  // forward interface
62  // virtual double operator()(int npar, double* params,int iflag = 4) const;
63 
64  /**
65  evaluate gradient hessian and function value needed by fumili
66  */
67  void EvaluateAll(const std::vector<double> &v);
68 
69 private:
70  // data member
71 
72  const Function &fFunc;
73  double fUp;
74 };
75 
76 template <class Function>
77 void FumiliFCNAdapter<Function>::EvaluateAll(const std::vector<double> &v)
78 {
79  MnPrint print("FumiliFCNAdaptor");
80 
81  // typedef FumiliFCNAdapter::Function Function;
82 
83  // evaluate all elements
84  unsigned int npar = Dimension();
85  if (npar != v.size())
86  print.Error("npar", npar, "v.size()", v.size());
87  assert(npar == v.size());
88  // must distinguish case of likelihood or LS
89 
90  std::vector<double> &grad = Gradient();
91  std::vector<double> &hess = Hessian();
92  // reset
93  assert(grad.size() == npar);
94  grad.assign(npar, 0.0);
95  hess.assign(hess.size(), 0.0);
96 
97  double sum = 0;
98  unsigned int ndata = fFunc.NPoints();
99 
100  std::vector<double> gf(npar);
101 
102  // loop on the data points
103 
104  // assume for now least-square
105  if (fFunc.Type() == Function::kLeastSquare) {
106 
107  for (unsigned int i = 0; i < ndata; ++i) {
108  // calculate data element and gradient
109  double fval = fFunc.DataElement(&v.front(), i, &gf[0]);
110 
111  // t.b.d should protect for bad values of fval
112  sum += fval * fval;
113 
114  for (unsigned int j = 0; j < npar; ++j) {
115  grad[j] += 2. * fval * gf[j];
116  for (unsigned int k = j; k < npar; ++k) {
117  int idx = j + k * (k + 1) / 2;
118  hess[idx] += 2.0 * gf[j] * gf[k];
119  }
120  }
121  }
122  } else if (fFunc.Type() == Function::kLogLikelihood) {
123 
124  for (unsigned int i = 0; i < ndata; ++i) {
125 
126  // calculate data element and gradient
127  // return value is log of pdf and derivative of the log(Pdf)
128  double fval = fFunc.DataElement(&v.front(), i, &gf[0]);
129 
130  sum -= fval;
131 
132  for (unsigned int j = 0; j < npar; ++j) {
133  double gfj = gf[j];
134  grad[j] -= gfj;
135  for (unsigned int k = j; k < npar; ++k) {
136  int idx = j + k * (k + 1) / 2;
137  hess[idx] += gfj * gf[k];
138  }
139  }
140  }
141  } else {
142  print.Error("Type of fit method is not supported, it must be chi2 or log-likelihood");
143  }
144 }
145 
146 } // end namespace Minuit2
147 
148 } // end namespace ROOT
149 
150 #endif // ROOT_Minuit2_FCNAdapter
f
#define f(i)
Definition: RSha256.hxx:104
ROOT::Minuit2::FumiliFCNAdapter::EvaluateAll
void EvaluateAll(const std::vector< double > &v)
evaluate gradient hessian and function value needed by fumili
Definition: FumiliFCNAdapter.h:77
ROOT::Minuit2::FumiliFCNBase
Extension of the FCNBase for the Fumili method.
Definition: FumiliFCNBase.h:46
ROOT::Minuit2::FumiliFCNAdapter
template wrapped class for adapting to FumiliFCNBase signature
Definition: FumiliFCNAdapter.h:43
ROOT::Minuit2::FumiliFCNAdapter::Up
double Up() const
Error definition of the function.
Definition: FumiliFCNAdapter.h:55
sum
static uint64_t sum(uint64_t i)
Definition: Factory.cxx:2345
Function
Double_t(* Function)(Double_t)
Definition: Functor.C:4
FitMethodFunction.h
ROOT::Minuit2::MnPrint::Error
void Error(const Ts &... args)
Definition: MnPrint.h:120
v
@ v
Definition: rootcling_impl.cxx:3635
ROOT::Minuit2::FumiliFCNAdapter::operator()
double operator()(const std::vector< double > &v) const
The meaning of the vector of parameters is of course defined by the user, who uses the values of thos...
Definition: FumiliFCNAdapter.h:53
ROOT::Minuit2::FumiliFCNAdapter::operator()
double operator()(const double *v) const
Definition: FumiliFCNAdapter.h:54
ROOT::Minuit2::FumiliFCNAdapter::SetErrorDef
void SetErrorDef(double up)
add interface to set dynamically a new error definition Re-implement this function if needed.
Definition: FumiliFCNAdapter.h:57
ROOT::Minuit2::FumiliFCNAdapter::fFunc
const Function & fFunc
Definition: FumiliFCNAdapter.h:72
ROOT::Minuit2::FumiliFCNAdapter::~FumiliFCNAdapter
~FumiliFCNAdapter()
Definition: FumiliFCNAdapter.h:51
ROOT::Minuit2::FumiliFCNAdapter::fUp
double fUp
Definition: FumiliFCNAdapter.h:73
FumiliFCNBase.h
MnPrint.h
ROOT
tbb::task_arena is an alias of tbb::interface7::task_arena, which doesn't allow to forward declare tb...
Definition: EExecutionPolicy.hxx:4
ROOT::Minuit2::FumiliFCNAdapter::Type_t
Function::Type_t Type_t
Definition: FumiliFCNAdapter.h:47
ROOT::Minuit2::FumiliFCNAdapter::FumiliFCNAdapter
FumiliFCNAdapter(const Function &f, unsigned int ndim, double up=1.)
Definition: FumiliFCNAdapter.h:49
ROOT::Minuit2::MnPrint
Definition: MnPrint.h:73