ROOT
Version master
master
v6.34
v6.32
v6.30
v6.28
v6.26
v6.24
v6.22
v6.20
v6.18
v6.16
v6.14
v6.12
v6.10
v6.08
v6.06
v6.04
Reference Guide
►
ROOT
•
All
Classes
Namespaces
Files
Functions
Variables
Typedefs
Enumerations
Enumerator
Properties
Friends
Macros
Modules
Pages
Loading...
Searching...
No Matches
rf315_projectpdf.py
Go to the documentation of this file.
1
## \file
2
## \ingroup tutorial_roofit_main
3
## \notebook
4
## Multidimensional models: marginizalization of multi-dimensional pdfs through integration
5
##
6
## \macro_image
7
## \macro_code
8
## \macro_output
9
##
10
## \date February 2018
11
## \authors Clemens Lange, Wouter Verkerke (C++ version)
12
13
import
ROOT
14
15
16
# Create pdf m(x,y) = gx(x|y) * g(y)
17
# --------------------------------------------------------------
18
19
# Increase default precision of numeric integration
20
# as self exercise has high sensitivity to numeric integration precision
21
ROOT.RooAbsPdf.defaultIntegratorConfig
().setEpsRel(1e-8)
22
ROOT.RooAbsPdf.defaultIntegratorConfig
().setEpsAbs(1e-8)
23
24
# Create observables
25
x =
ROOT.RooRealVar
(
"x"
,
"x"
, -5, 5)
26
y =
ROOT.RooRealVar
(
"y"
,
"y"
, -2, 2)
27
28
# Create function f(y) = a0 + a1*y
29
a0 =
ROOT.RooRealVar
(
"a0"
,
"a0"
, 0)
30
a1 =
ROOT.RooRealVar
(
"a1"
,
"a1"
, -1.5, -3, 1)
31
fy =
ROOT.RooPolyVar
(
"fy"
,
"fy"
, y, [a0, a1])
32
33
# Create gaussx(x,f(y),sx)
34
sigmax =
ROOT.RooRealVar
(
"sigmax"
,
"width of gaussian"
, 0.5)
35
gaussx =
ROOT.RooGaussian
(
"gaussx"
,
"Gaussian in x with shifting mean in y"
, x, fy, sigmax)
36
37
# Create gaussy(y,0,2)
38
gaussy =
ROOT.RooGaussian
(
"gaussy"
,
"Gaussian in y"
, y, 0.0, 2.0)
39
40
# Create gaussx(x,sx|y) * gaussy(y)
41
model =
ROOT.RooProdPdf
(
42
"model"
,
43
"gaussx(x|y)*gaussy(y)"
,
44
{gaussy},
45
Conditional=({gaussx}, {x}),
46
)
47
48
# Marginalize m(x,y) to m(x)
49
# ----------------------------------------------------
50
51
# modelx(x) = Int model(x,y) dy
52
modelx =
model.createProjection
({y})
53
54
# Use marginalized pdf as regular 1D pdf
55
# -----------------------------------------------
56
57
# Sample 1000 events from modelx
58
data =
modelx.generateBinned
({x}, 1000)
59
60
# Fit modelx to toy data
61
modelx.fitTo
(data, Verbose=
True
, PrintLevel=-1)
62
63
# Plot modelx over data
64
frame =
x.frame
(40)
65
data.plotOn
(frame)
66
modelx.plotOn
(frame)
67
68
# Make 2D histogram of model(x,y)
69
hh =
model.createHistogram
(
"x,y"
)
70
hh.SetLineColor
(
"kBlue"
)
71
72
c =
ROOT.TCanvas
(
"rf315_projectpdf"
,
"rf315_projectpdf"
, 800, 400)
73
c.Divide
(2)
74
c.cd
(1)
75
ROOT.gPad.SetLeftMargin
(0.15)
76
frame.GetYaxis
().SetTitleOffset(1.4)
77
frame.Draw
()
78
c.cd
(2)
79
ROOT.gPad.SetLeftMargin
(0.20)
80
hh.GetZaxis
().SetTitleOffset(2.5)
81
hh.Draw
(
"surf"
)
82
c.SaveAs
(
"rf315_projectpdf.png"
)
TRangeDynCast
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
Definition
TCollection.h:358
tutorials
roofit
roofit
rf315_projectpdf.py
ROOT master - Reference Guide Generated on Fri Mar 7 2025 21:57:58 (GVA Time) using Doxygen 1.10.0