ROOT
master
Reference Guide
Loading...
Searching...
No Matches
rf102_dataimport.py
Go to the documentation of this file.
1
## \file
2
## \ingroup tutorial_roofit_main
3
## \notebook
4
## 'BASIC FUNCTIONALITY' RooFit tutorial macro #102
5
## Importing data from ROOT TTrees and THx histograms
6
##
7
## \macro_image
8
## \macro_code
9
## \macro_output
10
##
11
## \date February 2018
12
## \authors Clemens Lange, Wouter Verkerke (C version)
13
14
import
ROOT
15
from
array
import
array
16
import
numpy
as
np
17
18
19
def
makeTH1
(trnd):
20
21
# Create ROOT ROOT.TH1 filled with a Gaussian distribution
22
23
hh =
ROOT.TH1D
(
"hh"
,
"hh"
, 25, -10, 10)
24
hh.Fill
(
np.array
([
trnd.Gaus
(0, 3)
for
_
in
range
(100)]))
25
return
hh
26
27
28
def
makeTTree
(trnd):
29
# Create ROOT ROOT.TTree filled with a Gaussian distribution in x and a
30
# uniform distribution in y
31
32
tree =
ROOT.TTree
(
"tree"
,
"tree"
)
33
px = array(
"d"
, [0])
34
py = array(
"d"
, [0])
35
tree.Branch
(
"x"
, px,
"x/D"
)
36
tree.Branch
(
"y"
, py,
"y/D"
)
37
for
i
in
range
(100):
38
px[0] =
trnd.Gaus
(0, 3)
39
py[0] =
trnd.Uniform
() * 30 - 15
40
tree.Fill
()
41
return
tree
42
43
trnd =
ROOT.TRandom3
()
44
45
############################
46
# Importing ROOT histograms
47
############################
48
# Import ROOT TH1 into a RooDataHist
49
# ---------------------------------------------------------
50
# Create a ROOT TH1 histogram
51
hh =
makeTH1
(trnd)
52
53
# Declare observable x
54
x =
ROOT.RooRealVar
(
"x"
,
"x"
, -10, 10)
55
56
# Create a binned dataset that imports contents of ROOT.TH1 and associates
57
# its contents to observable 'x'
58
dh =
ROOT.RooDataHist
(
"dh"
,
"dh"
, [x], Import=hh)
59
60
# Plot and fit a RooDataHist
61
# ---------------------------------------------------
62
# Make plot of binned dataset showing Poisson error bars (RooFit default)
63
frame =
x.frame
(Title=
"Imported ROOT.TH1 with Poisson error bars"
)
64
dh.plotOn
(frame)
65
66
# Fit a Gaussian p.d.f to the data
67
mean =
ROOT.RooRealVar
(
"mean"
,
"mean"
, 0, -10, 10)
68
sigma =
ROOT.RooRealVar
(
"sigma"
,
"sigma"
, 3, 0.1, 10)
69
gauss =
ROOT.RooGaussian
(
"gauss"
,
"gauss"
, x, mean, sigma)
70
gauss.fitTo
(dh, PrintLevel=-1)
71
gauss.plotOn
(frame)
72
73
# Plot and fit a RooDataHist with internal errors
74
# ---------------------------------------------------------------------------------------------
75
76
# If histogram has custom error (i.e. its contents is does not originate from a Poisson process
77
# but e.g. is a sum of weighted events) you can data with symmetric 'sum-of-weights' error instead
78
# (same error bars as shown by ROOT)
79
frame2 =
x.frame
(Title=
"Imported ROOT.TH1 with internal errors"
)
80
dh.plotOn
(frame2, DataError=
"SumW2"
)
81
gauss.plotOn
(frame2)
82
83
# Please note that error bars shown (Poisson or SumW2) are for visualization only, the are NOT used
84
# in a maximum likelihood fit
85
#
86
# A (binned) ML fit will ALWAYS assume the Poisson error interpretation of data (the mathematical definition
87
# of likelihood does not take any external definition of errors). Data with non-unit weights can only be correctly
88
# fitted with a chi^2 fit (see rf602_chi2fit.py)
89
#
90
# Importing ROOT TTrees
91
# -----------------------------------------------------------
92
# Import ROOT TTree into a RooDataSet
93
94
tree =
makeTTree
(trnd)
95
96
# Define 2nd observable y
97
y =
ROOT.RooRealVar
(
"y"
,
"y"
, -10, 10)
98
99
# Construct unbinned dataset importing tree branches x and y matching between branches and ROOT.RooRealVars
100
# is done by name of the branch/RRV
101
#
102
# Note that ONLY entries for which x,y have values within their allowed ranges as defined in
103
# ROOT.RooRealVar x and y are imported. Since the y values in the import tree are in the range [-15,15]
104
# and RRV y defines a range [-10,10] this means that the ROOT.RooDataSet
105
# below will have less entries than the ROOT.TTree 'tree'
106
107
ds =
ROOT.RooDataSet
(
"ds"
,
"ds"
, {x, y}, Import=tree)
108
109
# Use ascii import/export for datasets
110
# ------------------------------------------------------------------------------------
111
112
113
def
write_dataset
(ds, filename):
114
# Write data to output stream
115
outstream =
ROOT.std.ofstream
(filename)
116
# Optionally, adjust the stream here (e.g. std::setprecision)
117
ds.write
(outstream)
118
outstream.close
()
119
120
121
write_dataset
(ds,
"rf102_testData.txt"
)
122
123
# Read data from input stream. The variables of the dataset need to be supplied
124
# to the RooDataSet::read() function.
125
print(
"\n-----------------------\nReading data from ASCII"
)
126
dataReadBack =
ROOT.RooDataSet.read
(
127
"rf102_testData.txt"
,
128
[x, y],
# variables to be read. If the file has more fields, these are ignored.
129
"D"
,
# Prints if a RooFit message stream listens for debug messages. Use Q for quiet.
130
)
131
132
dataReadBack.Print
(
"V"
)
133
134
print(
"\nOriginal data, line 20:"
)
135
ds.get
(20).
Print
(
"V"
)
136
137
print(
"\nRead-back data, line 20:"
)
138
dataReadBack.get
(20).
Print
(
"V"
)
139
140
141
# Plot data set with multiple binning choices
142
# ------------------------------------------------------------------------------------
143
# Print number of events in dataset
144
ds.Print
()
145
146
# Print unbinned dataset with default frame binning (100 bins)
147
frame3 =
y.frame
(Title=
"Unbinned data shown in default frame binning"
)
148
ds.plotOn
(frame3)
149
150
# Print unbinned dataset with custom binning choice (20 bins)
151
frame4 =
y.frame
(Title=
"Unbinned data shown with custom binning"
)
152
ds.plotOn
(frame4, Binning=20)
153
154
frame5 =
y.frame
(Title=
"Unbinned data read back from ASCII file"
)
155
ds.plotOn
(frame5, Binning=20)
156
dataReadBack.plotOn
(frame5, Binning=20, MarkerColor=
"r"
, MarkerStyle=5)
157
158
# Draw all frames on a canvas
159
c =
ROOT.TCanvas
(
"rf102_dataimport"
,
"rf102_dataimport"
, 800, 800)
160
c.Divide
(3, 2)
161
c.cd
(1)
162
ROOT.gPad.SetLeftMargin
(0.15)
163
frame.GetYaxis
().SetTitleOffset(1.4)
164
frame.Draw
()
165
c.cd
(2)
166
ROOT.gPad.SetLeftMargin
(0.15)
167
frame2.GetYaxis
().SetTitleOffset(1.4)
168
frame2.Draw
()
169
c.cd
(4)
170
ROOT.gPad.SetLeftMargin
(0.15)
171
frame3.GetYaxis
().SetTitleOffset(1.4)
172
frame3.Draw
()
173
c.cd
(5)
174
ROOT.gPad.SetLeftMargin
(0.15)
175
frame4.GetYaxis
().SetTitleOffset(1.4)
176
frame4.Draw
()
177
c.cd
(6)
178
ROOT.gPad.SetLeftMargin
(0.15)
179
frame4.GetYaxis
().SetTitleOffset(1.4)
180
frame5.Draw
()
181
182
c.SaveAs
(
"rf102_dataimport.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
Print
void Print(GNN_Data &d, std::string txt="")
Definition
TMVA_SOFIE_GNN_Application.C:59
ROOT::Detail::TRangeCast
Definition
TCollection.h:311
tutorials
roofit
roofit
rf102_dataimport.py
ROOT master - Reference Guide Generated on Fri Feb 14 2025 04:17:31 (GVA Time) using Doxygen 1.10.0