Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
vectorizedFit.C
Go to the documentation of this file.
1/// \file
2/// \ingroup tutorial_fit
3/// \notebook
4/// Tutorial for creating a Vectorized TF1 function using a formula expression and
5/// use it for fitting an histogram
6///
7/// To create a vectorized function (if ROOT has been compiled with support for vectorization)
8/// is very easy. One needs to create the TF1 object with the option "VEC" or call the method
9/// TF1::SetVectorized
10///
11/// \macro_image
12/// \macro_output
13/// \macro_code
14///
15/// \author Lorenzo Moneta
16
18#include <TCanvas.h>
19#include <TF1.h>
20#include <TH1D.h>
21#include <TStopwatch.h>
22#include <TStyle.h>
23
24#include <iostream>
25
26void vectorizedFit() {
27
28 gStyle->SetOptFit(111111);
29
30
32
33 int nbins = 40000;
34 auto h1 = new TH1D("h1","h1",nbins,-3,3);
35 h1->FillRandom("gaus",nbins*50);
36 auto c1 = new TCanvas("Fit","Fit",800,1000);
37 c1->Divide(1,2);
38 c1->cd(1);
40
41 std::cout << "Doing Serial Gaussian Fit " << std::endl;
42 auto f1 = new TF1("f1","gaus");
43 f1->SetNpx(nbins*10);
44 w.Start();
45 h1->Fit(f1);
46 h1->Fit(f1,"L+");
47 w.Print();
48
49 std::cout << "Doing Vectorized Gaussian Fit " << std::endl;
50 auto f2 = new TF1("f2","gaus",-3,3,"VEC");
51 // alternatively you can also use the TF1::SetVectorized function
52 //f2->SetVectorized(true);
53 w.Start();
54 h1->Fit(f2);
55 h1->Fit(f2,"L+");
56 w.Print();
57 // rebin histograms and scale it back to the function
58 h1->Rebin(nbins/100);
59 h1->Scale(100./nbins);
60 ((TF1 *)h1->GetListOfFunctions()->At(0))->SetTitle("Chi2 Fit");
61 ((TF1 *)h1->GetListOfFunctions()->At(1))->SetTitle("Likelihood Fit");
63 //c1->cd(1)->BuildLegend();
64
65 /// Do a polynomial fit now
66 c1->cd(2);
67 auto f3 = new TF1("f3","[A]*x^2+[B]*x+[C]",0,10);
68 f3->SetParameters(0.5,3,2);
69 f3->SetNpx(nbins*10);
70 // generate the events
71 auto h2 = new TH1D("h2","h2",nbins,0,10);
72 h2->FillRandom("f3",10*nbins);
73 std::cout << "Doing Serial Polynomial Fit " << std::endl;
74 f3->SetParameters(2,2,2);
75 w.Start();
76 h2->Fit(f3);
77 h2->Fit(f3,"L+");
78 w.Print();
79
80 std::cout << "Doing Vectorized Polynomial Fit " << std::endl;
81 auto f4 = new TF1("f4","[A]*x*x+[B]*x+[C]",0,10);
82 f4->SetVectorized(true);
83 f4->SetParameters(2,2,2);
84 w.Start();
85 h2->Fit(f4);
86 h2->Fit(f4,"L+");
87 w.Print();
88
89 // rebin histograms and scale it back to the function
90 h2->Rebin(nbins/100);
91 h2->Scale(100./nbins);
92 ((TF1 *)h2->GetListOfFunctions()->At(0))->SetTitle("Chi2 Fit");
93 ((TF1 *)h2->GetListOfFunctions()->At(1))->SetTitle("Likelihood Fit");
94 ((TF1 *)h2->GetListOfFunctions()->At(1))->SetLineColor(kBlue);
95 //c1->cd(2)->BuildLegend();
96
97}
@ kBlue
Definition Rtypes.h:66
Option_t Option_t SetLineColor
R__EXTERN TStyle * gStyle
Definition TStyle.h:433
static void SetDefaultMinimizer(const char *type, const char *algo=nullptr)
Set the default Minimizer type and corresponding algorithms.
The Canvas class.
Definition TCanvas.h:23
1-Dim function class
Definition TF1.h:233
virtual void SetNpx(Int_t npx=100)
Set the number of points used to draw the function.
Definition TF1.cxx:3433
1-D histogram with a double per channel (see TH1 documentation)
Definition TH1.h:664
virtual void FillRandom(const char *fname, Int_t ntimes=5000, TRandom *rng=nullptr)
Fill histogram following distribution in function fname.
Definition TH1.cxx:3515
virtual TFitResultPtr Fit(const char *formula, Option_t *option="", Option_t *goption="", Double_t xmin=0, Double_t xmax=0)
Fit histogram with function fname.
Definition TH1.cxx:3894
virtual TH1 * Rebin(Int_t ngroup=2, const char *newname="", const Double_t *xbins=nullptr)
Rebin this histogram.
Definition TH1.cxx:6239
TList * GetListOfFunctions() const
Definition TH1.h:244
virtual void Scale(Double_t c1=1, Option_t *option="")
Multiply this histogram by a constant c1.
Definition TH1.cxx:6568
TObject * At(Int_t idx) const override
Returns the object at position idx. Returns 0 if idx is out of range.
Definition TList.cxx:355
Stopwatch class.
Definition TStopwatch.h:28
void SetOptFit(Int_t fit=1)
The type of information about fit parameters printed in the histogram statistics box can be selected ...
Definition TStyle.cxx:1589
return c1
Definition legend1.C:41
TH1F * h1
Definition legend1.C:5
TF1 * f1
Definition legend1.C:11