Logo ROOT   6.12/07
Reference Guide
FumiliGradientCalculator.cxx
Go to the documentation of this file.
1 // @(#)root/minuit2:$Id$
2 // Authors: M. Winkler, F. James, L. Moneta, A. Zsenei 2003-2005
3 
4 /**********************************************************************
5  * *
6  * Copyright (c) 2005 LCG ROOT Math team, CERN/PH-SFT *
7  * *
8  **********************************************************************/
9 
11 #include "Minuit2/FumiliFCNBase.h"
15 #include "Minuit2/FumiliChi2FCN.h"
17 
18 //to compare with N2P calculator
19 //#define DEBUG 1
20 #ifdef DEBUG
21 #include "Minuit2/MnPrint.h"
23 #include "Minuit2/MnStrategy.h"
24 #include "Minuit2/MnUserFcn.h"
25 #endif
26 
27 namespace ROOT {
28 
29  namespace Minuit2 {
30 
31 
33 
34  // Calculate gradient for Fumili using the gradient and Hessian provided by the FCN Fumili function
35  // applying the external to int trasformation.
36 
37  int nvar = par.Vec().size();
38  std::vector<double> extParam = fTransformation( par.Vec() );
39  // std::vector<double> deriv;
40  // deriv.reserve( nvar );
41  // for (int i = 0; i < nvar; ++i) {
42  // unsigned int ext = fTransformation.ExtOfInt(i);
43  // if ( fTransformation.Parameter(ext).HasLimits())
44  // deriv.push_back( fTransformation.DInt2Ext( i, par.Vec()(i) ) );
45  // else
46  // deriv.push_back(1.0);
47  // }
48 
49  // eval Gradient
50  FumiliFCNBase & fcn = const_cast<FumiliFCNBase &>(fFcn);
51 
52  fcn.EvaluateAll(extParam);
53 
54 
55  MnAlgebraicVector v(nvar);
56  MnAlgebraicSymMatrix h(nvar);
57 
58 
59  const std::vector<double> & fcn_gradient = fFcn.Gradient();
60  assert( fcn_gradient.size() == extParam.size() );
61 
62 
63  // for (int i = 0; i < nvar; ++i) {
64  // unsigned int iext = fTransformation.ExtOfInt(i);
65  // double ideriv = 1.0;
66  // if ( fTransformation.Parameter(iext).HasLimits())
67  // ideriv = fTransformation.DInt2Ext( i, par.Vec()(i) ) ;
68 
69 
70  // // v(i) = fcn_gradient[iext]*deriv;
71  // v(i) = ideriv*fcn_gradient[iext];
72 
73  // for (int j = i; j < nvar; ++j) {
74  // unsigned int jext = fTransformation.ExtOfInt(j);
75  // double jderiv = 1.0;
76  // if ( fTransformation.Parameter(jext).HasLimits())
77  // jderiv = fTransformation.DInt2Ext( j, par.Vec()(j) ) ;
78 
79  // // h(i,j) = deriv[i]*deriv[j]*fFcn.Hessian(iext,jext);
80  // h(i,j) = ideriv*jderiv*fFcn.Hessian(iext,jext);
81  // }
82  // }
83 
84 
85  // cache deriv and Index values .
86  // in large Parameter limit then need to re-optimize and see if better not caching
87 
88  std::vector<double> deriv(nvar);
89  std::vector<unsigned int> extIndex(nvar);
90  for (int i = 0; i < nvar; ++i) {
91  extIndex[i] = fTransformation.ExtOfInt(i);
92  deriv[i] = 1;
93  if ( fTransformation.Parameter(extIndex[i]).HasLimits())
94  deriv[i] = fTransformation.DInt2Ext( i, par.Vec()(i) ) ;
95 
96  v(i) = fcn_gradient[extIndex[i]]*deriv[i];
97 
98  for (int j = 0; j <= i; ++j) {
99  h(i,j) = deriv[i]*deriv[j]*fFcn.Hessian(extIndex[i],extIndex[j]);
100  }
101  }
102 
103 #ifdef DEBUG
104  // compare with other Gradient
105  // // calculate Gradient using Minuit method
106 
108  FunctionGradient g2 = gc(par);
109 
110  std::cout << "Fumili Gradient " << v << std::endl;
111  std::cout << "Minuit Gradient " << g2.Vec() << std::endl;
112 #endif
113 
114  // store calculated Hessian
115  fHessian = h;
116  return FunctionGradient(v);
117 }
118 
120  const FunctionGradient&) const
121 
122 {
123  // Needed for interface of base class.
124  return this->operator()(par);
125 
126 }
127 
128  } // namespace Minuit2
129 
130 } // namespace ROOT
Namespace for new ROOT classes and functions.
Definition: StringConv.hxx:21
const MnAlgebraicVector & Vec() const
FunctionGradient operator()(const MinimumParameters &) const
TH1 * h
Definition: legend2.C:5
Class describing a symmetric matrix of size n.
Definition: LASymMatrix.h:51
double DInt2Ext(unsigned int, double) const
const MnAlgebraicVector & Vec() const
unsigned int ExtOfInt(unsigned int internal) const
virtual double Hessian(unsigned int row, unsigned int col) const
Return Value of the i-th j-th element of the Hessian matrix estimated previously using the FumiliFCNB...
class performing the numerical gradient calculation
SVector< double, 2 > v
Definition: Dict.h:5
unsigned int size() const
Definition: LAVector.h:198
virtual void EvaluateAll(const std::vector< double > &par)=0
Evaluate function Value, Gradient and Hessian using Fumili approximation, for values of parameters p ...
Wrapper used by Minuit of FCN interface containing a reference to the transformation object...
Definition: MnUserFcn.h:26
const MinuitParameter & Parameter(unsigned int) const
virtual const std::vector< double > & Gradient() const
Return cached Value of function Gradient estimated previously using the FumiliFCNBase::EvaluateAll me...
Extension of the FCNBase for the Fumili method.
Definition: FumiliFCNBase.h:47
API class for defining three levels of strategies: low (0), medium (1), high (>=2); acts on: Migrad (...
Definition: MnStrategy.h:27