Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
rf212_plottingInRanges_blinding.C
Go to the documentation of this file.
1/// \file
2/// \ingroup tutorial_roofit
3/// \notebook -js
4/// Plot a PDF in disjunct ranges, and get normalisation right.
5///
6/// Usually, when comparing a fit to data, one should first plot the data, and then the PDF.
7/// In this case, the PDF is automatically normalised to match the number of data events in the plot.
8/// However, when plotting only a sub-range, when e.g. a signal region has to be blinded,
9/// one has to exclude the blinded region from the computation of the normalisation.
10///
11/// In this tutorial, we show how to explicitly choose the normalisation when plotting using `NormRange()`.
12///
13/// Thanks to Marc Escalier for asking how to do this correctly.
14///
15/// \macro_image
16/// \macro_code
17/// \macro_output
18///
19/// \date March 2020
20/// \author Stephan Hageboeck
21
22#include <RooDataSet.h>
23#include <RooExponential.h>
24#include <RooPlot.h>
25#include <RooRealVar.h>
26#include <TCanvas.h>
27
28using namespace RooFit;
29
30void rf212_plottingInRanges_blinding()
31{
32 // Make a fit model
33 RooRealVar x("x", "The observable", 1, 30);
34 RooRealVar tau("tau", "The exponent", -0.1337, -10., -0.1);
35 RooExponential exp("exp", "A falling exponential function", x, tau);
36
37 // Define the sidebands (e.g. background regions)
38 x.setRange("full", 1, 30);
39 x.setRange("left", 1, 10);
40 x.setRange("right", 20, 30);
41
42 // Generate toy data, and cut out the blinded region.
43 RooDataSet* data = exp.generate(x, 1000);
44 auto blindedData = data->reduce(CutRange("left,right"));
45
46 // Kick tau a bit, and run an unbinned fit where the blinded data are missing.
47 // ----------------------------------------------------------------------------------------------------------
48 tau.setVal(-2.);
49 exp.fitTo(*blindedData);
50
51
52 // Here we will plot the results
53 TCanvas *canvas=new TCanvas("canvas","canvas",800,600);
54 canvas->Divide(2,1);
55
56
57 // Wrong:
58 // Plotting each slice on its own normalises the PDF over its plotting range. For the full curve, that means
59 // that the blinded region where data is missing is included in the normalisation calculation. The PDF therefore
60 // comes out too low, and doesn't match up with the slices in the side bands, which are normalised to "their" data.
61 // ----------------------------------------------------------------------------------------------------------
62
63 std::cout << "Now plotting with unique normalisation for each slice." << std::endl;
64 canvas->cd(1);
65 RooPlot* plotFrame = x.frame(RooFit::Title("Wrong: Each slice normalised over its plotting range"));
66
67 // Plot only the blinded data, and then plot the PDF over the full range as well as both sidebands
68 blindedData->plotOn(plotFrame);
69 exp.plotOn(plotFrame, LineColor(kRed), Range("full"));
70 exp.plotOn(plotFrame, LineColor(kBlue), Range("left"));
71 exp.plotOn(plotFrame, LineColor(kGreen), Range("right"));
72
73 plotFrame->Draw();
74
75 // Right:
76 // Make the same plot, but normalise each piece with respect to the regions "left" AND "right". This requires setting
77 // a "NormRange", which tells RooFit over which range the PDF has to be integrated to normalise.
78 // This is means that the normalisation of the blue and green curves is slightly different from the left plot,
79 // because they get a common scale factor.
80 // ----------------------------------------------------------------------------------------------------------
81
82 std::cout << "\n\nNow plotting with correct norm ranges:" << std::endl;
83 canvas->cd(2);
84 RooPlot* plotFrameWithNormRange = x.frame(RooFit::Title("Right: All slices have common normalisation"));
85
86 // Plot only the blinded data, and then plot the PDF over the full range as well as both sidebands
87 blindedData->plotOn(plotFrameWithNormRange);
88 exp.plotOn(plotFrameWithNormRange, LineColor(kBlue), Range("left"), RooFit::NormRange("left,right"));
89 exp.plotOn(plotFrameWithNormRange, LineColor(kGreen), Range("right"), RooFit::NormRange("left,right"));
90 exp.plotOn(plotFrameWithNormRange, LineColor(kRed), Range("full"), RooFit::NormRange("left,right"), LineStyle(10));
91
92 plotFrameWithNormRange->Draw();
93
94 canvas->Draw();
95
96}
@ kRed
Definition Rtypes.h:66
@ kGreen
Definition Rtypes.h:66
@ kBlue
Definition Rtypes.h:66
double exp(double)
RooAbsData * reduce(const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg(), const RooCmdArg &arg3=RooCmdArg(), const RooCmdArg &arg4=RooCmdArg(), const RooCmdArg &arg5=RooCmdArg(), const RooCmdArg &arg6=RooCmdArg(), const RooCmdArg &arg7=RooCmdArg(), const RooCmdArg &arg8=RooCmdArg())
Create a reduced copy of this dataset.
RooDataSet is a container class to hold unbinned data.
Definition RooDataSet.h:33
Exponential PDF.
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
void Draw(Option_t *option="") override
Draw a canvas.
Definition TCanvas.cxx:845
TVirtualPad * cd(Int_t subpadnumber=0) override
Set current canvas & pad.
Definition TCanvas.cxx:708
void Divide(Int_t nx=1, Int_t ny=1, Float_t xmargin=0.01, Float_t ymargin=0.01, Int_t color=0) override
Automatic pad generation by division.
Definition TPad.cxx:1177
RooCmdArg Title(const char *name)
RooCmdArg NormRange(const char *rangeNameList)
RooCmdArg CutRange(const char *rangeName)
RooCmdArg LineColor(Color_t color)
RooCmdArg LineStyle(Style_t style)
Double_t x[n]
Definition legend1.C:17
The namespace RooFit contains mostly switches that change the behaviour of functions of PDFs (or othe...
Ta Range(0, 0, 1, 1)