Logo ROOT  
Reference Guide
FittingDemo.C
Go to the documentation of this file.
1/// \file
2/// \ingroup tutorial_fit
3/// \notebook -js
4/// Example for fitting signal/background.
5/// This example can be executed with:
6///
7/// ~~~{.cpp}
8/// root > .x FittingDemo.C (using the CINT interpreter)
9/// root > .x FittingDemo.C+ (using the native complier via ACLIC)
10/// ~~~
11///
12/// \macro_image
13/// \macro_output
14/// \macro_code
15///
16/// \author Rene Brun
17
18#include "TH1.h"
19#include "TMath.h"
20#include "TF1.h"
21#include "TLegend.h"
22#include "TCanvas.h"
23
24// Quadratic background function
25Double_t background(Double_t *x, Double_t *par) {
26 return par[0] + par[1]*x[0] + par[2]*x[0]*x[0];
27}
28
29
30// Lorenzian Peak function
31Double_t lorentzianPeak(Double_t *x, Double_t *par) {
32 return (0.5*par[0]*par[1]/TMath::Pi()) /
33 TMath::Max( 1.e-10,(x[0]-par[2])*(x[0]-par[2])
34 + .25*par[1]*par[1]);
35}
36
37// Sum of background and peak function
38Double_t fitFunction(Double_t *x, Double_t *par) {
39 return background(x,par) + lorentzianPeak(x,&par[3]);
40}
41
42void FittingDemo() {
43 //Bevington Exercise by Peter Malzacher, modified by Rene Brun
44
45 const int nBins = 60;
46
47 Double_t data[nBins] = { 6, 1,10,12, 6,13,23,22,15,21,
48 23,26,36,25,27,35,40,44,66,81,
49 75,57,48,45,46,41,35,36,53,32,
50 40,37,38,31,36,44,42,37,32,32,
51 43,44,35,33,33,39,29,41,32,44,
52 26,39,29,35,32,21,21,15,25,15};
53 TCanvas *c1 = new TCanvas("c1","Fitting Demo",10,10,700,500);
54 c1->SetFillColor(33);
55 c1->SetFrameFillColor(41);
56 c1->SetGrid();
57
58 TH1F *histo = new TH1F("histo",
59 "Lorentzian Peak on Quadratic Background",60,0,3);
60 histo->SetMarkerStyle(21);
61 histo->SetMarkerSize(0.8);
62 histo->SetStats(0);
63
64 for(int i=0; i < nBins; i++) histo->SetBinContent(i+1,data[i]);
65
66 // create a TF1 with the range from 0 to 3 and 6 parameters
67 TF1 *fitFcn = new TF1("fitFcn",fitFunction,0,3,6);
68 fitFcn->SetNpx(500);
69 fitFcn->SetLineWidth(4);
70 fitFcn->SetLineColor(kMagenta);
71
72 // first try without starting values for the parameters
73 // This defaults to 1 for each param.
74 // this results in an ok fit for the polynomial function
75 // however the non-linear part (lorenzian) does not
76 // respond well.
77 fitFcn->SetParameters(1,1,1,1,1,1);
78 histo->Fit("fitFcn","0");
79
80 // second try: set start values for some parameters
81 fitFcn->SetParameter(4,0.2); // width
82 fitFcn->SetParameter(5,1); // peak
83
84 histo->Fit("fitFcn","V+","ep");
85
86 // improve the picture:
87 TF1 *backFcn = new TF1("backFcn",background,0,3,3);
88 backFcn->SetLineColor(kRed);
89 TF1 *signalFcn = new TF1("signalFcn",lorentzianPeak,0,3,3);
90 signalFcn->SetLineColor(kBlue);
91 signalFcn->SetNpx(500);
92 Double_t par[6];
93
94 // writes the fit results into the par array
95 fitFcn->GetParameters(par);
96
97 backFcn->SetParameters(par);
98 backFcn->Draw("same");
99
100 signalFcn->SetParameters(&par[3]);
101 signalFcn->Draw("same");
102
103 // draw the legend
104 TLegend *legend=new TLegend(0.6,0.65,0.88,0.85);
105 legend->SetTextFont(72);
106 legend->SetTextSize(0.04);
107 legend->AddEntry(histo,"Data","lpe");
108 legend->AddEntry(backFcn,"Background fit","l");
109 legend->AddEntry(signalFcn,"Signal fit","l");
110 legend->AddEntry(fitFcn,"Global Fit","l");
111 legend->Draw();
112
113}
double Double_t
Definition: RtypesCore.h:57
@ kRed
Definition: Rtypes.h:64
@ kMagenta
Definition: Rtypes.h:64
@ kBlue
Definition: Rtypes.h:64
virtual void SetLineWidth(Width_t lwidth)
Set the line width.
Definition: TAttLine.h:43
virtual void SetLineColor(Color_t lcolor)
Set the line color.
Definition: TAttLine.h:40
virtual void SetMarkerStyle(Style_t mstyle=1)
Set the marker style.
Definition: TAttMarker.h:40
virtual void SetMarkerSize(Size_t msize=1)
Set the marker size.
Definition: TAttMarker.h:41
virtual void SetTextFont(Font_t tfont=62)
Set the text font.
Definition: TAttText.h:45
virtual void SetTextSize(Float_t tsize=1)
Set the text size.
Definition: TAttText.h:46
The Canvas class.
Definition: TCanvas.h:27
1-Dim function class
Definition: TF1.h:210
virtual void SetNpx(Int_t npx=100)
Set the number of points used to draw the function.
Definition: TF1.cxx:3435
virtual Double_t * GetParameters() const
Definition: TF1.h:514
virtual void Draw(Option_t *option="")
Draw this function with its current attributes.
Definition: TF1.cxx:1320
virtual void SetParameters(const Double_t *params)
Definition: TF1.h:638
virtual void SetParameter(Int_t param, Double_t value)
Definition: TF1.h:628
1-D histogram with a float per channel (see TH1 documentation)}
Definition: TH1.h:571
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:3808
virtual void SetBinContent(Int_t bin, Double_t content)
Set bin content see convention for numbering bins in TH1::GetBin In case the bin number is greater th...
Definition: TH1.cxx:8678
virtual void SetStats(Bool_t stats=kTRUE)
Set statistics option on/off.
Definition: TH1.cxx:8446
This class displays a legend box (TPaveText) containing several legend entries.
Definition: TLegend.h:23
TLegendEntry * AddEntry(const TObject *obj, const char *label="", Option_t *option="lpf")
Add a new entry to this legend.
Definition: TLegend.cxx:330
virtual void Draw(Option_t *option="")
Draw this legend with its current attributes.
Definition: TLegend.cxx:423
return c1
Definition: legend1.C:41
Double_t x[n]
Definition: legend1.C:17
Short_t Max(Short_t a, Short_t b)
Definition: TMathBase.h:212
constexpr Double_t Pi()
Definition: TMath.h:38