Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
rf102_dataimport.py File Reference

Detailed Description

View in nbviewer Open in SWAN
'BASIC FUNCTIONALITY' RooFit tutorial macro #102 Importing data from ROOT TTrees and THx histograms

import ROOT
from array import array
import numpy as np
def makeTH1(trnd):
# Create ROOT ROOT.TH1 filled with a Gaussian distribution
hh = ROOT.TH1D("hh", "hh", 25, -10, 10)
hh.Fill(np.array([trnd.Gaus(0, 3) for _ in range(100)]))
return hh
def makeTTree(trnd):
# Create ROOT ROOT.TTree filled with a Gaussian distribution in x and a
# uniform distribution in y
tree = ROOT.TTree("tree", "tree")
px = array("d", [0])
py = array("d", [0])
tree.Branch("x", px, "x/D")
tree.Branch("y", py, "y/D")
for i in range(100):
px[0] = trnd.Gaus(0, 3)
py[0] = trnd.Uniform() * 30 - 15
tree.Fill()
return tree
trnd = ROOT.TRandom3()
############################
# Importing ROOT histograms
############################
# Import ROOT TH1 into a RooDataHist
# ---------------------------------------------------------
# Create a ROOT TH1 histogram
hh = makeTH1(trnd)
# Declare observable x
x = ROOT.RooRealVar("x", "x", -10, 10)
# Create a binned dataset that imports contents of ROOT.TH1 and associates
# its contents to observable 'x'
dh = ROOT.RooDataHist("dh", "dh", [x], Import=hh)
# Plot and fit a RooDataHist
# ---------------------------------------------------
# Make plot of binned dataset showing Poisson error bars (RooFit default)
frame = x.frame(Title="Imported ROOT.TH1 with Poisson error bars")
dh.plotOn(frame)
# Fit a Gaussian p.d.f to the data
mean = ROOT.RooRealVar("mean", "mean", 0, -10, 10)
sigma = ROOT.RooRealVar("sigma", "sigma", 3, 0.1, 10)
gauss = ROOT.RooGaussian("gauss", "gauss", x, mean, sigma)
gauss.fitTo(dh, PrintLevel=-1)
gauss.plotOn(frame)
# Plot and fit a RooDataHist with internal errors
# ---------------------------------------------------------------------------------------------
# If histogram has custom error (i.e. its contents is does not originate from a Poisson process
# but e.g. is a sum of weighted events) you can data with symmetric 'sum-of-weights' error instead
# (same error bars as shown by ROOT)
frame2 = x.frame(Title="Imported ROOT.TH1 with internal errors")
dh.plotOn(frame2, DataError="SumW2")
gauss.plotOn(frame2)
# Please note that error bars shown (Poisson or SumW2) are for visualization only, the are NOT used
# in a maximum likelihood fit
#
# A (binned) ML fit will ALWAYS assume the Poisson error interpretation of data (the mathematical definition
# of likelihood does not take any external definition of errors). Data with non-unit weights can only be correctly
# fitted with a chi^2 fit (see rf602_chi2fit.py)
#
# Importing ROOT TTrees
# -----------------------------------------------------------
# Import ROOT TTree into a RooDataSet
tree = makeTTree(trnd)
# Define 2nd observable y
y = ROOT.RooRealVar("y", "y", -10, 10)
# Construct unbinned dataset importing tree branches x and y matching between branches and ROOT.RooRealVars
# is done by name of the branch/RRV
#
# Note that ONLY entries for which x,y have values within their allowed ranges as defined in
# ROOT.RooRealVar x and y are imported. Since the y values in the import tree are in the range [-15,15]
# and RRV y defines a range [-10,10] this means that the ROOT.RooDataSet
# below will have less entries than the ROOT.TTree 'tree'
ds = ROOT.RooDataSet("ds", "ds", {x, y}, Import=tree)
# Use ascii import/export for datasets
# ------------------------------------------------------------------------------------
def write_dataset(ds, filename):
# Write data to output stream
outstream = ROOT.std.ofstream(filename)
# Optionally, adjust the stream here (e.g. std::setprecision)
ds.write(outstream)
outstream.close()
write_dataset(ds, "rf102_testData.txt")
# Read data from input stream. The variables of the dataset need to be supplied
# to the RooDataSet::read() function.
print("\n-----------------------\nReading data from ASCII")
dataReadBack = ROOT.RooDataSet.read(
"rf102_testData.txt",
[x, y], # variables to be read. If the file has more fields, these are ignored.
"D", # Prints if a RooFit message stream listens for debug messages. Use Q for quiet.
)
dataReadBack.Print("V")
print("\nOriginal data, line 20:")
ds.get(20).Print("V")
print("\nRead-back data, line 20:")
dataReadBack.get(20).Print("V")
# Plot data set with multiple binning choices
# ------------------------------------------------------------------------------------
# Print number of events in dataset
ds.Print()
# Print unbinned dataset with default frame binning (100 bins)
frame3 = y.frame(Title="Unbinned data shown in default frame binning")
ds.plotOn(frame3)
# Print unbinned dataset with custom binning choice (20 bins)
frame4 = y.frame(Title="Unbinned data shown with custom binning")
ds.plotOn(frame4, Binning=20)
frame5 = y.frame(Title="Unbinned data read back from ASCII file")
ds.plotOn(frame5, Binning=20)
dataReadBack.plotOn(frame5, Binning=20, MarkerColor="r", MarkerStyle=5)
# Draw all frames on a canvas
c = ROOT.TCanvas("rf102_dataimport", "rf102_dataimport", 800, 800)
c.Divide(3, 2)
c.cd(1)
ROOT.gPad.SetLeftMargin(0.15)
frame.GetYaxis().SetTitleOffset(1.4)
frame.Draw()
c.cd(2)
ROOT.gPad.SetLeftMargin(0.15)
frame2.GetYaxis().SetTitleOffset(1.4)
frame2.Draw()
c.cd(4)
ROOT.gPad.SetLeftMargin(0.15)
frame3.GetYaxis().SetTitleOffset(1.4)
frame3.Draw()
c.cd(5)
ROOT.gPad.SetLeftMargin(0.15)
frame4.GetYaxis().SetTitleOffset(1.4)
frame4.Draw()
c.cd(6)
ROOT.gPad.SetLeftMargin(0.15)
frame4.GetYaxis().SetTitleOffset(1.4)
frame5.Draw()
c.SaveAs("rf102_dataimport.png")
void Print(GNN_Data &d, std::string txt="")
[#1] INFO:Fitting -- RooAbsPdf::fitTo(gauss_over_gauss_Int[x]) fixing normalization set for coefficient determination to observables in data
[#1] INFO:Fitting -- using CPU computation library compiled with -mavx2
[#1] INFO:Fitting -- RooAddition::defaultErrorLevel(nll_gauss_over_gauss_Int[x]_dh) Summation contains a RooNLLVar, using its error level
[#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: activating const optimization
[#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: deactivating const optimization
[#1] INFO:DataHandling -- RooTreeDataStore::loadValues(ds) Skipping event #0 because y cannot accommodate the value 14.424
[#1] INFO:DataHandling -- RooTreeDataStore::loadValues(ds) Skipping event #3 because y cannot accommodate the value -12.0022
[#1] INFO:DataHandling -- RooTreeDataStore::loadValues(ds) Skipping event #5 because y cannot accommodate the value 13.8261
[#1] INFO:DataHandling -- RooTreeDataStore::loadValues(ds) Skipping event #6 because y cannot accommodate the value -14.9925
[#1] INFO:DataHandling -- RooTreeDataStore::loadValues(ds) Skipping ...
[#0] WARNING:DataHandling -- RooTreeDataStore::loadValues(ds) Ignored 36 out-of-range events
[#1] INFO:DataHandling -- RooDataSet::read: reading file rf102_testData.txt
[#1] INFO:DataHandling -- RooDataSet::read: read 64 events (ignored 0 out of range events)
DataStore dataset (rf102_testData.txt)
Contains 64 entries
Observables:
1) x = 0.0174204 L(-10 - 10) "x"
2) y = 9.46654 L(-10 - 10) "y"
3) blindState = Normal(idx = 0)
"Blinding State"
1) RooRealVar:: x = -0.79919
2) RooRealVar:: y = 0.0106407
1) RooRealVar:: x = -0.79919
2) RooRealVar:: y = 0.0106407
3) RooCategory:: blindState = Normal(idx = 0)
RooDataSet::ds[x,y] = 64 entries
-----------------------
Reading data from ASCII
Original data, line 20:
Read-back data, line 20:
Date
February 2018
Authors
Clemens Lange, Wouter Verkerke (C version)

Definition in file rf102_dataimport.py.