ROOT  6.06/09
Reference Guide
timeGA.cxx
Go to the documentation of this file.
1 #include "TH1.h"
2 #include "TF1.h"
3 #include "Fit/BinData.h"
4 #include "Fit/Chi2FCN.h"
5 
6 #include "Math/WrappedMultiTF1.h"
7 #include "Math/Minimizer.h"
9 #include "Math/Factory.h"
10 #include "HFitInterface.h"
11 
12 #include "TMath.h"
13 
14 #include <TStopwatch.h>
15 
16 //#define DEBUG
17 
18 // Quadratic background function
20  return par[0] + par[1]*x[0];
21 }
22 
23 // Gaussian Peak function
25  return par[0]*TMath::Gaus(x[0],par[1], par[2]);
26 }
27 
28 // Sum of background and peak function
30  return background(x,par) + gaussianPeak(x,&par[2]) + gaussianPeak(x,&par[5]);
31 }
32 
33 // We'll look for the minimum at 2 and 7
34 double par0[8] = { 1, 0.05, 10 , 2, 0.5 , 10 , 7 , 1. };
35 const int ndata = 10000;
36 const double gAbsTolerance = 0.1;
37 int gVerbose = 0;
38 
39 using std::cout;
40 using std::endl;
41 
42 int GAMinimize(ROOT::Math::IMultiGenFunction& chi2Func, double& xm1, double& xm2)
43 {
44  // minimize the function
46  if (min == 0) {
47  std::cout << "Error creating minimizer " << std::endl;
48  return -1;
49  }
50 
51  // Set the parameters for the minimization.
52  min->SetMaxFunctionCalls(1000000);
53  min->SetMaxIterations(100000);
55  min->SetPrintLevel(gVerbose);
56  min->SetFunction(chi2Func);
57  ROOT::Math::GeneticMinimizerParameters params; // construct with default values
58  params.fNsteps = 100; // increset number of steps top 100 (default is 40)
59  min->SetParameters(params);
60 
61 
62  // initial values of the function
63  double x0[8];
64  std::copy(par0, par0 + 8, x0);
65  x0[3] = xm1;
66  x0[6] = xm2;
67 
68  for (unsigned int i = 0; i < chi2Func.NDim(); ++i) {
69 #ifdef DEBUG
70  std::cout << "set variable " << i << " to value " << x0[i] << std::endl;
71 #endif
72  if ( i == 3 || i == 6 )
73  min->SetLimitedVariable(i,"x" + ROOT::Math::Util::ToString(i),x0[i], 0.1,2,8);
74  else
75  min->SetLimitedVariable(i,"x" + ROOT::Math::Util::ToString(i),x0[i], 0.1,x0[i]-2,x0[i]+2);
76  }
77 
78 
79  // minimize
80  if ( !min->Minimize() ) return 1;
81  min->MinValue();
82 
83  // show the results
84  std::cout << "Min values by GeneticMinimizer: " << min->X()[3] << " " << min->X()[6] << std::endl;
85  xm1 = min->X()[3];
86  xm2 = min->X()[6];
87 
88  return 0;
89 }
90 
91 
93 {
94  double x1=0, x2=0;
95 
96  // create a TF1 from fit function and generate histogram
97  TF1 * fitFunc = new TF1("fitFunc",fitFunction,0,10,8);
98  fitFunc->SetParameters(par0);
99 
100  // Create a histogram filled with random data from TF1
101  TH1D * h1 = new TH1D("h1","h1",100, 0, 10);
102  for (int i = 0; i < ndata; ++i) {
103  h1->Fill(fitFunc->GetRandom() );
104  }
105 
106  // perform the fit
108  ROOT::Fit::FillData(d,h1);
109  ROOT::Math::WrappedMultiTF1 f(*fitFunc);
110 
111  // Create the function for fitting.
113 
114  // Look for an approximation with a Genetic Algorithm
115  TStopwatch t;
116  t.Start();
117  GAMinimize(chi2Func, x1,x2);
118  t.Stop();
119  std::cout << "Time :\t " << t.RealTime() << " " << t.CpuTime() << std::endl;
120 
121  return 0;
122 }
123 
124 int main(int argc, char **argv)
125 {
126  // Parse command line arguments
127  for (Int_t i=1 ; i<argc ; i++) {
128  std::string arg = argv[i] ;
129  if (arg == "-v") {
130  gVerbose = 1;
131  }
132  if (arg == "-vv") {
133  gVerbose = 3;
134  }
135  if (arg == "-h") {
136  std::cout << "Usage: " << argv[0] << " [-v] [-vv]\n";
137  std::cout << " where:\n";
138  std::cout << " -v : verbose mode\n";
139  std::cout << " -vv : very verbose mode\n";
140  std::cout << std::endl;
141  return -1;
142  }
143  }
144 
145  int status = GAMinTutorial();
146 
147  return status;
148 }
double par[1]
Definition: unuranDistr.cxx:38
Double_t fitFunction(Double_t *x, Double_t *par)
Definition: timeGA.cxx:29
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
Definition: TH1.cxx:3159
virtual void SetParameters(const Double_t *params)
Definition: TF1.h:439
virtual Double_t GetRandom()
Return a random number following this function shape.
Definition: TF1.cxx:1785
void SetMaxIterations(unsigned int maxiter)
set maximum iterations (one iteration can have many function calls)
Definition: Minimizer.h:459
virtual bool Minimize()
method to perform the minimization
static Vc_ALWAYS_INLINE int_v min(const int_v &x, const int_v &y)
Definition: vector.h:433
Double_t RealTime()
Stop the stopwatch (if it is running) and return the realtime (in seconds) passed between the start a...
Definition: TStopwatch.cxx:108
void Start(Bool_t reset=kTRUE)
Start the stopwatch.
Definition: TStopwatch.cxx:56
Class to Wrap a ROOT Function class (like TF1) in a IParamMultiFunction interface of multi-dimensions...
Double_t background(Double_t *x, Double_t *par)
Definition: timeGA.cxx:19
Double_t gaussianPeak(Double_t *x, Double_t *par)
Definition: timeGA.cxx:24
double par0[8]
Definition: timeGA.cxx:34
double fitFunc(double *x, double *p)
Double_t CpuTime()
Stop the stopwatch (if it is running) and return the cputime (in seconds) passed between the start an...
Definition: TStopwatch.cxx:123
TAlienJobStatus * status
Definition: TAlienJob.cxx:51
int Int_t
Definition: RtypesCore.h:41
THist< 1, double > TH1D
Definition: THist.h:314
static const double x2[5]
void Stop()
Stop the stopwatch.
Definition: TStopwatch.cxx:75
Double_t x[n]
Definition: legend1.C:17
TH1F * h1
Definition: legend1.C:5
int gVerbose
Definition: timeGA.cxx:37
Chi2FCN class for binnned fits using the least square methods.
Definition: Chi2FCN.h:66
int main(int argc, char **argv)
Definition: timeGA.cxx:124
void FillData(BinData &dv, const TH1 *hist, TF1 *func=0)
fill the data vector from a TH1.
virtual unsigned int NDim() const =0
Retrieve the dimension of the function.
const int ndata
Definition: timeGA.cxx:35
Class describing the binned data sets : vectors of x coordinates, y values and optionally error on y ...
Definition: BinData.h:61
void SetParameters(const GeneticMinimizerParameters &params)
1-D histogram with a double per channel (see TH1 documentation)}
Definition: TH1.h:613
int GAMinimize(ROOT::Math::IMultiGenFunction &chi2Func, double &xm1, double &xm2)
Definition: timeGA.cxx:42
Double_t Gaus(Double_t x, Double_t mean=0, Double_t sigma=1, Bool_t norm=kFALSE)
Calculate a gaussian function with mean and sigma.
Definition: TMath.cxx:453
void SetMaxFunctionCalls(unsigned int maxfcn)
set maximum of function calls
Definition: Minimizer.h:456
static const double x1[5]
double f(double x)
virtual const double * X() const
return pointer to X values at the minimum
double Double_t
Definition: RtypesCore.h:55
void SetTolerance(double tol)
set the tolerance
Definition: Minimizer.h:462
const double gAbsTolerance
Definition: timeGA.cxx:36
virtual void SetFunction(const ROOT::Math::IMultiGenFunction &func)
set the function to minimize
int GAMinTutorial()
Definition: timeGA.cxx:92
virtual double MinValue() const
return minimum function value
std::string ToString(const T &val)
Utility function for conversion to strings.
Definition: Util.h:42
1-Dim function class
Definition: TF1.h:149
void SetPrintLevel(int level)
set print level
Definition: Minimizer.h:453
Documentation for the abstract class IBaseFunctionMultiDim.
Definition: IFunction.h:63
virtual bool SetLimitedVariable(unsigned int, const std::string &, double, double, double, double)
set a new upper/lower limited variable (override if minimizer supports them ) otherwise as default se...
Stopwatch class.
Definition: TStopwatch.h:30