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