ROOT
master
Reference Guide
Loading...
Searching...
No Matches
rf504_simwstool.py
Go to the documentation of this file.
1
## \file
2
## \ingroup tutorial_roofit_main
3
## \notebook -nodraw
4
## Organization and simultaneous fits: using RooSimWSTool to construct a simultaneous pdf
5
## that is built of variations of an input pdf
6
##
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 master pdf
17
# ---------------------------------
18
19
# Construct gauss(x,m,s)
20
x =
ROOT.RooRealVar
(
"x"
,
"x"
, -10, 10)
21
m =
ROOT.RooRealVar
(
"m"
,
"m"
, 0, -10, 10)
22
s =
ROOT.RooRealVar
(
"s"
,
"s"
, 1, -10, 10)
23
gauss =
ROOT.RooGaussian
(
"g"
,
"g"
, x, m, s)
24
25
# Construct poly(x,p0)
26
p0 =
ROOT.RooRealVar
(
"p0"
,
"p0"
, 0.01, 0.0, 1.0)
27
poly =
ROOT.RooPolynomial
(
"p"
,
"p"
, x, [p0])
28
29
# model = f*gauss(x) + (1-f)*poly(x)
30
f =
ROOT.RooRealVar
(
"f"
,
"f"
, 0.5, 0.0, 1.0)
31
model =
ROOT.RooAddPdf
(
"model"
,
"model"
, [gauss, poly], [f])
32
33
# Create category observables for splitting
34
# ----------------------------------------------------------------------------------
35
36
# Define two categories that can be used for splitting
37
c =
ROOT.RooCategory
(
"c"
,
"c"
)
38
c.defineType
(
"run1"
)
39
c.defineType
(
"run2"
)
40
41
d =
ROOT.RooCategory
(
"d"
,
"d"
)
42
d.defineType
(
"foo"
)
43
d.defineType
(
"bar"
)
44
45
# Set up SimWSTool
46
# -----------------------------
47
48
# Import ingredients in a workspace
49
w =
ROOT.RooWorkspace
(
"w"
,
"w"
)
50
w.Import
({model, c, d})
51
52
# Make Sim builder tool
53
sct =
ROOT.RooSimWSTool
(w)
54
55
# Build a simultaneous model with one split
56
# ---------------------------------------------------------------------------------
57
58
# Construct a simultaneous pdf with the following form
59
#
60
# model_run1(x) = f*gauss_run1(x,m_run1,s) + (1-f)*poly
61
# model_run2(x) = f*gauss_run2(x,m_run2,s) + (1-f)*poly
62
# simpdf(x,c) = model_run1(x) if c=="run1"
63
# = model_run2(x) if c=="run2"
64
#
65
# Returned pdf is owned by the workspace
66
model_sim =
sct.build
(
"model_sim"
,
"model"
, SplitParam=(
"m"
,
"c"
))
67
68
# Print tree structure of model
69
model_sim.Print
(
"t"
)
70
71
# Adjust model_sim parameters in workspace
72
w.var
(
"m_run1"
).setVal(-3)
73
w.var
(
"m_run2"
).setVal(+3)
74
75
# Print contents of workspace
76
w.Print
(
"v"
)
77
78
# Build a simultaneous model with product split
79
# -----------------------------------------------------------------------------------------
80
81
# Build another simultaneous pdf using a composite split in states c X d
82
model_sim2 =
sct.build
(
"model_sim2"
,
"model"
, SplitParam=(
"p0"
,
"c,d"
))
83
84
# Print tree structure of self model
85
model_sim2.Print
(
"t"
)
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
rf504_simwstool.py
ROOT master - Reference Guide Generated on Sat Feb 8 2025 04:21:36 (GVA Time) using Doxygen 1.10.0