Logo ROOT  
Reference Guide
rf402_datahandling.py
Go to the documentation of this file.
1 ## \file
2 ## \ingroup tutorial_roofit
3 ## \notebook
4 ##
5 ## Data and categories: tools for manipulation of (un)binned datasets
6 ##
7 ## \macro_code
8 ##
9 ## \date February 2018
10 ## \authors Clemens Lange, Wouter Verkerke (C++ version)
11 
12 from __future__ import print_function
13 import ROOT
14 import math
15 
16 # WVE Add reduction by range
17 
18 # Binned (RooDataHist) and unbinned datasets (RooDataSet) share
19 # many properties and inherit from a common abstract base class
20 # (RooAbsData), provides an interface for all operations
21 # that can be performed regardless of the data format
22 
23 x = ROOT.RooRealVar("x", "x", -10, 10)
24 y = ROOT.RooRealVar("y", "y", 0, 40)
25 c = ROOT.RooCategory("c", "c")
26 c.defineType("Plus", +1)
27 c.defineType("Minus", -1)
28 
29 # Basic operations on unbinned datasetss
30 # --------------------------------------------------------------
31 
32 # ROOT.RooDataSet is an unbinned dataset (a collection of points in
33 # N-dimensional space)
34 d = ROOT.RooDataSet("d", "d", ROOT.RooArgSet(x, y, c))
35 
36 # Unlike ROOT.RooAbsArgs (ROOT.RooAbsPdf, ROOT.RooFormulaVar,....) datasets are not attached to
37 # the variables they are constructed from. Instead they are attached to an internal
38 # clone of the supplied set of arguments
39 
40 # Fill d with dummy values
41 for i in range(1000):
42  x.setVal(i / 50 - 10)
43  y.setVal(math.sqrt(1.0 * i))
44  if (i % 2):
45  c.setLabel("Plus")
46  else:
47  c.setLabel("Minus")
48 
49  # We must explicitly refer to x,y, here to pass the values because
50  # d is not linked to them (as explained above)
51  if i < 3:
52  print(x, y, c)
53  print(type(x))
54  d.add(ROOT.RooArgSet(x, y, c))
55 
56 d.Print("v")
57 print("")
58 
59 # The get() function returns a pointer to the internal copy of the RooArgSet(x,y,c)
60 # supplied in the constructor
61 row = d.get()
62 row.Print("v")
63 print("")
64 
65 # Get with an argument loads a specific data point in row and returns
66 # a pointer to row argset. get() always returns the same pointer, unless
67 # an invalid row number is specified. In that case a null ptr is returned
68 d.get(900).Print("v")
69 print("")
70 
71 # Reducing, appending and merging
72 # -------------------------------------------------------------
73 
74 # The reduce() function returns a dataset which is a subset of the
75 # original
76 print("\n >> d1 has only columns x,c")
77 d1 = d.reduce(ROOT.RooArgSet(x, c))
78 d1.Print("v")
79 
80 print("\n >> d2 has only column y")
81 d2 = d.reduce(ROOT.RooArgSet(y))
82 d2.Print("v")
83 
84 print("\n >> d3 has only the points with y>5.17")
85 d3 = d.reduce("y>5.17")
86 d3.Print("v")
87 
88 print("\n >> d4 has only columns x, for data points with y>5.17")
89 d4 = d.reduce(ROOT.RooArgSet(x, c), "y>5.17")
90 d4.Print("v")
91 
92 # The merge() function adds two data set column-wise
93 print("\n >> merge d2(y) with d1(x,c) to form d1(x,c,y)")
94 d1.merge(d2)
95 d1.Print("v")
96 
97 # The append() function addes two datasets row-wise
98 print("\n >> append data points of d3 to d1")
99 d1.append(d3)
100 d1.Print("v")
101 
102 # Operations on binned datasets
103 # ---------------------------------------------------------
104 
105 # A binned dataset can be constructed empty, an unbinned dataset, or
106 # from a ROOT native histogram (TH1,2,3)
107 
108 print(">> construct dh (binned) from d(unbinned) but only take the x and y dimensions, ")
109 print(">> the category 'c' will be projected in the filling process")
110 
111 # The binning of real variables (like x,y) is done using their fit range
112 # 'get/setRange()' and number of specified fit bins 'get/setBins()'.
113 # Category dimensions of binned datasets get one bin per defined category
114 # state
115 x.setBins(10)
116 y.setBins(10)
117 dh = ROOT.RooDataHist("dh", "binned version of d", ROOT.RooArgSet(x, y), d)
118 dh.Print("v")
119 
120 yframe = y.frame(ROOT.RooFit.Bins(10), ROOT.RooFit.Title(
121  "Operations on binned datasets"))
122 dh.plotOn(yframe) # plot projection of 2D binned data on y
123 
124 # Examine the statistics of a binned dataset
125 print(">> number of bins in dh : ", dh.numEntries())
126 print(">> sum of weights in dh : ", dh.sum(ROOT.kFALSE))
127 # accounts for bin volume
128 print(">> integral over histogram: ", dh.sum(ROOT.kTRUE))
129 
130 # Locate a bin from a set of coordinates and retrieve its properties
131 x.setVal(0.3)
132 y.setVal(20.5)
133 print(">> retrieving the properties of the bin enclosing coordinate (x,y) = (0.3,20.5) bin center:")
134 # load bin center coordinates in internal buffer
135 dh.get(ROOT.RooArgSet(x, y)).Print("v")
136 print(" weight = ", dh.weight()) # return weight of last loaded coordinates
137 
138 # Reduce the 2-dimensional binned dataset to a 1-dimensional binned dataset
139 #
140 # All reduce() methods are interfaced in RooAbsData. All reduction techniques
141 # demonstrated on unbinned datasets can be applied to binned datasets as
142 # well.
143 print(">> Creating 1-dimensional projection on y of dh for bins with x>0")
144 dh2 = dh.reduce(ROOT.RooArgSet(y), "x>0")
145 dh2.Print("v")
146 
147 # Add dh2 to yframe and redraw
148 dh2.plotOn(yframe, ROOT.RooFit.LineColor(ROOT.kRed),
149  ROOT.RooFit.MarkerColor(ROOT.kRed))
150 
151 # Saving and loading from file
152 # -------------------------------------------------------
153 
154 # Datasets can be persisted with ROOT I/O
155 print("\n >> Persisting d via ROOT I/O")
156 f = ROOT.TFile("rf402_datahandling.root", "RECREATE")
157 d.Write()
158 f.ls()
159 
160 # To read back in future session:
161 # > ROOT.TFile f("rf402_datahandling.root")
162 # > d = (ROOT.RooDataSet*) f.FindObject("d")
163 
164 c = ROOT.TCanvas("rf402_datahandling", "rf402_datahandling", 600, 600)
165 ROOT.gPad.SetLeftMargin(0.15)
166 yframe.GetYaxis().SetTitleOffset(1.4)
167 yframe.Draw()
168 
169 c.SaveAs("rf402_datahandling.png")
ROOT::Math::IntegOptionsUtil::Print
void Print(std::ostream &os, const OptionType &opt)
Definition: IntegratorOptions.cxx:101
type
int type
Definition: TGX11.cxx:121