Logo ROOT  
Reference Guide
rf105_funcbinding.py
Go to the documentation of this file.
1## \file
2## \ingroup tutorial_roofit
3## \notebook
4## 'BASIC FUNCTIONALITY' RooFit tutorial macro #105
5## Demonstration of binding ROOT Math functions as RooFit functions
6## and pdfs
7##
8## \macro_code
9##
10## \date February 2018
11## \authors Clemens Lange, Wouter Verkerke (C version)
12
13import ROOT
14
15# Bind ROOT TMath::Erf C function
16# ---------------------------------------------------
17
18# Bind one-dimensional ROOT.TMath.Erf function as ROOT.RooAbsReal function
19x = ROOT.RooRealVar("x", "x", -3, 3)
20erf = ROOT.RooFit.bindFunction("erf", ROOT.TMath.Erf, x)
21
22# Print erf definition
23erf.Print()
24
25# Plot erf on frame
26frame1 = x.frame(Title="TMath.Erf bound as ROOT.RooFit function")
27erf.plotOn(frame1)
28
29# Bind ROOT::Math::beta_pdf C function
30# -----------------------------------------------------------------------
31
32# Bind pdf ROOT.Math.Beta with three variables as ROOT.RooAbsPdf function
33x2 = ROOT.RooRealVar("x2", "x2", 0, 0.999)
34a = ROOT.RooRealVar("a", "a", 5, 0, 10)
35b = ROOT.RooRealVar("b", "b", 2, 0, 10)
36beta = ROOT.RooFit.bindPdf("beta", ROOT.Math.beta_pdf, x2, a, b)
37
38# Perf beta definition
39beta.Print()
40
41# Generate some events and fit
42data = beta.generate({x2}, 10000)
43beta.fitTo(data)
44
45# Plot data and pdf on frame
46frame2 = x2.frame(Title="ROOT.Math.Beta bound as ROOT.RooFit pdf")
47data.plotOn(frame2)
48beta.plotOn(frame2)
49
50# Bind ROOT TF1 as RooFit function
51# ---------------------------------------------------------------
52
53# Create a ROOT TF1 function
54fa1 = ROOT.TF1("fa1", "sin(x)/x", 0, 10)
55
56# Create an observable
57x3 = ROOT.RooRealVar("x3", "x3", 0.01, 20)
58
59# Create binding of TF1 object to above observable
60rfa1 = ROOT.RooFit.bindFunction(fa1, x3)
61
62# Print rfa1 definition
63rfa1.Print()
64
65# Make plot frame in observable, TF1 binding function
66frame3 = x3.frame(Title="TF1 bound as ROOT.RooFit function")
67rfa1.plotOn(frame3)
68
69c = ROOT.TCanvas("rf105_funcbinding", "rf105_funcbinding", 1200, 400)
70c.Divide(3)
71c.cd(1)
72ROOT.gPad.SetLeftMargin(0.15)
73frame1.GetYaxis().SetTitleOffset(1.6)
74frame1.Draw()
75c.cd(2)
76ROOT.gPad.SetLeftMargin(0.15)
77frame2.GetYaxis().SetTitleOffset(1.6)
78frame2.Draw()
79c.cd(3)
80ROOT.gPad.SetLeftMargin(0.15)
81frame3.GetYaxis().SetTitleOffset(1.6)
82frame3.Draw()
83
84c.SaveAs("rf105_funcbinding.png")