Logo ROOT  
Reference Guide
rf103_interprfuncs.py
Go to the documentation of this file.
1 ## \file
2 ## \ingroup tutorial_roofit
3 ## \notebook
4 ##
5 ## Basic functionality: interpreted functions and p.d.f.s
6 ##
7 ## \macro_code
8 ##
9 ## \date February 2018
10 ## \authors Clemens Lange, Wouter Verkerke (C++ version)
11 
12 import ROOT
13 
14 # Generic interpreted p.d.f.
15 # ------------------------------
16 
17 # Declare observable x
18 x = ROOT.RooRealVar("x", "x", -20, 20)
19 
20 # Construct generic pdf from interpreted expression
21 # ------------------------------------------------------
22 
23 # ROOT.To construct a proper p.d.f, the formula expression is explicitly normalized internally by dividing
24 # it by a numeric integral of the expresssion over x in the range [-20,20]
25 #
26 alpha = ROOT.RooRealVar("alpha", "alpha", 5, 0.1, 10)
27 genpdf = ROOT.RooGenericPdf(
28  "genpdf",
29  "genpdf",
30  "(1+0.1*abs(x)+sin(sqrt(abs(x*alpha+0.1))))",
31  ROOT.RooArgList(
32  x,
33  alpha))
34 
35 # Sample, fit and plot generic pdf
36 # ---------------------------------------------------------------
37 
38 # Generate a toy dataset from the interpreted p.d.f
39 data = genpdf.generate(ROOT.RooArgSet(x), 10000)
40 
41 # Fit the interpreted p.d.f to the generated data
42 genpdf.fitTo(data)
43 
44 # Make a plot of the data and the p.d.f overlaid
45 xframe = x.frame(ROOT.RooFit.Title("Interpreted expression pdf"))
46 data.plotOn(xframe)
47 genpdf.plotOn(xframe)
48 
49 # Standard p.d.f. adjust with interpreted helper function
50 # ------------------------------------------------------------------------------------------------------------
51 # Make a gauss(x,sqrt(mean2),sigma) from a standard ROOT.RooGaussian #
52 #
53 # Construct standard pdf with formula replacing parameter
54 # ------------------------------------------------------------------------------------------------------------
55 
56 # Construct parameter mean2 and sigma
57 mean2 = ROOT.RooRealVar("mean2", "mean^2", 10, 0, 200)
58 sigma = ROOT.RooRealVar("sigma", "sigma", 3, 0.1, 10)
59 
60 # Construct interpreted function mean = sqrt(mean^2)
61 mean = ROOT.RooFormulaVar(
62  "mean", "mean", "sqrt(mean2)", ROOT.RooArgList(mean2))
63 
64 # Construct a gaussian g2(x,sqrt(mean2),sigma)
65 g2 = ROOT.RooGaussian("g2", "h2", x, mean, sigma)
66 
67 # Generate toy data
68 # ---------------------------------
69 
70 # Construct a separate gaussian g1(x,10,3) to generate a toy Gaussian
71 # dataset with mean 10 and width 3
72 g1 = ROOT.RooGaussian("g1", "g1", x, ROOT.RooFit.RooConst(
73  10), ROOT.RooFit.RooConst(3))
74 data2 = g1.generate(ROOT.RooArgSet(x), 1000)
75 
76 # Fit and plot tailored standard pdf
77 # -------------------------------------------------------------------
78 
79 # Fit g2 to data from g1
80 r = g2.fitTo(data2, ROOT.RooFit.Save()) # ROOT.RooFitResult
81 r.Print()
82 
83 # Plot data on frame and overlay projection of g2
84 xframe2 = x.frame(ROOT.RooFit.Title("Tailored Gaussian pdf"))
85 data2.plotOn(xframe2)
86 g2.plotOn(xframe2)
87 
88 # Draw all frames on a canvas
89 c = ROOT.TCanvas("rf103_interprfuncs", "rf103_interprfuncs", 800, 400)
90 c.Divide(2)
91 c.cd(1)
92 ROOT.gPad.SetLeftMargin(0.15)
93 xframe.GetYaxis().SetTitleOffset(1.4)
94 xframe.Draw()
95 c.cd(2)
96 ROOT.gPad.SetLeftMargin(0.15)
97 xframe2.GetYaxis().SetTitleOffset(1.4)
98 xframe2.Draw()
99 
100 c.SaveAs("rf103_interprfuncs.png")