Usually, when comparing a fit to data, one should first plot the data, and then the PDF. In this case, the PDF is automatically normalised to match the number of data events in the plot. However, when plotting only a sub-range, when e.g. a signal region has to be blinded, one has to exclude the blinded region from the computation of the normalisation.
In this tutorial, we show how to explicitly choose the normalisation when plotting using NormRange()
.
import ROOT
x = ROOT.RooRealVar("x", "The observable", 1, 30)
tau = ROOT.RooRealVar("tau", "The exponent", -0.1337, -10.0, -0.1)
exp = ROOT.RooExponential("exp", "A falling exponential function", x, tau)
x.setRange("full", 1, 30)
x.setRange("left", 1, 10)
x.setRange("right", 20, 30)
data = exp.generate(x, 1000)
blindedData = data.reduce(CutRange="left,right")
tau.setVal(-2.0)
exp.fitTo(blindedData)
canvas = ROOT.TCanvas("canvas", "canvas", 800, 600)
canvas.Divide(2, 1)
print("Now plotting with unique normalisation for each slice.\n")
canvas.cd(1)
plotFrame = x.frame(Title="Wrong: Each slice normalised over its plotting range")
blindedData.plotOn(plotFrame)
exp.plotOn(plotFrame, LineColor="r", Range="full")
exp.plotOn(plotFrame, LineColor="b", Range="left")
exp.plotOn(plotFrame, LineColor="g", Range="right")
plotFrame.Draw()
print("\n\nNow plotting with correct norm ranges:\n")
canvas.cd(2)
plotFrameWithNormRange = x.frame(Title="Right: All slices have common normalisation")
blindedData.plotOn(plotFrameWithNormRange)
exp.plotOn(plotFrameWithNormRange, LineColor="b", Range="left", NormRange="left,right")
exp.plotOn(plotFrameWithNormRange, LineColor="g", Range="right", NormRange="left,right")
exp.plotOn(plotFrameWithNormRange, LineColor="r", Range="full", NormRange="left,right", LineStyle=10)
plotFrameWithNormRange.Draw()
canvas.Draw()
canvas.SaveAs("rf212_plottingInRanges_blinding.png")