Logo ROOT   6.10/09
Reference Guide
testFitPerf.cxx
Go to the documentation of this file.
1 #include "TH1.h"
2 #include "TF1.h"
3 #include "TF2.h"
4 #include "TMath.h"
5 #include "TSystem.h"
6 #include "TRandom3.h"
7 #include "TTree.h"
8 #include "TROOT.h"
9 
10 #include "Fit/BinData.h"
11 #include "Fit/UnBinData.h"
12 //#include "Fit/BinPoint.h"
13 #include "Fit/Fitter.h"
14 #include "HFitInterface.h"
15 
16 #include "Math/IParamFunction.h"
17 #include "Math/WrappedTF1.h"
18 #include "Math/WrappedMultiTF1.h"
21 
22 #include "TGraphErrors.h"
23 
24 #include "TStyle.h"
25 
26 #include "TSeqCollection.h"
27 
28 #include "Math/Polynomial.h"
29 #include "Math/DistFunc.h"
30 
31 #include <string>
32 #include <iostream>
33 
34 #include "TStopwatch.h"
35 
36 #include "TVirtualFitter.h"
37 //#include "TFitterMinuit.h"
38 // #include "TFitterFumili.h"
39 // #include "TFumili.h"
40 
41 #include "GaussFunction.h"
42 
43 #include "RooDataHist.h"
44 #include "RooDataSet.h"
45 #include "RooRealVar.h"
46 #include "RooGaussian.h"
47 #include "RooMinuit.h"
48 #include "RooChi2Var.h"
49 #include "RooGlobalFunc.h"
50 #include "RooFitResult.h"
51 #include "RooProdPdf.h"
52 
53 #include <cassert>
54 
55 #include "MinimizerTypes.h"
56 
57 //#define DEBUG
58 
59 int nfit;
60 const int N = 20;
61 double iniPar[2*N];
62 
63 
65  for (unsigned int i = 0; i < data.Size(); ++i) {
66  std::cout << data.Coords(i)[0] << "\t";
67  }
68  std::cout << "\ndata size is " << data.Size() << std::endl;
69 }
70 
71 void printResult(int iret) {
72  std::cout << "\n************************************************************\n";
73  std::cout << "Test\t\t\t\t";
74  if (iret == 0) std::cout << "OK";
75  else std::cout << "FAILED";
76  std::cout << "\n************************************************************\n";
77 }
78 
79 bool USE_BRANCH = false;
80 ROOT::Fit::UnBinData * FillUnBinData(TTree * tree, bool copyData = true, unsigned int dim = 1 ) {
81 
82  // fill the unbin data set from a TTree
83  ROOT::Fit::UnBinData * d = 0;
84  // for the large tree
85  if (std::string(tree->GetName()) == "t2") {
86  d = new ROOT::Fit::UnBinData();
87  // large tree
88  unsigned int n = tree->GetEntries();
89 #ifdef DEBUG
90  std::cout << "number of unbin data is " << n << " of dim " << N << std::endl;
91 #endif
92  d->Append(n,N);
93  TBranch * bx = tree->GetBranch("x");
94  double vx[N];
95  bx->SetAddress(vx);
96  std::vector<double> m(N);
97  for (int unsigned i = 0; i < n; ++i) {
98  bx->GetEntry(i);
99  d->Add(vx);
100  for (int j = 0; j < N; ++j)
101  m[j] += vx[j];
102  }
103 
104 #ifdef DEBUG
105  std::cout << "average values of means :\n";
106  for (int j = 0; j < N; ++j)
107  std::cout << m[j]/n << " ";
108  std::cout << "\n";
109 #endif
110 
111  return d;
112  }
113  if (USE_BRANCH)
114  {
115  d = new ROOT::Fit::UnBinData();
116  unsigned int n = tree->GetEntries();
117  //std::cout << "number of unbin data is " << n << std::endl;
118 
119  if (dim == 2) {
120  d->Append(n,2);
121  TBranch * bx = tree->GetBranch("x");
122  TBranch * by = tree->GetBranch("y");
123  double v[2];
124  bx->SetAddress(&v[0]);
125  by->SetAddress(&v[1]);
126  for (int unsigned i = 0; i < n; ++i) {
127  bx->GetEntry(i);
128  by->GetEntry(i);
129  d->Add(v);
130  }
131  }
132  else if (dim == 1) {
133  d->Append(n,1);
134  TBranch * bx = tree->GetBranch("x");
135  double v[1];
136  bx->SetAddress(&v[0]);
137  for (int unsigned i = 0; i < n; ++i) {
138  bx->GetEntry(i);
139  d->Add(v);
140  }
141  }
142 
143  return d;
144 
145  //printData(d);
146  }
147  else {
148  tree->SetEstimate(tree->GetEntries());
149 
150  // use TTREE::Draw
151  if (dim == 2) {
152  tree->Draw("x:y",0,"goff"); // goff is used to turn off the graphics
153  double * x = tree->GetV1();
154  double * y = tree->GetV2();
155 
156  if (x == 0 || y == 0) {
157  USE_BRANCH= true;
158  return FillUnBinData(tree, true, dim);
159  }
160 
161  // use array pre-allocated in tree->Draw . This is faster
162  //assert(x != 0);
163  unsigned int n = tree->GetSelectedRows();
164 
165  if (copyData) {
166  d = new ROOT::Fit::UnBinData(n,2);
167  double vx[2];
168  for (int unsigned i = 0; i < n; ++i) {
169  vx[0] = x[i];
170  vx[1] = y[i];
171  d->Add(vx);
172  }
173  }
174  else // use data pointers directly
175  d = new ROOT::Fit::UnBinData(n,x,y);
176 
177  }
178  else if ( dim == 1) {
179 
180  tree->Draw("x",0,"goff"); // goff is used to turn off the graphics
181  double * x = tree->GetV1();
182 
183  if (x == 0) {
184  USE_BRANCH= true;
185  return FillUnBinData(tree, true, dim);
186  }
187  unsigned int n = tree->GetSelectedRows();
188 
189  if (copyData) {
190  d = new ROOT::Fit::UnBinData(n,1);
191  for (int unsigned i = 0; i < n; ++i) {
192  d->Add(x[i]);
193  }
194  }
195  else
196  d = new ROOT::Fit::UnBinData(n,x);
197  }
198  return d;
199  }
200 
201  //std::copy(x,x+n, d.begin() );
202  return 0;
203 }
204 
205 
206 
207 
208 // print the data
209 template <class T>
210 void printData(const T & data) {
211  for (typename T::const_iterator itr = data.begin(); itr != data.end(); ++itr) {
212  std::cout << itr->Coords()[0] << " " << itr->Value() << " " << itr->Error() << std::endl;
213  }
214  std::cout << "\ndata size is " << data.Size() << std::endl;
215 }
216 
217 
218 
219 // fitting using new fitter
221 template <class MinType, class T>
222 int DoBinFit(T * hist, Func & func, bool debug = false, bool useGrad = false) {
223 
224  //std::cout << "Fit histogram " << std::endl;
225 
227  ROOT::Fit::FillData(d,hist);
228 
229  //printData(d);
230 
231  // create the fitter
232 
233  ROOT::Fit::Fitter fitter;
234  fitter.Config().SetMinimizer(MinType::name().c_str(),MinType::name2().c_str());
235 
236  if (debug)
237  fitter.Config().MinimizerOptions().SetPrintLevel(3);
238 
239 
240  // create the function
241  if (!useGrad) {
242 
243  // use simply TF1 wrapper
244  //ROOT::Math::WrappedMultiTF1 f(*func);
245  //ROOT::Math::WrappedTF1 f(*func);
246  fitter.SetFunction(func);
247 
248  } else { // only for gaus fits
249  // use function gradient
250 #ifdef USE_MATHMORE_FUNC
251  // use mathmore for polynomial
252  ROOT::Math::Polynomial pol(2);
253  assert(pol.NPar() == func->GetNpar());
254  pol.SetParameters(func->GetParameters() );
255  ROOT::Math::WrappedParamFunction<ROOT::Math::Polynomial> f(pol,1,func->GetParameters(),func->GetParameters()+func->GetNpar() );
256 #endif
258  f.SetParameters(func.Parameters());
259  fitter.SetFunction(f);
260  }
261 
262 
263  bool ret = fitter.Fit(d);
264  if (!ret) {
265  std::cout << " Fit Failed " << std::endl;
266  return -1;
267  }
268  if (debug)
269  fitter.Result().Print(std::cout);
270  return 0;
271 }
272 
273 // unbin fit
274 template <class MinType, class T>
275 int DoUnBinFit(T * tree, Func & func, bool debug = false, bool copyData = false ) {
276 
277  ROOT::Fit::UnBinData * d = FillUnBinData(tree, copyData, func.NDim() );
278  // need to have done Tree->Draw() before fit
279  //FillUnBinData(d,tree);
280 
281  //std::cout << "data size type and size is " << typeid(*d).name() << " " << d->Size() << std::endl;
282  if (copyData)
283  std::cout << "\tcopy data in FitData\n";
284  else
285  std::cout << "\tre-use original data \n";
286 
287 
288 
289  //printData(d);
290 
291  // create the fitter
292  //std::cout << "Fit parameter 2 " << f.Parameters()[2] << std::endl;
293 
294  ROOT::Fit::Fitter fitter;
295  fitter.Config().SetMinimizer(MinType::name().c_str(),MinType::name2().c_str());
296 
297  if (debug)
298  fitter.Config().MinimizerOptions().SetPrintLevel(3);
299 
300  // set tolerance 1 for tree to be same as in TTTreePlayer::UnBinFIt
301  fitter.Config().MinimizerOptions().SetTolerance(1);
302 
303 
304  // create the function
305 
306  fitter.SetFunction(func);
307  // need to fix param 0 , normalization in the unbinned fits
308  //fitter.Config().ParSettings(0).Fix();
309 
310  bool ret = fitter.Fit(*d);
311  if (!ret) {
312  std::cout << " Fit Failed " << std::endl;
313  return -1;
314  }
315  if (debug)
316  fitter.Result().Print(std::cout);
317 
318  delete d;
319 
320  return 0;
321 
322 }
323 
324 template <class MinType>
325 int DoFit(TTree * tree, Func & func, bool debug = false, bool copyData = false ) {
326  return DoUnBinFit<MinType, TTree>(tree, func, debug, copyData);
327 }
328 template <class MinType>
329 int DoFit(TH1 * h1, Func & func, bool debug = false, bool copyData = false ) {
330  return DoBinFit<MinType, TH1>(h1, func, debug, copyData);
331 }
332 template <class MinType>
333 int DoFit(TGraph * gr, Func & func, bool debug = false, bool copyData = false ) {
334  return DoBinFit<MinType, TGraph>(gr, func, debug, copyData);
335 }
336 
337 template <class MinType, class FitObj>
338 int FitUsingNewFitter(FitObj * fitobj, Func & func, bool useGrad=false) {
339 
340  std::cout << "\n************************************************************\n";
341  std::cout << "\tFit using new Fit::Fitter " << typeid(*fitobj).name() << std::endl;
342  std::cout << "\tMinimizer is " << MinType::name() << " " << MinType::name2() << " func dim = " << func.NDim() << std::endl;
343 
344  int iret = 0;
345  TStopwatch w; w.Start();
346 
347 #ifdef DEBUG
348  // std::cout << "initial Parameters " << iniPar << " " << *iniPar << " " << *(iniPar+1) << std::endl;
349  func.SetParameters(iniPar);
350  iret |= DoFit<MinType>(fitobj,func,true, useGrad);
351  if (iret != 0) {
352  std::cout << "Fit failed " << std::endl;
353  }
354 
355 #else
356  for (int i = 0; i < nfit; ++i) {
357  func.SetParameters(iniPar);
358  iret = DoFit<MinType>(fitobj,func, false, useGrad);
359  if (iret != 0) {
360  std::cout << "Fit failed " << std::endl;
361  break;
362  }
363  }
364 #endif
365  w.Stop();
366  std::cout << "\nTime: \t" << w.RealTime() << " , " << w.CpuTime() << std::endl;
367  std::cout << "\n************************************************************\n";
368 
369  return iret;
370 }
371 
372 
373 //------------old fit methods
374 
375 // fit using Fit method
376 template <class T, class MinType>
377 int FitUsingTFit(T * hist, TF1 * func) {
378 
379  std::cout << "\n************************************************************\n";
380  std::cout << "\tFit using " << hist->ClassName() << "::Fit\n";
381  std::cout << "\tMinimizer is " << MinType::name() << std::endl;
382 
383  std::string opt = "BFQ0";
384  if (MinType::name() == "Linear")
385  opt = "Q0";
386 
387  // check gminuit
388 // TSeqCollection * l = gROOT->GetListOfSpecials();
389 // for (int i = 0; i < l->GetEntries(); ++i) {
390 // TObject * obj = l->At(i);
391 // std::cout << " entries for i " << i << " of " << l->GetEntries() << " is " << obj << std::endl;
392 // if (obj != 0) std::cout << " object name is " << obj->GetName() << std::endl;
393 // }
394 
395  int iret = 0;
397 
398  TStopwatch w; w.Start();
399  for (int i = 0; i < nfit; ++i) {
400  func->SetParameters(iniPar);
401  iret |= hist->Fit(func,opt.c_str());
402  if (iret != 0) {
403  std::cout << "Fit failed " << std::endl;
404  return iret;
405  }
406  }
407  // std::cout << "iret " << iret << std::endl;
408 #ifdef DEBUG
409  func->SetParameters(iniPar);
410 // if (fitter == "Minuit2") {
411 // // increase print level
412 // TVirtualFitter * tvf = TVirtualFitter::GetFitter();
413 // TFitterMinuit * minuit2 = dynamic_cast<TFitterMinuit * >(tvf);
414 // assert (minuit2 != 0);
415 // minuit2->SetPrintLevel(3);
416 // }
417  if (MinType::name() != "Linear")
418  iret |= hist->Fit(func,"BFV0");
419  else
420  iret |= hist->Fit(func,"V0");
421  // get precice value of minimum
422  int pr = std::cout.precision(18);
423  std::cout << "Chi2 value = " << func->GetChisquare() << std::endl;
424  std::cout.precision(pr);
425 
426 #endif
427  w.Stop();
428  std::cout << "\nTime: \t" << w.RealTime() << " , " << w.CpuTime() << std::endl;
429  std::cout << "\n************************************************************\n";
430 
431  return iret;
432 }
433 
434 
435 // unbinned fit using Fit of trees
436 template <class MinType>
437 int FitUsingTTreeFit(TTree * tree, TF1 * func, const std::string & vars = "x") {
438 
439  std::cout << "\n************************************************************\n";
440  std::cout << "\tFit using TTree::UnbinnedFit\n";
441  std::cout << "\tMinimizer is " << MinType::name() << std::endl;
442 
443 
444  std::string sel = "";
445 
446  int iret = 0;
448 
449  TStopwatch w; w.Start();
450  for (int i = 0; i < nfit; ++i) {
451  func->SetParameters(iniPar);
452  iret |= tree->UnbinnedFit(func->GetName(),vars.c_str(),sel.c_str(),"Q");
453  if (iret != 0) {
454  std::cout << "Fit failed : iret = " << iret << std::endl;
455  return iret;
456  }
457  }
458  // std::cout << "iret " << iret << std::endl;
459 #ifdef DEBUG
460  func->SetParameters(iniPar);
461 
462  iret |= tree->UnbinnedFit(func->GetName(),vars.c_str(),sel.c_str(),"V");
463 
464 #endif
465  w.Stop();
466  std::cout << "\nTime: \t" << w.RealTime() << " , " << w.CpuTime() << std::endl;
467  std::cout << "\n************************************************************\n";
468 
469  return iret;
470 }
471 
472 
473 // int FitUsingTVirtualFitter(TH1 * hist, TF1 * func, const std::string & fitter) {
474 
475 // TVirtualFitter::SetDefaultFitter(fitter.c_str());
476 // std::cout << "Fit using TVirtualFitter and " << fitter << "\t Time: \t" << w.RealTime() << " , " << w.CpuTime() < std::endl;
477 // }
478 
479 
480 // ROoFit
481 
482 
483 //binned roo fit
484 int FitUsingRooFit(TH1 * hist, TF1 * func) {
485 
486  int iret = 0;
487  std::cout << "\n************************************************************\n";
488  std::cout << "\tFit using RooFit (Chi2 Fit)\n";
489  std::cout << "\twith function " << func->GetName() << "\n";
490 
491 
492  // start counting the time
493  TStopwatch w; w.Start();
494 
495  for (int i = 0; i < nfit; ++i) {
496 
497  RooRealVar x("x","x",-5,5) ;
498 
499  RooDataHist data("bindata","bin dataset with x",x,hist) ;
500  // define params
501  RooAbsPdf * pdf = 0;
502  RooRealVar * mean = 0;
503  RooRealVar * sigma = 0;
504 
505  func->SetParameters(iniPar);
506  std::string fname = func->GetName();
507  if (fname == "gaussian") {
508  double val,pmin,pmax;
509  val = func->GetParameter(1); //func->GetParLimits(1,-100,100);
510  mean = new RooRealVar("mean","Mean of Gaussian",val) ;
511  val = func->GetParameter(2); func->GetParLimits(1,pmin,pmax);
512  sigma = new RooRealVar("sigma","Width of Gaussian",val) ;
513 
514  pdf = new RooGaussian("gauss","gauss(x,mean,sigma)",x,*mean,*sigma) ;
515  }
516 
517  assert(pdf != 0);
518 #define USE_CHI2_FIT
519 #ifdef USE_CHI2_FIT
520  RooChi2Var chi2("chi2","chi2",*pdf,data) ;
521  RooMinuit m(chi2) ;
522  m.setPrintLevel(-1);
523  m.fit("mh") ;
524 #else
525  pdf->fitTo(data);
526 #endif
527 // if (iret != 0) return iret;
528  delete pdf;
529  delete mean; delete sigma;
530  }
531 
532  w.Stop();
533  std::cout << "\nTime: \t" << w.RealTime() << " , " << w.CpuTime() << std::endl;
534  std::cout << "\n************************************************************\n";
535  return iret;
536 }
537 
538 //unbinned roo fit
540 
541  int iret = 0;
542  std::cout << "\n************************************************************\n";
543  std::cout << "\tFit using RooFit (Likelihood Fit)\n";
544  std::cout << "\twith function " << func->GetName() << "\n";
545 
546 
547  // start counting the time
548  TStopwatch w; w.Start();
549 
550  for (int i = 0; i < nfit; ++i) {
551 
552  RooRealVar x("x","x",-100,100) ;
553  RooRealVar y("y","y",-100,100);
554 
555  RooDataSet data("unbindata","unbin dataset with x",tree,RooArgSet(x,y)) ;
556 
557 
558  RooRealVar mean("mean","Mean of Gaussian",iniPar[0], -100,100) ;
559  RooRealVar sigma("sigma","Width of Gaussian",iniPar[1], -100, 100) ;
560 
561  RooGaussian pdfx("gaussx","gauss(x,mean,sigma)",x,mean,sigma);
562 
563  // for 2d data
564  RooRealVar meany("meanx","Mean of Gaussian",iniPar[2], -100,100) ;
565  RooRealVar sigmay("sigmay","Width of Gaussian",iniPar[3], -100, 100) ;
566  RooGaussian pdfy("gaussy","gauss(y,meanx,sigmay)",y,meany,sigmay);
567 
568  RooProdPdf pdf("gausxy","gausxy",RooArgSet(pdfx,pdfy) );
569 
570 
571 #ifdef DEBUG
572  int level = 3;
573  std::cout << "num entries = " << data.numEntries() << std::endl;
574  bool save = true;
575  (pdf.getVariables())->Print("v"); // print the parameters
576 #else
577  int level = -1;
578  bool save = false;
579 #endif
580 
581 //#ifndef _WIN32 // until a bug 30762 is fixed
582  RooFitResult * result = pdf.fitTo(data, RooFit::Minos(0), RooFit::Hesse(1) , RooFit::PrintLevel(level), RooFit::Save(save) );
583 // #else
584 // RooFitResult * result = pdf.fitTo(data );
585 // #endif
586 
587 #ifdef DEBUG
588  mean.Print();
589  sigma.Print();
590  assert(result != 0);
591  std::cout << " Roofit status " << result->status() << std::endl;
592  result->Print();
593 #endif
594  if (save) iret |= (result == 0);
595 
596  if (iret != 0) {
597  std::cout << "Fit failed " << std::endl;
598  return iret;
599  }
600 
601  }
602 
603  w.Stop();
604  std::cout << "\nTime: \t" << w.RealTime() << " , " << w.CpuTime() << std::endl;
605  std::cout << "\n************************************************************\n";
606  return iret;
607 }
608 
609 //unbinned roo fit (large tree)
611 
612  int iret = 0;
613  std::cout << "\n************************************************************\n";
614  std::cout << "\tFit using RooFit (Likelihood Fit)\n";
615 
616 
617  // start counting the time
618  TStopwatch w; w.Start();
619 
620  for (int i = 0; i < nfit; ++i) {
621 
622  RooArgSet xvars;
623  std::vector<RooRealVar *> x(N);
624  std::vector<RooRealVar *> m(N);
625  std::vector<RooRealVar *> s(N);
626 
627  std::vector<RooGaussian *> g(N);
628  std::vector<RooProdPdf *> pdf(N);
629 
630  for (int j = 0; j < N; ++j) {
631  std::string xname = "x_" + ROOT::Math::Util::ToString(j);
632  x[j] = new RooRealVar(xname.c_str(),xname.c_str(),-10000,10000) ;
633  xvars.add( *x[j] );
634  }
635 
636 
637  RooDataSet data("unbindata","unbin dataset with x",tree,xvars) ;
638 
639  // create the gaussians
640  for (int j = 0; j < N; ++j) {
641  std::string mname = "m_" + ROOT::Math::Util::ToString(j);
642  std::string sname = "s_" + ROOT::Math::Util::ToString(j);
643 
644 
645  m[j] = new RooRealVar(mname.c_str(),mname.c_str(),iniPar[2*j],-100,100) ;
646  s[j] = new RooRealVar(sname.c_str(),sname.c_str(),iniPar[2*j+1],-100,100) ;
647 
648  std::string gname = "g_" + ROOT::Math::Util::ToString(j);
649  g[j] = new RooGaussian(gname.c_str(),"gauss(x,mean,sigma)",*x[j],*m[j],*s[j]);
650 
651  std::string pname = "prod_" + ROOT::Math::Util::ToString(j);
652  if (j == 1) {
653  pdf[1] = new RooProdPdf(pname.c_str(),pname.c_str(),RooArgSet(*g[1],*g[0]) );
654  }
655  else if (j > 1) {
656  pdf[j] = new RooProdPdf(pname.c_str(),pname.c_str(),RooArgSet(*g[j],*pdf[j-1]) );
657  }
658 // else
659 // pdf[0] = g[0];
660 
661  }
662 
663 
664 
665 
666 #ifdef DEBUG
667  int level = 3;
668  std::cout << "num entries = " << data.numEntries() << std::endl;
669  bool save = true;
670  (pdf[N-1]->getVariables())->Print("v"); // print the parameters
671  std::cout << "\n\nDo the fit now \n\n";
672 #else
673  int level = -1;
674  bool save = false;
675 #endif
676 
677 
678 #ifndef _WIN32 // until a bug 30762 is fixed
679  RooFitResult * result = pdf[N-1]->fitTo(data, RooFit::Minos(0), RooFit::Hesse(1) , RooFit::PrintLevel(level), RooFit::Save(save) );
680 #else
681  RooFitResult * result = pdf[N-1]->fitTo(data);
682 #endif
683 
684 #ifdef DEBUG
685  assert(result != 0);
686  std::cout << " Roofit status " << result->status() << std::endl;
687  result->Print();
688 #endif
689 
690 
691  iret |= (result != 0);
692 
693  // free
694  for (int j = 0; j < N; ++j) {
695  delete x[j];
696  delete m[j];
697  delete s[j];
698  delete g[j];
699  if (j> 0) delete pdf[j];
700  }
701 
702  if (iret != 0) return iret;
703  //assert(iret == 0);
704 
705 
706  }
707 
708  w.Stop();
709  std::cout << "\nTime: \t" << w.RealTime() << " , " << w.CpuTime() << std::endl;
710  std::cout << "\n************************************************************\n";
711  return iret;
712 }
713 
714 
715 double poly2(const double *x, const double *p) {
716  return p[0] + (p[1]+p[2]*x[0] ) * x[0];
717 }
718 
719 int testPolyFit() {
720 
721  int iret = 0;
722 
723 
724  std::cout << "\n\n************************************************************\n";
725  std::cout << "\t POLYNOMIAL FIT\n";
726  std::cout << "************************************************************\n";
727 
728  std::string fname("pol2");
729  //TF1 * func = (TF1*)gROOT->GetFunction(fname.c_str());
730  TF1 * f1 = new TF1("pol2",fname.c_str(),-5,5.);
731 
732  f1->SetParameter(0,1);
733  f1->SetParameter(1,0.0);
734  f1->SetParameter(2,1.0);
735 
736 
737  // fill an histogram
738  TH1D * h1 = new TH1D("h1","h1",20,-5.,5.);
739 // h1->FillRandom(fname.c_str(),100);
740  for (int i = 0; i <1000; ++i)
741  h1->Fill( f1->GetRandom() );
742 
743  //h1->Print();
744  //h1->Draw();
745  iniPar[0] = 2.; iniPar[1] = 2.; iniPar[2] = 2.;
746 
747 
748  iret |= FitUsingTFit<TH1,TMINUIT>(h1,f1);
749  iret |= FitUsingTFit<TH1,MINUIT2>(h1,f1);
750  // dummy for testing
751  //iret |= FitUsingNewFitter<DUMMY>(h1,f1);
752 
753  // use simply TF1 wrapper
754  //ROOT::Math::WrappedMultiTF1 f2(*f1);
756 
757 
758  // if Minuit2 is later than TMinuit on Interl is much slower , why ??
759  iret |= FitUsingNewFitter<MINUIT2>(h1,f2);
760  iret |= FitUsingNewFitter<TMINUIT>(h1,f2);
761 
762  // test with linear fitter
763  // for this test need to pass a multi-dim function
764  ROOT::Math::WrappedTF1 wf(*f1);
766  iret |= FitUsingNewFitter<LINEAR>(h1,lfunc,true);
767  iret |= FitUsingTFit<TH1,LINEAR>(h1,f1);
768 
769  // test with a graph
770 
771  gStyle->SetErrorX(0.); // to seto zero error on X
772  TGraphErrors * gr = new TGraphErrors(h1);
773  iret |= FitUsingTFit<TGraph,TMINUIT>(gr,f1);
774 
775  iret |= FitUsingTFit<TGraph,MINUIT2>(gr,f1);
776 
777  iret |= FitUsingNewFitter<MINUIT2>(gr,f2);
778 
779 
780  std::cout << "\n-----> test now TGraphErrors with errors in X coordinates\n\n";
781  // try with error in X
782  gStyle->SetErrorX(0.5); // to set zero error on X
783  TGraphErrors * gr2 = new TGraphErrors(h1);
784  iret |= FitUsingTFit<TGraph,TMINUIT>(gr2,f1);
785  iret |= FitUsingTFit<TGraph,MINUIT2>(gr2,f1);
786 
787  iret |= FitUsingNewFitter<MINUIT2>(gr2,f2);
788 
789  printResult(iret);
790 
791  return iret;
792 }
793 
794 double gaussian(const double *x, const double *p) {
795  //return p[0]*TMath::Gaus(x[0],p[1],p[2]);
796  double tmp = (x[0]-p[1])/p[2];
797  return p[0] * std::exp(-tmp*tmp/2);
798 }
799 
800 double gausnorm(const double *x, const double *p) {
801  //return p[0]*TMath::Gaus(x[0],p[1],p[2]);
802  double invsig = 1./p[1];
803  double tmp = (x[0]-p[0]) * invsig;
804  const double sqrt_2pi = 1./std::sqrt(2.* 3.14159 );
805  return std::exp(-0.5 * tmp*tmp ) * sqrt_2pi * invsig;
806 }
807 double gausnorm2D(const double *x, const double *p) {
808  //return p[0]*TMath::Gaus(x[0],p[1],p[2]);
809  return gausnorm(x,p)*gausnorm(x+1,p+2);
810 }
811 double gausnormN(const double *x, const double *p) {
812  //return p[0]*TMath::Gaus(x[0],p[1],p[2]);
813  double f = 1;
814  for (int i = 0; i < N; ++i)
815  f *= gausnorm(x+i,p+2*i);
816 
817  return f;
818 }
819 
820 int testGausFit() {
821 
822  int iret = 0;
823 
824  std::cout << "\n\n************************************************************\n";
825  std::cout << "\t GAUSSIAN FIT\n";
826  std::cout << "************************************************************\n";
827 
828 
829 
830  //std::string fname = std::string("gaus");
831  //TF1 * func = (TF1*)gROOT->GetFunction(fname.c_str());
832  //TF1 * f1 = new TF1("gaus",fname.c_str(),-5,5.);
833  TF1 * f1 = new TF1("gaussian",gaussian,-5,5.,3);
834  //f2->SetParameters(0,1,1);
835 
836  // fill an histogram
837  int nbin = 10000;
838  TH1D * h2 = new TH1D("h2","h2",nbin,-5.,5.);
839 // h1->FillRandom(fname.c_str(),100);
840  for (int i = 0; i < 10000000; ++i)
841  h2->Fill( gRandom->Gaus(0,10) );
842 
843  iniPar[0] = 100.; iniPar[1] = 2.; iniPar[2] = 2.;
844 
845 
846  // use simply TF1 wrapper
847  //ROOT::Math::WrappedMultiTF1 f2(*f1);
849 
850 
851  iret |= FitUsingNewFitter<MINUIT2>(h2,f2);
852  iret |= FitUsingNewFitter<TMINUIT>(h2,f2);
853 
854 // iret |= FitUsingNewFitter<GSL_PR>(h2,f2);
855 
856 
857  iret |= FitUsingTFit<TH1,TMINUIT>(h2,f1);
858  iret |= FitUsingTFit<TH1,MINUIT2>(h2,f1);
859 
860 
861  iret |= FitUsingNewFitter<GSL_FR>(h2,f2);
862  iret |= FitUsingNewFitter<GSL_PR>(h2,f2);
863  iret |= FitUsingNewFitter<GSL_BFGS>(h2,f2);
864  iret |= FitUsingNewFitter<GSL_BFGS2>(h2,f2);
865 
866 
867  // test also fitting a TGraphErrors with histogram data
868  gStyle->SetErrorX(0.); // to seto zero error on X
869  TGraphErrors * gr = new TGraphErrors(h2);
870 
871  iret |= FitUsingTFit<TGraph,TMINUIT>(gr,f1);
872  iret |= FitUsingTFit<TGraph,MINUIT2>(gr,f1);
873 
874  iret |= FitUsingNewFitter<MINUIT2>(gr,f2);
875 
876  // try with error in X
877  gStyle->SetErrorX(0.5); // to seto zero error on X
878  TGraphErrors * gr2 = new TGraphErrors(h2);
879  iret |= FitUsingTFit<TGraph,TMINUIT>(gr2,f1);
880 
881  iret |= FitUsingNewFitter<MINUIT2>(gr2,f2);
882 
883 
884 
885 //#ifdef LATER
886  // test using grad function
887  std::cout << "\n\nTest Using pre-calculated gradients\n\n";
888  bool useGrad=true;
889  iret |= FitUsingNewFitter<MINUIT2>(h2,f2,useGrad);
890  iret |= FitUsingNewFitter<TMINUIT>(h2,f2,useGrad);
891  iret |= FitUsingNewFitter<GSL_FR>(h2,f2,useGrad);
892  iret |= FitUsingNewFitter<GSL_PR>(h2,f2,useGrad);
893  iret |= FitUsingNewFitter<GSL_BFGS>(h2,f2,useGrad);
894  iret |= FitUsingNewFitter<GSL_BFGS2>(h2,f2,useGrad);
895 
896 
897  // test LS algorithm
898  std::cout << "\n\nTest Least Square algorithms\n\n";
899  iret |= FitUsingNewFitter<GSL_NLS>(h2,f2);
900  iret |= FitUsingNewFitter<FUMILI2>(h2,f2);
901  iret |= FitUsingNewFitter<TFUMILI>(h2,f2);
902 
903 // iret |= FitUsingTFit<TH1,FUMILI2>(h2,f1);
904 // iret |= FitUsingTFit<TH1,TFUMILI>(h2,f1);
905 //#endif
906 
907  //iret |= FitUsingRooFit(h2,f1);
908 
909  printResult(iret);
910 
911  return iret;
912 }
913 
914 int testTreeFit() {
915 
916  std::cout << "\n\n************************************************************\n";
917  std::cout << "\t UNBINNED TREE (GAUSSIAN) FIT\n";
918  std::cout << "************************************************************\n";
919 
920 
921  TTree t1("t1","a simple Tree with simple variables");
922  double x, y;
923  Int_t ev;
924  t1.Branch("x",&x,"x/D");
925  t1.Branch("y",&y,"y/D");
926 // t1.Branch("pz",&pz,"pz/F");
927 // t1.Branch("random",&random,"random/D");
928  t1.Branch("ev",&ev,"ev/I");
929 
930  //fill the tree
931  int nrows = 10000;
932 #ifdef TREE_FIT2D
933  nrows = 10000;
934 #endif
935  for (Int_t i=0;i<nrows;i++) {
936  gRandom->Rannor(x,y);
937  x *= 2; x += 1.;
938  y *= 3; y -= 2;
939 
940  ev = i;
941  t1.Fill();
942 
943  }
944  //t1.Draw("x"); // to select fit variable
945 
946  TF1 * f1 = new TF1("gausnorm", gausnorm, -10,10, 2);
947  TF2 * f2 = new TF2("gausnorm2D", gausnorm2D, -10,10, -10,10, 4);
948 
951 
952 
953  iniPar[0] = 0;
954  iniPar[1] = 1;
955  iniPar[2] = 0;
956  iniPar[3] = 1;
957 
958  // use simply TF1 wrapper
959  //ROOT::Math::WrappedMultiTF1 f2(*f1);
960 
961  int iret = 0;
962 
963  // fit 1D first
964 
965 
966  iret |= FitUsingTTreeFit<MINUIT2>(&t1,f1,"x");
967  iret |= FitUsingTTreeFit<MINUIT2>(&t1,f1,"x");
968 
969  iret |= FitUsingTTreeFit<TMINUIT>(&t1,f1,"x");
970 
971  iret |= FitUsingNewFitter<MINUIT2>(&t1,wf1,false); // not copying the data
972  iret |= FitUsingNewFitter<TMINUIT>(&t1,wf1,false); // not copying the data
973  iret |= FitUsingNewFitter<MINUIT2>(&t1,wf1,true); // copying the data
974  iret |= FitUsingNewFitter<TMINUIT>(&t1,wf1,true); // copying the data
975 
976  // fit 2D
977 
978 
979 
980  iret |= FitUsingTTreeFit<MINUIT2>(&t1,f2,"x:y");
981  iret |= FitUsingTTreeFit<TMINUIT>(&t1,f2,"x:y");
982 
983  iret |= FitUsingNewFitter<MINUIT2>(&t1,wf2, true);
984  iret |= FitUsingNewFitter<MINUIT2>(&t1,wf2, false);
985 
986 
987 
988  iret |= FitUsingRooFit(&t1,f1);
989 
990  printResult(iret);
991  return iret;
992 
993 }
994 
995 int testLargeTreeFit(int nevt = 1000) {
996 
997  std::cout << "\n\n************************************************************\n";
998  std::cout << "\t UNBINNED TREE (GAUSSIAN MULTI-DIM) FIT\n";
999  std::cout << "************************************************************\n";
1000 
1001  TTree t1("t2","a large Tree with simple variables");
1002  double x[N];
1003  Int_t ev;
1004  t1.Branch("x",x,"x[20]/D");
1005  t1.Branch("ev",&ev,"ev/I");
1006 
1007  //fill the tree
1008  TRandom3 r;
1009  for (Int_t i=0;i<nevt;i++) {
1010  for (int j = 0; j < N; ++j) {
1011  double mu = double(j)/10.;
1012  double s = 1.0 + double(j)/10.;
1013  x[j] = r.Gaus(mu,s);
1014  }
1015 
1016  ev = i;
1017  t1.Fill();
1018 
1019  }
1020  //t1.Draw("x"); // to select fit variable
1021 
1022 
1023  for (int i = 0; i <N; ++i) {
1024  iniPar[2*i] = 0;
1025  iniPar[2*i+1] = 1;
1026  }
1027 
1028  // use simply TF1 wrapper
1029  //ROOT::Math::WrappedMultiTF1 f2(*f1);
1031 
1032  int iret = 0;
1033  iret |= FitUsingNewFitter<MINUIT2>(&t1,f2);
1034  iret |= FitUsingNewFitter<TMINUIT>(&t1,f2);
1035  iret |= FitUsingNewFitter<GSL_BFGS2>(&t1,f2);
1036 
1037 
1038 
1039  printResult(iret);
1040  return iret;
1041 
1042 }
1043 int testLargeTreeRooFit(int nevt = 1000) {
1044 
1045  int iret = 0;
1046 
1047  TTree t2("t2b","a large Tree with simple variables");
1048  double x[N];
1049  Int_t ev;
1050  for (int j = 0; j < N; ++j) {
1051  std::string xname = "x_" + ROOT::Math::Util::ToString(j);
1052  std::string xname2 = "x_" + ROOT::Math::Util::ToString(j) + "/D";
1053  t2.Branch(xname.c_str(),&x[j],xname2.c_str());
1054  }
1055  t2.Branch("ev",&ev,"ev/I");
1056  //fill the tree
1057  TRandom3 r;
1058  for (Int_t i=0;i<nevt;i++) {
1059  for (int j = 0; j < N; ++j) {
1060  double mu = double(j)/10.;
1061  double s = 1.0 + double(j)/10.;
1062  x[j] = r.Gaus(mu,s);
1063  }
1064 
1065  ev = i;
1066  t2.Fill();
1067  }
1068 
1069 
1070  for (int i = 0; i <N; ++i) {
1071  iniPar[2*i] = 0;
1072  iniPar[2*i+1] = 1;
1073  }
1074 
1075 
1076  //TF1 * f = new TF1("gausnormN", gausnormN, -100,100, 2*N);
1077 
1078  iret |= FitUsingRooFit2(&t2);
1079 
1080 
1081  printResult(iret);
1082 
1083  return iret;
1084 
1085 }
1086 
1088 
1089  int iret = 0;
1090 
1091 
1092 
1093 #ifndef DEBUG
1094  nfit = 10;
1095 #endif
1096  iret |= testGausFit();
1097 
1098 
1099 #ifdef DEBUG
1100  nfit = 1;
1101 #else
1102  nfit = 1;
1103 #endif
1104  iret |= testTreeFit();
1105 
1106 
1107 
1108  //return iret;
1109 
1110 #ifndef DEBUG
1111  nfit = 1000;
1112 #endif
1113  iret |= testPolyFit();
1114 
1115 
1116 
1117 
1118  nfit = 1;
1119  iret |= testLargeTreeRooFit(500);
1120  iret |= testLargeTreeFit(500);
1121 
1122 
1123  if (iret != 0)
1124  std::cerr << "testFitPerf :\t FAILED " << std::endl;
1125  else
1126  std::cerr << "testFitPerf :\t OK " << std::endl;
1127  return iret;
1128 }
1129 
1130 int main() {
1131  return testFitPerf();
1132 }
1133 
double poly2(const double *x, const double *p)
RooArgSet * getVariables(Bool_t stripDisconnected=kTRUE) const
Return RooArgSet with all variables (tree leaf nodes of expresssion tree)
Definition: RooAbsArg.cxx:2082
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
virtual const double * Parameters() const =0
Access the parameter values.
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
Definition: TH1.cxx:3126
virtual void SetParameters(const Double_t *params)
Definition: TF1.h:588
virtual Double_t GetRandom()
Return a random number following this function shape.
Definition: TF1.cxx:1815
double gaussian(const double *x, const double *p)
virtual void Rannor(Float_t &a, Float_t &b)
Return 2 numbers distributed following a gaussian with mean=0 and sigma=1.
Definition: TRandom.cxx:460
void SetPrintLevel(int level)
set print level
Random number generator class based on M.
Definition: TRandom3.h:27
bool USE_BRANCH
Definition: testFitPerf.cxx:79
Double_t RealTime()
Stop the stopwatch (if it is running) and return the realtime (in seconds) passed between the start a...
Definition: TStopwatch.cxx:110
virtual void SetAddress(void *add)
Set address of this branch.
Definition: TBranch.cxx:2148
virtual Bool_t add(const RooAbsCollection &col, Bool_t silent=kFALSE)
Add a collection of arguments to this collection by calling add() for each element in the source coll...
Definition: RooArgSet.h:86
void Start(Bool_t reset=kTRUE)
Start the stopwatch.
Definition: TStopwatch.cxx:58
void SetTolerance(double tol)
set the tolerance
double gausnorm(const double *x, const double *p)
int DoBinFit(T *hist, Func &func, bool debug=false, bool useGrad=false)
virtual Double_t Gaus(Double_t mean=0, Double_t sigma=1)
Samples a random number from the standard Normal (Gaussian) Distribution with the given mean and sigm...
Definition: TRandom.cxx:235
double T(double x)
Definition: ChebyshevPol.h:34
RooCmdArg PrintLevel(Int_t code)
RooProdPdf is an efficient implementation of a product of PDFs of the form.
Definition: RooProdPdf.h:31
RooCmdArg Minos(Bool_t flag=kTRUE)
R__EXTERN TStyle * gStyle
Definition: TStyle.h:402
virtual Int_t Fill()
Fill all branches.
Definition: TTree.cxx:4383
static void SetDefaultFitter(const char *name="")
static: set name of default fitter
Class to Wrap a ROOT Function class (like TF1) in a IParamFunction interface of one dimensions to be ...
Definition: WrappedTF1.h:37
Double_t CpuTime()
Stop the stopwatch (if it is running) and return the cputime (in seconds) passed between the start an...
Definition: TStopwatch.cxx:125
Class describing the unbinned data sets (just x coordinates values) of any dimensions.
Definition: UnBinData.h:42
int Int_t
Definition: RtypesCore.h:41
ROOT::Math::MinimizerOptions & MinimizerOptions()
access to the minimizer control parameter (non const method)
Definition: FitConfig.h:138
TLatex * t1
Definition: textangle.C:20
virtual void Print(Option_t *options=0) const
Print TNamed name and title.
Definition: RooFitResult.h:64
int testLargeTreeRooFit(int nevt=1000)
const double * Coords(unsigned int ipoint) const
return a pointer to the coordinates data for the given fit point
Definition: FitData.h:245
ROOT::Math::IParamMultiFunction Func
const FitResult & Result() const
get fit result
Definition: Fitter.h:337
int testFitPerf()
double sqrt(double)
double iniPar[2 *N]
Definition: testFitPerf.cxx:61
RooDataSet is a container class to hold N-dimensional binned data.
Definition: RooDataHist.h:40
void Stop()
Stop the stopwatch.
Definition: TStopwatch.cxx:77
Double_t x[n]
Definition: legend1.C:17
virtual Long64_t GetSelectedRows()
Definition: TTree.h:428
virtual Int_t UnbinnedFit(const char *funcname, const char *varexp, const char *selection="", Option_t *option="", Long64_t nentries=kMaxEntries, Long64_t firstentry=0)
Unbinned fit of one or more variable(s) from a tree.
Definition: TTree.cxx:9114
Int_t setPrintLevel(Int_t newLevel)
Change the MINUIT internal printing level.
Definition: RooMinuit.cxx:569
Plain Gaussian p.d.f.
Definition: RooGaussian.h:25
virtual void Print(Option_t *options=0) const
Print TNamed name and title.
Definition: RooAbsArg.h:227
void printResult(int iret)
Definition: testFitPerf.cxx:71
const Double_t sigma
const int N
Definition: testFitPerf.cxx:60
TH1F * h1
Definition: legend1.C:5
double gausnormN(const double *x, const double *p)
const FitConfig & Config() const
access to the fit configuration (const method)
Definition: Fitter.h:365
unsigned int Size() const
return number of fit points
Definition: FitData.h:302
RooRealVar represents a fundamental (non-derived) real valued object.
Definition: RooRealVar.h:36
virtual TBranch * GetBranch(const char *name)
Return pointer to the branch with the given name in this tree or its friends.
Definition: TTree.cxx:4979
ROOT::Fit::UnBinData * FillUnBinData(TTree *tree, bool copyData=true, unsigned int dim=1)
Definition: testFitPerf.cxx:80
bool Fit(const Data &data, const Function &func)
fit a data set using any generic model function If data set is binned a least square fit is performed...
Definition: Fitter.h:129
virtual void SetParameters(const double *p)
Set the parameter values.
void SetMinimizer(const char *type, const char *algo=0)
set minimizer type
Definition: FitConfig.h:152
void FillData(BinData &dv, const TH1 *hist, TF1 *func=0)
fill the data vector from a TH1.
virtual void SetEstimate(Long64_t nentries=1000000)
Set number of entries to estimate variable limits.
Definition: TTree.cxx:8580
int main()
TRandom2 r(17)
SVector< double, 2 > v
Definition: Dict.h:5
virtual Double_t * GetV2()
Definition: TTree.h:452
Fitter class, entry point for performing all type of fits.
Definition: Fitter.h:77
TMarker * m
Definition: textangle.C:8
int nfit
Definition: testFitPerf.cxx:59
int FitUsingNewFitter(FitObj *fitobj, Func &func, bool useGrad=false)
virtual Int_t GetEntry(Long64_t entry=0, Int_t getall=0)
Read all leaves of entry and return total number of bytes read.
Definition: TBranch.cxx:1291
Class describing the binned data sets : vectors of x coordinates, y values and optionally error on y ...
Definition: BinData.h:53
const std::string sname
Definition: testIO.cxx:45
A 2-Dim function with parameters.
Definition: TF2.h:29
R__EXTERN TRandom * gRandom
Definition: TRandom.h:62
tomato 1-D histogram with a double per channel (see TH1 documentation)}
Definition: TH1.h:594
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:29
void SetFunction(const IModelFunction &func, bool useGradient=false)
Set the fitted function (model function) from a parametric function interface.
Definition: Fitter.cxx:104
Parametric Function class describing polynomials of order n.
Definition: Polynomial.h:63
TGraphErrors * gr
Definition: legend1.C:25
IParamFunction interface (abstract class) describing multi-dimensional parameteric functions It is a ...
int FitUsingTTreeFit(TTree *tree, TF1 *func, const std::string &vars="x")
Double_t GetChisquare() const
Definition: TF1.h:398
virtual void SetParameters(const double *p)=0
Set the parameter values.
RooCmdArg Hesse(Bool_t flag=kTRUE)
double f(double x)
void Print(std::ostream &os, const OptionType &opt)
int DoUnBinFit(T *tree, Func &func, bool debug=false, bool copyData=false)
void Add(double x)
preallocate a data set given size and dimension of the coordinates if a vector already exists with co...
Definition: UnBinData.h:200
virtual void Draw(Option_t *opt)
Default Draw method for all objects.
Definition: TTree.h:355
Double_t y[n]
Definition: legend1.C:17
double func(double *x, double *p)
Definition: stressTF1.cxx:213
The TH1 histogram class.
Definition: TH1.h:56
void printData(const ROOT::Fit::UnBinData &data)
Definition: testFitPerf.cxx:64
double gausnorm2D(const double *x, const double *p)
void Append(unsigned int newPoints, unsigned int dim=1, bool isWeighted=false)
Definition: UnBinData.h:288
virtual Long64_t GetEntries() const
Definition: TTree.h:381
virtual Int_t Branch(TCollection *list, Int_t bufsize=32000, Int_t splitlevel=99, const char *name="")
Create one branch for each element in the collection.
Definition: TTree.cxx:1660
RooCmdArg Save(Bool_t flag=kTRUE)
virtual void GetParLimits(Int_t ipar, Double_t &parmin, Double_t &parmax) const
Return limits for parameter ipar.
Definition: TF1.cxx:1666
MultiDimParamGradFunctionAdapter class to wrap a one-dimensional parametric gradient function in a mu...
RooFitResult * fit(const char *options)
Parse traditional RooAbsPdf::fitTo driver options.
Definition: RooMinuit.cxx:273
int testGausFit()
RooAbsPdf is the abstract interface for all probability density functions The class provides hybrid a...
Definition: RooAbsPdf.h:41
void Print(std::ostream &os, bool covmat=false) const
print the result and optionaly covariance matrix and correlations
Definition: FitResult.cxx:446
std::string ToString(const T &val)
Utility function for conversion to strings.
Definition: Util.h:40
double f2(const double *x)
virtual Double_t GetParameter(Int_t ipar) const
Definition: TF1.h:466
1-Dim function class
Definition: TF1.h:150
A Graph is a graphics object made of two arrays X and Y with npoints each.
Definition: TGraph.h:41
TF1 * f1
Definition: legend1.C:11
A TGraphErrors is a TGraph with error bars.
Definition: TGraphErrors.h:26
int testPolyFit()
THist< 1, double, THistStatContent, THistStatUncertainty > TH1D
Definition: THist.hxx:310
bool debug
virtual RooFitResult * fitTo(RooAbsData &data, const RooCmdArg &arg1=RooCmdArg::none(), const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none())
Fit PDF to given dataset.
Definition: RooAbsPdf.cxx:1056
Definition: tree.py:1
int testTreeFit()
A TTree object has a header with a name and a title.
Definition: TTree.h:78
double result[121]
void SetParameters(const double *p)
Set the parameter values.
virtual void SetParameter(Int_t param, Double_t value)
Definition: TF1.h:578
int FitUsingTFit(T *hist, TF1 *func)
A TTree is a list of TBranches.
Definition: TBranch.h:57
double exp(double)
void SetErrorX(Float_t errorx=0.5)
Definition: TStyle.h:315
unsigned int NPar() const
Return the number of parameters.
int FitUsingRooFit2(TTree *tree)
const Int_t n
Definition: legend1.C:16
int FitUsingRooFit(TH1 *hist, TF1 *func)
virtual Double_t * GetV1()
Definition: TTree.h:450
virtual Int_t numEntries() const
Definition: RooAbsData.cxx:269
WrappedParamFunction class to wrap any multi-dimensional function pbject implementing the operator()(...
virtual unsigned int NDim() const =0
Retrieve the dimension of the function.
int DoFit(TTree *tree, Func &func, bool debug=false, bool copyData=false)
RooMinuit is a wrapper class around TFitter/TMinuit that provides a seamless interface between the MI...
Definition: RooMinuit.h:39
Stopwatch class.
Definition: TStopwatch.h:28
int testLargeTreeFit(int nevt=1000)
Int_t status() const
Definition: RooFitResult.h:75