Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
rf204a_extendedLikelihood.C
Go to the documentation of this file.
1/// \file
2/// \ingroup tutorial_roofit
3/// \notebook -js
4/// Extended maximum likelihood fit in multiple ranges.
5///
6/// When an extended pdf and multiple ranges are used, the
7/// RooExtendPdf cannot correctly interpret the coefficients
8/// used for extension.
9/// This can be solved by using a RooAddPdf for extending the model.
10///
11/// \macro_output
12/// \macro_code
13///
14/// \date December 2018
15/// \author Stephan Hageboeck
16
17
18#include "RooRealVar.h"
19#include "RooDataSet.h"
20#include "RooGaussian.h"
21#include "RooChebychev.h"
22#include "RooAddPdf.h"
23#include "RooExtendPdf.h"
24#include "RooFitResult.h"
25#include "TCanvas.h"
26#include "TAxis.h"
27#include "RooPlot.h"
28using namespace RooFit ;
29
30
31void rf204a_extendedLikelihood()
32{
33
34
35 // S e t u p c o m p o n e n t p d f s
36 // ---------------------------------------
37
38 // Declare observable x
39 RooRealVar x("x","x",0,11) ;
40
41 // Create two Gaussian PDFs g1(x,mean1,sigma) anf g2(x,mean2,sigma) and their parameters
42 RooRealVar mean("mean","mean of gaussians",5) ;
43 RooRealVar sigma1("sigma1","width of gaussians",0.5) ;
44 RooRealVar sigma2("sigma2","width of gaussians",1) ;
45
46 RooGaussian sig1("sig1","Signal component 1",x,mean,sigma1) ;
47 RooGaussian sig2("sig2","Signal component 2",x,mean,sigma2) ;
48
49 // Build Chebychev polynomial pdf
50 RooRealVar a0("a0","a0",0.5,0.,1.) ;
51 RooRealVar a1("a1","a1",0.2,0.,1.) ;
52 RooChebychev bkg("bkg","Background",x,RooArgSet(a0,a1)) ;
53
54 // Sum the signal components into a composite signal pdf
55 RooRealVar sig1frac("sig1frac","fraction of component 1 in signal",0.8,0.,1.) ;
56 RooAddPdf sig("sig","Signal",RooArgList(sig1,sig2),sig1frac) ;
57
58
59 // E x t e n d t h e p d f s
60 // -----------------------------
61
62
63 // Define signal range in which events counts are to be defined
64 x.setRange("signalRange",4,6) ;
65
66 // Associated nsig/nbkg as expected number of events with sig/bkg _in_the_range_ "signalRange"
67 RooRealVar nsig("nsig","number of signal events in signalRange",500,0.,10000) ;
68 RooRealVar nbkg("nbkg","number of background events in signalRange",500,0,10000) ;
69
70 // Use AddPdf to extend the model. Giving as many coefficients as pdfs switches
71 // on extension.
72 RooAddPdf model("model","(g1+g2)+a", RooArgList(bkg,sig), RooArgList(nbkg,nsig)) ;
73
74
75 // S a m p l e d a t a , f i t m o d e l
76 // -------------------------------------------
77
78 // Generate 1000 events from model so that nsig,nbkg come out to numbers <<500 in fit
79 RooDataSet *data = model.generate(x,1000) ;
80
81
82
83 auto canv = new TCanvas("Canvas", "Canvas", 1500, 600);
84 canv->Divide(3,1);
85
86 // Fit full range
87 // -------------------------------------------
88
89 canv->cd(1);
90
91 // Perform unbinned ML fit to data, full range
92
93 // IMPORTANT:
94 // The model needs to be copied when fitting with different ranges because
95 // the interpretation of the coefficients is tied to the fit range
96 // that's used in the first fit
97 RooAddPdf model1(model);
98 RooFitResult* r = model1.fitTo(*data,Save()) ;
99 r->Print() ;
100
101 RooPlot * frame = x.frame(Title("Full range fitted"));
102 data->plotOn(frame);
103 model1.plotOn(frame, VisualizeError(*r));
104 model1.plotOn(frame);
105 model1.paramOn(frame);
106 frame->Draw();
107
108
109 // Fit in two regions
110 // -------------------------------------------
111
112 canv->cd(2);
113 x.setRange("left", 0., 4.);
114 x.setRange("right", 6., 10.);
115
116 RooAddPdf model2(model);
117 RooFitResult* r2 = model2.fitTo(*data,
118 Range("left,right"),
119 Save()) ;
120 r2->Print();
121
122
123 RooPlot * frame2 = x.frame(Title("Fit in left/right sideband"));
124 data->plotOn(frame2);
125 model2.plotOn(frame2, VisualizeError(*r2));
126 model2.plotOn(frame2);
127 model2.paramOn(frame2);
128 frame2->Draw();
129
130
131 // Fit in one region
132 // -------------------------------------------
133 // Note how restricting the region to only the left tail increases
134 // the fit uncertainty
135
136 canv->cd(3);
137 x.setRange("leftToMiddle", 0., 5.);
138
139 RooAddPdf model3(model);
140 RooFitResult* r3 = model3.fitTo(*data,
141 Range("leftToMiddle"),
142 Save()) ;
143 r3->Print();
144
145
146 RooPlot * frame3 = x.frame(Title("Fit from left to middle"));
147 data->plotOn(frame3);
148 model3.plotOn(frame3, VisualizeError(*r3));
149 model3.plotOn(frame3);
150 model3.paramOn(frame3);
151 frame3->Draw();
152
153 canv->Draw();
154}
ROOT::R::TRInterface & r
Definition Object.C:4
virtual RooPlot * plotOn(RooPlot *frame, 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()) const
RooAddPdf is an efficient implementation of a sum of PDFs of the form.
Definition RooAddPdf.h:32
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition RooArgList.h:21
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition RooArgSet.h:29
Chebychev polynomial p.d.f.
RooDataSet is a container class to hold unbinned data.
Definition RooDataSet.h:33
RooFitResult is a container class to hold the input and output of a PDF fit to a dataset.
virtual void Print(Option_t *options=0) const
Print TNamed name and title.
Plain Gaussian p.d.f.
Definition RooGaussian.h:24
A RooPlot is a plot frame and a container for graphics objects within that frame.
Definition RooPlot.h:44
static RooPlot * frame(const RooAbsRealLValue &var, Double_t xmin, Double_t xmax, Int_t nBins)
Create a new frame for a given variable in x.
Definition RooPlot.cxx:249
virtual void Draw(Option_t *options=0)
Draw this plot and all of the elements it contains.
Definition RooPlot.cxx:691
RooRealVar represents a variable that can be changed from the outside.
Definition RooRealVar.h:39
The Canvas class.
Definition TCanvas.h:23
virtual void Print(Option_t *option="") const
This method must be overridden when a class wants to print itself.
Definition TObject.cxx:552
RooCmdArg Save(Bool_t flag=kTRUE)
RooCmdArg VisualizeError(const RooDataSet &paramData, Double_t Z=1)
Double_t x[n]
Definition legend1.C:17
The namespace RooFit contains mostly switches that change the behaviour of functions of PDFs (or othe...
const char * Title
Definition TXMLSetup.cxx:68
Ta Range(0, 0, 1, 1)