ROOT
Version v6.34
master
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
Reference Guide
►
ROOT
•
All
Classes
Namespaces
Files
Functions
Variables
Typedefs
Enumerations
Enumerator
Properties
Friends
Macros
Modules
Pages
Loading...
Searching...
No Matches
rf510_wsnamedsets.py
Go to the documentation of this file.
1
## \file
2
## \ingroup tutorial_roofit
3
## \notebook
4
## 'ORGANIZATION AND SIMULTANEOUS FITS' RooFit tutorial macro #510
5
##
6
## Working with named parameter sets and parameter snapshots in
7
## workspaces
8
##
9
## \macro_image
10
## \macro_code
11
## \macro_output
12
##
13
## \date February 2018
14
## \authors Clemens Lange, Wouter Verkerke (C version)
15
16
17
import
ROOT
18
19
20
def
fillWorkspace
(w):
21
# Create model
22
# -----------------------
23
24
# Declare observable x
25
x =
ROOT.RooRealVar
(
"x"
,
"x"
, 0, 10)
26
27
# Create two Gaussian PDFs g1(x,mean1,sigma) anf g2(x,mean2,sigma) and
28
# their parameters
29
mean =
ROOT.RooRealVar
(
"mean"
,
"mean of gaussians"
, 5, 0, 10)
30
sigma1 =
ROOT.RooRealVar
(
"sigma1"
,
"width of gaussians"
, 0.5)
31
sigma2 =
ROOT.RooRealVar
(
"sigma2"
,
"width of gaussians"
, 1)
32
33
sig1 =
ROOT.RooGaussian
(
"sig1"
,
"Signal component 1"
, x, mean, sigma1)
34
sig2 =
ROOT.RooGaussian
(
"sig2"
,
"Signal component 2"
, x, mean, sigma2)
35
36
# Build Chebychev polynomial p.d.f.
37
a0 =
ROOT.RooRealVar
(
"a0"
,
"a0"
, 0.5, 0.0, 1.0)
38
a1 =
ROOT.RooRealVar
(
"a1"
,
"a1"
, -0.2, 0.0, 1.0)
39
bkg =
ROOT.RooChebychev
(
"bkg"
,
"Background"
, x, [a0, a1])
40
41
# Sum the signal components into a composite signal p.d.f.
42
sig1frac =
ROOT.RooRealVar
(
"sig1frac"
,
"fraction of component 1 in signal"
, 0.8, 0.0, 1.0)
43
sig =
ROOT.RooAddPdf
(
"sig"
,
"Signal"
, [sig1, sig2], [sig1frac])
44
45
# Sum the composite signal and background
46
bkgfrac =
ROOT.RooRealVar
(
"bkgfrac"
,
"fraction of background"
, 0.5, 0.0, 1.0)
47
model =
ROOT.RooAddPdf
(
"model"
,
"g1+g2+a"
, [bkg, sig], [bkgfrac])
48
49
# Import model into p.d.f.
50
w.Import
(model)
51
52
# Encode definition of parameters in workspace
53
# ---------------------------------------------------------------------------------------
54
55
# Define named sets "parameters" and "observables", list which variables should be considered
56
# parameters and observables by the users convention
57
#
58
# Variables appearing in sets _must_ live in the workspace already, the autoImport flag
59
# of defineSet must be set to import them on the fly. Named sets contain only references
60
# to the original variables, the value of observables in named sets already
61
# reflect their 'current' value
62
params =
model.getParameters
({x})
63
w.defineSet
(
"parameters"
, params)
64
w.defineSet
(
"observables"
, {x})
65
66
# Encode reference value for parameters in workspace
67
# ---------------------------------------------------------------------------------------------------
68
69
# Define a parameter 'snapshot' in the p.d.f.
70
# Unlike a named set, parameter snapshot stores an independent set of values for
71
# a given set of variables in the workspace. The values can be stored and reloaded
72
# into the workspace variable objects using the loadSnapshot() and saveSnapshot()
73
# methods. A snapshot saves the value of each variable, errors that are stored
74
# with it as well as the 'Constant' flag that is used in fits to determine if a
75
# parameter is kept fixed or not.
76
77
# Do a dummy fit to a (supposedly) reference dataset here and store the results
78
# of that fit into a snapshot
79
refData =
model.generate
({x}, 10000)
80
model.fitTo
(refData, PrintLevel=-1)
81
82
# The kTRUE flag imports the values of the objects in (*params) into the workspace
83
# If not set, present values of the workspace parameters objects are stored
84
w.saveSnapshot
(
"reference_fit"
, params,
True
)
85
86
# Make another fit with the signal componentforced to zero
87
# and save those parameters too
88
89
bkgfrac.setVal
(1)
90
bkgfrac.setConstant
(
True
)
91
bkgfrac.removeError
()
92
model.fitTo
(refData, PrintLevel=-1)
93
94
w.saveSnapshot
(
"reference_fit_bkgonly"
, params,
True
)
95
96
97
# Create model and dataset
98
# -----------------------------------------------
99
100
w =
ROOT.RooWorkspace
(
"w"
)
101
fillWorkspace
(w)
102
103
# Exploit convention encoded in named set "parameters" and "observables"
104
# to use workspace contents w/o need for introspected
105
model = w[
"model"
]
106
107
# Generate data from p.d.f. in given observables
108
data =
model.generate
(
w.set
(
"observables"
), 1000)
109
110
# Fit model to data
111
model.fitTo
(data, PrintLevel=-1)
112
113
# Plot fitted model and data on frame of first (only) observable
114
frame = (
w.set
(
"observables"
).first()).frame()
115
data.plotOn
(frame)
116
model.plotOn
(frame)
117
118
# Overlay plot with model with reference parameters as stored in snapshots
119
w.loadSnapshot
(
"reference_fit"
)
120
model.plotOn
(frame, LineColor=
"r"
)
121
w.loadSnapshot
(
"reference_fit_bkgonly"
)
122
model.plotOn
(frame, LineColor=
"r"
, LineStyle=
"--"
)
123
124
# Draw the frame on the canvas
125
c =
ROOT.TCanvas
(
"rf510_wsnamedsets"
,
"rf503_wsnamedsets"
, 600, 600)
126
ROOT.gPad.SetLeftMargin
(0.15)
127
frame.GetYaxis
().SetTitleOffset(1.4)
128
frame.Draw
()
129
130
c.SaveAs
(
"rf510_wsnamedsets.png"
)
131
132
# Print workspace contents
133
w.Print
()
134
135
# Workspace will remain in memory after macro finishes
136
ROOT.gDirectory.Add
(w)
TRangeDynCast
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
Definition
TCollection.h:358
ROOT::Detail::TRangeCast
Definition
TCollection.h:311
tutorials
roofit
rf510_wsnamedsets.py
ROOT tags/6-34-04 - Reference Guide Generated on Wed Mar 19 2025 04:35:07 (GVA Time) using Doxygen 1.10.0