Logo ROOT  
Reference Guide
rf706_histpdf.py
Go to the documentation of this file.
1 ## \file
2 ## \ingroup tutorial_roofit
3 ## \notebook
4 ## Special pdf's: histogram based pdfs and functions
5 ##
6 ## \macro_code
7 ##
8 ## \date February 2018
9 ## \authors Clemens Lange, Wouter Verkerke (C++ version)
10 
11 import ROOT
12 
13 
14 # Create pdf for sampling
15 # ---------------------------------------------
16 
17 x = ROOT.RooRealVar("x", "x", 0, 20)
18 p = ROOT.RooPolynomial("p", "p", x, ROOT.RooArgList(ROOT.RooFit.RooConst(
19  0.01), ROOT.RooFit.RooConst(-0.01), ROOT.RooFit.RooConst(0.0004)))
20 
21 # Create low stats histogram
22 # ---------------------------------------------------
23 
24 # Sample 500 events from p
25 x.setBins(20)
26 data1 = p.generate(ROOT.RooArgSet(x), 500)
27 
28 # Create a binned dataset with 20 bins and 500 events
29 hist1 = data1.binnedClone()
30 
31 # Represent data in dh as pdf in x
32 histpdf1 = ROOT.RooHistPdf("histpdf1", "histpdf1", ROOT.RooArgSet(x), hist1, 0)
33 
34 # Plot unbinned data and histogram pdf overlaid
35 frame1 = x.frame(ROOT.RooFit.Title(
36  "Low statistics histogram pdf"), ROOT.RooFit.Bins(100))
37 data1.plotOn(frame1)
38 histpdf1.plotOn(frame1)
39 
40 # Create high stats histogram
41 # -----------------------------------------------------
42 
43 # Sample 100000 events from p
44 x.setBins(10)
45 data2 = p.generate(ROOT.RooArgSet(x), 100000)
46 
47 # Create a binned dataset with 10 bins and 100K events
48 hist2 = data2.binnedClone()
49 
50 # Represent data in dh as pdf in x, 2nd order interpolation
51 histpdf2 = ROOT.RooHistPdf("histpdf2", "histpdf2", ROOT.RooArgSet(x), hist2, 2)
52 
53 # Plot unbinned data and histogram pdf overlaid
54 frame2 = x.frame(ROOT.RooFit.Title(
55  "High stats histogram pdf with interpolation"), ROOT.RooFit.Bins(100))
56 data2.plotOn(frame2)
57 histpdf2.plotOn(frame2)
58 
59 c = ROOT.TCanvas("rf706_histpdf", "rf706_histpdf", 800, 400)
60 c.Divide(2)
61 c.cd(1)
62 ROOT.gPad.SetLeftMargin(0.15)
63 frame1.GetYaxis().SetTitleOffset(1.4)
64 frame1.Draw()
65 c.cd(2)
66 ROOT.gPad.SetLeftMargin(0.15)
67 frame2.GetYaxis().SetTitleOffset(1.8)
68 frame2.Draw()
69 
70 c.SaveAs("rf706_histpdf.png")