Logo ROOT   6.14/05
Reference Guide
RDataFrame.cxx
Go to the documentation of this file.
1 // Author: Enrico Guiraud, Danilo Piparo CERN 12/2016
2 
3 /*************************************************************************
4  * Copyright (C) 1995-2016, Rene Brun and Fons Rademakers. *
5  * All rights reserved. *
6  * *
7  * For the licensing terms see $ROOTSYS/LICENSE. *
8  * For the list of contributors see $ROOTSYS/README/CREDITS. *
9  *************************************************************************/
10 
11 #include <algorithm>
12 #include <stdexcept>
13 
14 #include "ROOT/RDataFrame.hxx"
15 #include "TChain.h"
16 #include "TDirectory.h"
17 
18 // clang-format off
19 /**
20 * \class ROOT::RDataFrame
21 * \ingroup dataframe
22 * \brief ROOT's RDataFrame offers a high level interface for analyses of data stored in `TTree`s, CSV's and other data formats.
23 
24 In addition, multi-threading and other low-level optimisations allow users to exploit all the resources available
25 on their machines completely transparently.<br>
26 Skip to the [class reference](#reference) or keep reading for the user guide.
27 
28 In a nutshell:
29 ~~~{.cpp}
30 ROOT::EnableImplicitMT(); // Tell ROOT you want to go parallel
31 ROOT::RDataFrame d("myTree", "file_*.root"); // Interface to TTree and TChain
32 auto myHisto = d.Histo1D("Branch_A"); // This happens in parallel!
33 myHisto->Draw();
34 ~~~
35 
36 Calculations are expressed in terms of a type-safe *functional chain of actions and transformations*, `RDataFrame` takes
37 care of their execution. The implementation automatically puts in place several low level optimisations such as
38 multi-thread parallelisation and caching.
39 
40 \htmlonly
41 <a href="https://doi.org/10.5281/zenodo.260230"><img src="https://zenodo.org/badge/DOI/10.5281/zenodo.260230.svg"
42 alt="DOI"></a>
43 \endhtmlonly
44 
45 ## Table of Contents
46 - [Cheat sheet](#cheatsheet)
47 - [Introduction](#introduction)
48 - [Crash course](#crash-course)
49 - [More features](#more-features)
50 - [Transformations](#transformations) -- manipulating data
51 - [Actions](#actions) -- getting results
52 - [Parallel execution](#parallel-execution) -- how to use it and common pitfalls
53 - [Class reference](#reference) -- most methods are implemented in the RInterface base class
54 
55 ## <a name="cheatsheet"></a>Cheat sheet
56 These are the operations which can be performed with RDataFrame
57 
58 ### Transformations
59 Transformations are a way to manipulated the data.
60 
61 | **Transformation** | **Description** |
62 |------------------|-----------------|
63 | [Define](classROOT_1_1RDF_1_1RInterface.html#a7d48eb23b4378e99ebccb35e94ad025a) | Creates a new column in the dataset. |
64 | [DefineSlot](classROOT_1_1RDF_1_1RInterface.html#acaacf727b8a41d27c6bb4513348ac892) | Same as `Define`, but the user-defined function must take an extra `unsigned int slot` as its first parameter. `slot` will take a different value, `0` to `nThreads - 1`, for each thread of execution. This is meant as a helper in writing thread-safe `Define` transformation when using `RDataFrame` after `ROOT::EnableImplicitMT()`. `DefineSlot` works just as well with single-thread execution: in that case `slot` will always be `0`. |
65 | [DefineSlotEntry](classROOT_1_1RDF_1_1RInterface.html#a4f17074d5771916e3df18f8458186de7) | Same as `DefineSlot`, but the entry number is passed in addition to the slot number. This is meant as a helper in case some dependency on the entry number needs to be honoured. |
66 | [Filter](classROOT_1_1RDF_1_1RInterface.html#a70284a3bedc72b19610aaa91b5007ebd) | Filter the rows of the dataset. |
67 | [Range](classROOT_1_1RDF_1_1RInterface.html#a1b36b7868831de2375e061bb06cfc225) | Creates a node that filters entries based on range of entries |
68 
69 ### Actions
70 Actions are a way to produce a result out of the data. Each one is described in more detail in the reference guide.
71 
72 In the following, whenever we say an action "returns" something, we always mean it returns a smart pointer to it. Also
73 note that all actions are only executed for events that pass all preceding filters.
74 
75 Lazy actions only trigger the event loop when one of the results is accessed for the first time, making it easy to
76 produce several different results in one event loop. Instant actions trigger the event loop instantly.
77 
78 
79 | **Lazy action** | **Description** |
80 |------------------|-----------------|
81 | [Aggregate](classROOT_1_1RDF_1_1RInterface.html#ae540b00addc441f9b504cbae0ef0a24d) | Execute a user-defined accumulation operation on the processed column values. |
82 | [Book](classROOT_1_1RDF_1_1RInterface.html#a9b2f61f3333d1669e57055b9ae8be9d9) | Book execution of a custom action using a user-defined helper object. |
83 | [Cache](classROOT_1_1RDF_1_1RInterface.html#aaaa0a7bb8eb21315d8daa08c3e25f6c9) | Caches in contiguous memory columns' entries. Custom columns can be cached as well, filtered entries are not cached. Users can specify which columns to save (default is all). |
84 | [Count](classROOT_1_1RDF_1_1RInterface.html#a37f9e00c2ece7f53fae50b740adc1456) | Return the number of events processed. |
85 | [Fill](classROOT_1_1RDF_1_1RInterface.html#a0cac4d08297c23d16de81ff25545440a) | Fill a user-defined object with the values of the specified branches, as if by calling `Obj.Fill(branch1, branch2, ...). |
86 | [Histo{1D,2D,3D}](classROOT_1_1RDF_1_1RInterface.html#a247ca3aeb7ce5b95015b7fae72983055) | Fill a {one,two,three}-dimensional histogram with the processed branch values. |
87 | [Max](classROOT_1_1RDF_1_1RInterface.html#a057179b1e77599466a0b02200d5cd8c3) | Return the maximum of processed branch values. If the type of the column is inferred, the return type is `double`, the type of the column otherwise.|
88 | [Mean](classROOT_1_1RDF_1_1RInterface.html#ade6b020284f2f4fe9d3b09246b5f376a) | Return the mean of processed branch values.|
89 | [Min](classROOT_1_1RDF_1_1RInterface.html#a7005702189e601972b6d19ecebcdc80c) | Return the minimum of processed branch values. If the type of the column is inferred, the return type is `double`, the type of the column otherwise.|
90 | [Profile{1D,2D}](classROOT_1_1RDF_1_1RInterface.html#a8ef7dc16b0e9f7bc9cfbe2d9e5de0cef) | Fill a {one,two}-dimensional profile with the branch values that passed all filters. |
91 | [Reduce](classROOT_1_1RDF_1_1RInterface.html#a118e723ae29834df8f2a992ded347354) | Reduce (e.g. sum, merge) entries using the function (lambda, functor...) passed as argument. The function must have signature `T(T,T)` where `T` is the type of the branch. Return the final result of the reduction operation. An optional parameter allows initialization of the result object to non-default values. |
92 | [Report](classROOT_1_1RDF_1_1RInterface.html#a94f322531dcb25beb8f53a602e5d6332) | Obtains statistics on how many entries have been accepted and rejected by the filters. See the section on [named filters](#named-filters-and-cutflow-reports) for a more detailed explanation. The method returns a RCutFlowReport instance which can be queried programmatically to get information about the effects of the individual cuts. |
93 | [Sum](classROOT_1_1RDF_1_1RInterface.html#a61d03407459120df6749af43ed506891) | Return the sum of the values in the column. If the type of the column is inferred, the return type is `double`, the type of the column otherwise. |
94 | [Take](classROOT_1_1RDF_1_1RInterface.html#a4fd694773a2931b6b07737ddcd1e73b4) | Extract a column from the dataset as a collection of values. If the type of the column is a C-style array, the type stored in the return container is a `ROOT::VecOps::RVec<T>` to guarantee the lifetime of the data involved. |
95 
96 | **Instant action** | **Description** |
97 |---------------------|-----------------|
98 | [Foreach](classROOT_1_1RDF_1_1RInterface.html#ad2822a7ccb8a9afdf3e5b2ea321886ca) | Execute a user-defined function on each entry. Users are responsible for the thread-safety of this lambda when executing with implicit multi-threading enabled. |
99 | [ForeachSlot](classROOT_1_1RDF_1_1RInterface.html#a3650ca30aae1ccd0d92bf3d680314129) | Same as `Foreach`, but the user-defined function must take an extra `unsigned int slot` as its first parameter. `slot` will take a different value, `0` to `nThreads - 1`, for each thread of execution. This is meant as a helper in writing thread-safe `Foreach` actions when using `RDataFrame` after `ROOT::EnableImplicitMT()`. `ForeachSlot` works just as well with single-thread execution: in that case `slot` will always be `0`. |
100 | [Snapshot](classROOT_1_1RDF_1_1RInterface.html#a233b7723e498967f4340705d2c4db7f8) | Writes processed data-set to disk, in a new `TTree` and `TFile`. Custom columns can be saved as well, filtered entries are not saved. Users can specify which columns to save (default is all). Snapshot, by default, overwrites the output file if it already exists. `Snapshot` can be made *lazy* setting the appropriate flage in the snapshot options.|
101 
102 
103 ### Other Operations
104 
105 | **Operation** | **Description** |
106 |---------------------|-----------------|
107 | [Alias](classROOT_1_1RDF_1_1RInterface.html#a31ca327e4a192dcc05a4aac240e1a725) | Introduce an alias for a particular column name |
108 | [GetColumnNames](classROOT_1_1RDF_1_1RInterface.html#a951fe60b74d3a9fda37df59fd1dac186) | Get all the available columns of the dataset |
109 
110 
111 ## <a name="introduction"></a>Introduction
112 Users define their analysis as a sequence of operations to be performed on the data-frame object; the framework
113 takes care of the management of the loop over entries as well as low-level details such as I/O and parallelisation.
114 `RDataFrame` provides methods to perform most common operations required by ROOT analyses;
115 at the same time, users can just as easily specify custom code that will be executed in the event loop.
116 
117 `RDataFrame` is built with a *modular* and *flexible* workflow in mind, summarised as follows:
118 
119 1. **build a data-frame** object by specifying your data-set
120 2. **apply a series of transformations** to your data
121  1. **filter** (e.g. apply some cuts) or
122  2. **define** a new column (e.g. the result of an expensive computation on branches)
123 3. **apply actions** to the transformed data to produce results (e.g. fill a histogram)
124 
125 The following table shows how analyses based on `TTreeReader` and `TTree::Draw` translate to `RDataFrame`. Follow the
126 [crash course](#crash-course) to discover more idiomatic and flexible ways to express analyses with `RDataFrame`.
127 <table>
128 <tr>
129  <td>
130  <b>TTreeReader</b>
131  </td>
132  <td>
133  <b>ROOT::RDataFrame</b>
134  </td>
135 </tr>
136 <tr>
137  <td>
138 ~~~{.cpp}
139 TTreeReader reader("myTree", file);
140 TTreeReaderValue<A_t> a(reader, "A");
141 TTreeReaderValue<B_t> b(reader, "B");
142 TTreeReaderValue<C_t> c(reader, "C");
143 while(reader.Next()) {
144  if(IsGoodEvent(a, b, c))
145  DoStuff(a, b, c);
146 }
147 ~~~
148  </td>
149  <td>
150 ~~~{.cpp}
151 ROOT::RDataFrame d("myTree", file, {"A", "B", "C"});
152 d.Filter(IsGoodEvent).Foreach(DoStuff);
153 ~~~
154  </td>
155 </tr>
156 <tr>
157  <td>
158  <b>TTree::Draw</b>
159  </td>
160  <td>
161  <b>ROOT::RDataFrame</b>
162  </td>
163 </tr>
164 <tr>
165  <td>
166 ~~~{.cpp}
167 TTree *t = static_cast<TTree*>(
168  file->Get("myTree")
169 );
170 t->Draw("var", "var > 2");
171 ~~~
172  </td>
173  <td>
174 ~~~{.cpp}
175 ROOT::RDataFrame d("myTree", file);
176 auto h = d.Filter("var > 2").Histo1D("var");
177 ~~~
178  </td>
179 </tr>
180 </table>
181 
182 ## <a name="crash-course"></a> Crash course
183 All snippets of code presented in the crash course can be executed in the ROOT interpreter. Simply precede them with
184 ~~~{.cpp}
185 using namespace ROOT; // RDataFrame's namespace
186 ~~~
187 which is omitted for brevity. The terms "column" and "branch" are used interchangeably.
188 
189 ### Creating a RDataFrame
190 RDataFrame's constructor is where the user specifies the dataset and, optionally, a default set of columns that
191 operations should work with. Here are the most common methods to construct a RDataFrame object:
192 ~~~{.cpp}
193 // single file -- all ctors are equivalent
194 RDataFrame d1("treeName", "file.root");
195 TFile *f = TFile::Open("file.root");
196 RDataFrame d2("treeName", f); // same as TTreeReader
197 TTree *t = nullptr;
198 f.GetObject("treeName", t);
199 RDataFrame d3(*t); // TTreeReader takes a pointer, RDF takes a reference
200 
201 // multiple files -- all ctors are equivalent
202 RDataFrame d3("myTree", {"file1.root", "file2.root"});
203 std::vector<std::string> files = {"file1.root", "file2.root"};
204 RDataFrame d3("myTree", files);
205 RDataFrame d4("myTree", "file*.root"); // see TRegexp's documentation for a list of valid regexes
206 TChain chain("myTree");
207 chain.Add("file1.root");
208 chain.Add("file2.root");
209 RDataFrame d3(chain);
210 ~~~
211 Additionally, users can construct a RDataFrame specifying just an integer number. This is the number of "events" that
212 will be generated by this RDataFrame.
213 ~~~{.cpp}
214 RDataFrame d(10); // a RDF with 10 entries (and no columns/branches, for now)
215 d.Foreach([] { static int i = 0; std::cout << i++ << std::endl; }); // silly example usage: count to ten
216 ~~~
217 This is useful to generate simple data-sets on the fly: the contents of each event can be specified via the `Define`
218 transformation (explained below). For example, we have used this method to generate Pythia events (with a `Define`
219 transformation) and write them to disk in parallel (with the `Snapshot` action).
220 
221 ### Programmatically get the list of column names
222 The list of column names available in the dataset can be obtained with the `GetColumnsNames` method:
223 ~~~{.cpp}
224 RDataFrame d("myTree", "file.root");
225 auto colNames = d.GetColumnNames();
226 for (auto &&colName : colNames) {
227  std::cout << colName << std::endl;
228  }
229 ~~~
230 
231 ### Filling a histogram
232 Let's now tackle a very common task, filling a histogram:
233 ~~~{.cpp}
234 // Fill a TH1D with the "MET" branch
235 RDataFrame d("myTree", "file.root");
236 auto h = d.Histo1D("MET");
237 h->Draw();
238 ~~~
239 The first line creates a `RDataFrame` associated to the `TTree` "myTree". This tree has a branch named "MET".
240 
241 `Histo1D` is an *action*; it returns a smart pointer (a `RResultPtr` to be precise) to a `TH1D` histogram filled
242 with the `MET` of all events. If the quantity stored in the branch is a collection (e.g. a vector or an array), the
243 histogram is filled with its elements.
244 
245 You can use the objects returned by actions as if they were pointers to the desired results. There are many other
246 possible [actions](#overview), and all their results are wrapped in smart pointers; we'll see why in a minute.
247 
248 ### Applying a filter
249 Let's say we want to cut over the value of branch "MET" and count how many events pass this cut. This is one way to do it:
250 ~~~{.cpp}
251 RDataFrame d("myTree", "file.root");
252 auto c = d.Filter("MET > 4.").Count();
253 std::cout << *c << std::endl;
254 ~~~
255 The filter string (which must contain a valid c++ expression) is applied to the specified branches for each event;
256 the name and types of the columns are inferred automatically. The string expression is required to return a `bool`
257 which signals whether the event passes the filter (`true`) or not (`false`).
258 
259 You can think of your data as "flowing" through the chain of calls, being transformed, filtered and finally used to
260 perform actions. Multiple `Filter` calls can be chained one after another.
261 
262 Using string filters is nice for simple things, but they are limited to specifying the equivalent of a single return
263 statement or the body of a lambda, so it's cumbersome to use strings with more complex filters. They also add a small
264 runtime overhead, as ROOT needs to just-in-time compile the string into C++ code. When more freedom is required or
265 runtime performance is very important, a C++ callable can be specified instead (a lambda in the following snippet,
266 but it can be any kind of function or even a functor class), together with a list of branch names.
267 This snippet is analogous to the one above:
268 ~~~{.cpp}
269 RDataFrame d("myTree", "file.root");
270 auto metCut = [](double x) { return x > 4.; }; // a c++11 lambda function checking "x > 4"
271 auto c = d.Filter(metCut, {"MET"}).Count();
272 std::cout << *c << std::endl;
273 ~~~
274 
275 An example of a more complex filter expressed as a string containing C++ code is shown below
276 
277 ~~~{.cpp}
278 RDataFrame d("myTree", "file.root");
279 auto df = d.Define("p", "std::array<double, 4> p{px, py, pz}; return p;")
280  .Filter("double p2 = 0.0; for (auto&& x : p) p2 += x*x; return sqrt(p2) < 10.0;");
281 ~~~
282 
283 The code snippet above defines a column `p` that is a fixed-size array using the component column names and then
284 filters on its magnitude by looping over its elements. It must be noted that the usage of strings to define columns
285 like the one above is a major advantage when using PyROOT. However, only constants and data coming from other columns
286 in the dataset can be involved in the code passed as a string. Local variables and functions cannot be used, since
287 the interpreter will not know how to find them. When capturing local state is necessary, a C++ callable can be used.
288 
289 More information on filters and how to use them to automatically generate cutflow reports can be found [below](#Filters).
290 
291 ### Defining custom columns
292 Let's now consider the case in which "myTree" contains two quantities "x" and "y", but our analysis relies on a derived
293 quantity `z = sqrt(x*x + y*y)`. Using the `Define` transformation, we can create a new column in the data-set containing
294 the variable "z":
295 ~~~{.cpp}
296 RDataFrame d("myTree", "file.root");
297 auto sqrtSum = [](double x, double y) { return sqrt(x*x + y*y); };
298 auto zMean = d.Define("z", sqrtSum, {"x","y"}).Mean("z");
299 std::cout << *zMean << std::endl;
300 ~~~
301 `Define` creates the variable "z" by applying `sqrtSum` to "x" and "y". Later in the chain of calls we refer to
302 variables created with `Define` as if they were actual tree branches/columns, but they are evaluated on demand, at most
303 once per event. As with filters, `Define` calls can be chained with other transformations to create multiple custom
304 columns. `Define` and `Filter` transformations can be concatenated and intermixed at will.
305 
306 As with filters, it is possible to specify new columns as string expressions. This snippet is analogous to the one above:
307 ~~~{.cpp}
308 RDataFrame d("myTree", "file.root");
309 auto zMean = d.Define("z", "sqrt(x*x + y*y)").Mean("z");
310 std::cout << *zMean << std::endl;
311 ~~~
312 Again the names of the branches used in the expression and their types are inferred automatically. The string must be
313 valid c++ and is just-in-time compiled by the ROOT interpreter, cling -- the process has a small runtime overhead.
314 
315 Previously, when showing the different ways a RDataFrame can be created, we showed a constructor that only takes a
316 number of entries a parameter. In the following example we show how to combine such an "empty" `RDataFrame` with `Define`
317 transformations to create a data-set on the fly. We then save the generated data on disk using the `Snapshot` action.
318 ~~~{.cpp}
319 RDataFrame d(100); // a RDF that will generate 100 entries (currently empty)
320 int x = -1;
321 auto d_with_columns = d.Define("x", [&x] { return ++x; })
322  .Define("xx", [&x] { return x*x; });
323 d_with_columns.Snapshot("myNewTree", "newfile.root");
324 ~~~
325 This example is slightly more advanced than what we have seen so far: for starters, it makes use of lambda captures (a
326 simple way to make external variables available inside the body of c++ lambdas) to act on the same variable `x` from
327 both `Define` transformations. Secondly we have *stored* the transformed data-frame in a variable. This is always
328 possible: at each point of the transformation chain, users can store the status of the data-frame for further use (more
329 on this [below](#callgraphs)).
330 
331 You can read more about defining new columns [here](#custom-columns).
332 
333 \image html RDF_Graph.png "A graph composed of two branches, one starting with a filter and one with a define. The end point of a branch is always an action."
334 
335 ### Running on a range of entries
336 It is sometimes necessary to limit the processing of the dataset to a range of entries. For this reason, the RDataFrame
337 offers the concept of ranges as a node of the RDataFrame chain of transformations; this means that filters, columns and
338 actions can be concatenated to and intermixed with `Range`s. If a range is specified after a filter, the range will act
339 exclusively on the entries passing the filter -- it will not even count the other entries! The same goes for a `Range`
340 hanging from another `Range`. Here are some commented examples:
341 ~~~{.cpp}
342 RDataFrame d("myTree", "file.root");
343 // Here we store a data-frame that loops over only the first 30 entries in a variable
344 auto d30 = d.Range(30);
345 // This is how you pick all entries from 15 onwards
346 auto d15on = d.Range(15, 0);
347 // We can specify a stride too, in this case we pick an event every 3
348 auto d15each3 = d.Range(0, 15, 3);
349 ~~~
350 Note that ranges are not available when multi-threading is enabled. More information on ranges is available
351 [here](#ranges).
352 
353 ### Executing multiple actions in the same event loop
354 As a final example let us apply two different cuts on branch "MET" and fill two different histograms with the "pt\_v" of
355 the filtered events.
356 By now, you should be able to easily understand what's happening:
357 ~~~{.cpp}
358 RDataFrame d("treeName", "file.root");
359 auto h1 = d.Filter("MET > 10").Histo1D("pt_v");
360 auto h2 = d.Histo1D("pt_v");
361 h1->Draw(); // event loop is run once here
362 h2->Draw("SAME"); // no need to run the event loop again
363 ~~~
364 `RDataFrame` executes all above actions by **running the event-loop only once**. The trick is that actions are not
365 executed at the moment they are called, but they are **lazy**, i.e. delayed until the moment one of their results is
366 accessed through the smart pointer. At that time, the event loop is triggered and *all* results are produced
367 simultaneously.
368 
369 It is therefore good practice to declare all your transformations and actions *before* accessing their results, allowing
370 `RDataFrame` to run the loop once and produce all results in one go.
371 
372 ### Going parallel
373 Let's say we would like to run the previous examples in parallel on several cores, dividing events fairly between cores.
374 The only modification required to the snippets would be the addition of this line *before* constructing the main
375 data-frame object:
376 ~~~{.cpp}
377 ROOT::EnableImplicitMT();
378 ~~~
379 Simple as that. More details are given [below](#parallel-execution).
380 
381 ## <a name="more-features"></a>More features
382 Here is a list of the most important features that have been omitted in the "Crash course" for brevity.
383 You don't need to read all these to start using `RDataFrame`, but they are useful to save typing time and runtime.
384 
385 ### Treatment of columns holding collections
386 When using RDataFrame to read data from a ROOT file, users can specify that the type of a branch is `RVec<T>`
387 to indicate the branch is a c-style array, an STL array or any other collection type associated to a
388 contiguous storage in memory.
389 
390 Column values of type `RVec<T>` perform no copy of the underlying array data, it's in some sense a view,
391 and offer a rich interface to operate on the array elements in a vectorised fashion.
392 
393 The `RVec<T>` type signals to RDataFrame that a special behaviour needs to be adopted when snapshotting
394 a dataset on disk. Indeed, if columns which are variable size C arrays are treated via the `RVec<T>`,
395 RDataFrame will correctly persistify them - if anything else is adopted, for example `std::span`, only
396 the first element of the array will be written.
397 
398 ### Callbacks
399 Acting on a RResultPtr, it is possible to register a callback that RDataFrame will execute "everyNEvents"
400 on a partial result.
401 
402 The callback must be a callable that takes a reference to the result type as argument and returns nothing.
403 RDataFrame, acting as a full fledged data processing framework, will invoke registered callbacks passing
404 partial action results as arguments to them (e.g. a histogram filled with a part of the selected events).
405 
406 Callbacks can be used e.g. to inspect partial results of the analysis while the event loop is running. For
407 example one can draw an up-to-date version of a result histogram every 100 entries like this:
408 ~~~{.cpp}
409 auto h = tdf.Histo1D("x");
410 TCanvas c("c","x hist");
411 h.OnPartialResult(100, [&c](TH1D &h_) { c.cd(); h_.Draw(); c.Update(); });
412 h->Draw(); // event loop runs here, this `Draw` is executed after the event loop is finished
413 ~~~
414 
415 ### Default branch lists
416 When constructing a `RDataFrame` object, it is possible to specify a **default column list** for your analysis, in the
417 usual form of a list of strings representing branch/column names. The default column list will be used as a fallback
418 whenever a list specific to the transformation/action is not present. RDataFrame will take as many of these columns as
419 needed, ignoring trailing extra names if present.
420 ~~~{.cpp}
421 // use "b1" and "b2" as default branches
422 RDataFrame d1("myTree", "file.root", {"b1","b2"});
423 auto h = d1.Filter([](int b1, int b2) { return b1 > b2; }) // will act on "b1" and "b2"
424  .Histo1D(); // will act on "b1"
425 
426 // just one default branch this time
427 RDataFrame d2("myTree", "file.root", {"b1"});
428 auto min = d2.Filter([](double b2) { return b2 > 0; }, {"b2"}) // we can still specify non-default branch lists
429  .Min(); // returns the minimum value of "b1" for the filtered entries
430 ~~~
431 
432 ### Branch type guessing and explicit declaration of branch types
433 C++ is a statically typed language: all types must be known at compile-time. This includes the types of the `TTree`
434 branches we want to work on. For filters, temporary columns and some of the actions, **branch types are deduced from the
435 signature** of the relevant filter function/temporary column expression/action function:
436 ~~~{.cpp}
437 // here b1 is deduced to be `int` and b2 to be `double`
438 dataFrame.Filter([](int x, double y) { return x > 0 && y < 0.; }, {"b1", "b2"});
439 ~~~
440 If we specify an incorrect type for one of the branches, an exception with an informative message will be thrown at
441 runtime, when the branch value is actually read from the `TTree`: `RDataFrame` detects type mismatches. The same would
442 happen if we swapped the order of "b1" and "b2" in the branch list passed to `Filter`.
443 
444 Certain actions, on the other hand, do not take a function as argument (e.g. `Histo1D`), so we cannot deduce the type of
445 the branch at compile-time. In this case **`RDataFrame` infers the type of the branch** from the `TTree` itself. This
446 is why we never needed to specify the branch types for all actions in the above snippets.
447 
448 When the branch type is not a common one such as `int`, `double`, `char` or `float` it is nonetheless good practice to
449 specify it as a template parameter to the action itself, like this:
450 ~~~{.cpp}
451 dataFrame.Histo1D("b1"); // OK, the type of "b1" is deduced at runtime
452 dataFrame.Min<MyNumber_t>("myObject"); // OK, "myObject" is deduced to be of type `MyNumber_t`
453 ~~~
454 
455 Deducing types at runtime requires the just-in-time compilation of the relevant actions, which has a small runtime
456 overhead, so specifying the type of the columns as template parameters to the action is good practice when performance is a goal.
457 
458 ### Generic actions
459 `RDataFrame` strives to offer a comprehensive set of standard actions that can be performed on each event. At the same
460 time, it **allows users to execute arbitrary code (i.e. a generic action) inside the event loop** through the `Foreach`
461 and `ForeachSlot` actions.
462 
463 `Foreach(f, columnList)` takes a function `f` (lambda expression, free function, functor...) and a list of columns, and
464 executes `f` on those columns for each event. The function passed must return nothing (i.e. `void`). It can be used to
465 perform actions that are not already available in the interface. For example, the following snippet evaluates the root
466 mean square of column "b":
467 ~~~{.cpp}
468 // Single-thread evaluation of RMS of column "b" using Foreach
469 double sumSq = 0.;
470 unsigned int n = 0;
471 RDataFrame d("bTree", bFilePtr);
472 d.Foreach([&sumSq, &n](double b) { ++n; sumSq += b*b; }, {"b"});
473 std::cout << "rms of b: " << std::sqrt(sumSq / n) << std::endl;
474 ~~~
475 When executing on multiple threads, users are responsible for the thread-safety of the expression passed to `Foreach`.
476 The code above would need to employ some resource protection mechanism to ensure non-concurrent writing of `rms`; but
477 this is probably too much head-scratch for such a simple operation.
478 
479 `ForeachSlot` can help in this situation. It is an alternative version of `Foreach` for which the function takes an
480 additional parameter besides the columns it should be applied to: an `unsigned int slot` parameter, where `slot` is a
481 number indicating which thread (0, 1, 2 , ..., poolSize - 1) the function is being run in. We can take advantage of
482 `ForeachSlot` to evaluate a thread-safe root mean square of branch "b":
483 ~~~{.cpp}
484 // Thread-safe evaluation of RMS of branch "b" using ForeachSlot
485 ROOT::EnableImplicitMT();
486 const unsigned int nSlots = ROOT::GetImplicitMTPoolSize();
487 std::vector<double> sumSqs(nSlots, 0.);
488 std::vector<unsigned int> ns(nSlots, 0);
489 
490 RDataFrame d("bTree", bFilePtr);
491 d.ForeachSlot([&sumSqs, &ns](unsigned int slot, double b) { sumSqs[slot] += b*b; ns[slot] += 1; }, {"b"});
492 double sumSq = std::accumulate(sumSqs.begin(), sumSqs.end(), 0.); // sum all squares
493 unsigned int n = std::accumulate(ns.begin(), ns.end(), 0); // sum all counts
494 std::cout << "rms of b: " << std::sqrt(sumSq / n) << std::endl;
495 ~~~
496 You see how we created one `double` variable for each thread in the pool, and later merged their results via
497 `std::accumulate`.
498 
499 ### Friend trees
500 Friend trees are supported by RDataFrame.
501 In order to deal with friend trees with RDataFrame, the user is required to build
502 the tree and its friends and instantiate a RDataFrame with it.
503 Two caveats are presents when using jitted `Define`s and `Filter`s:
504 1) the only columns which can be used in the strings passed to the aforementioned transformations are the top level branches of the friend trees.
505 2) the "friend columns" cannot be written with the notation involving a dot. For example, if a tree is created like this:
506 ~~~{.cpp}
507 TTree t([...]);
508 TTree ft([...]);
509 t.AddFriend(t,"myFriend");
510 ~~~
511 in order to access a certain column `col` of the tree ft, it will be necessary to alias it before. To continue the example:
512 ~~~{.cpp}
513 RDataFrame d(t);
514 d.Alias("myFriend_MyCol", "myFriend.MyCol");
515 auto f = d.Filter("myFriend_MyCol == 42");
516 ~~~
517 
518 ### Reading file formats different from ROOT's
519 RDataFrame can be interfaced with RDataSources. The RDataSource interface defines an API that RDataFrame can use to read arbitrary data formats.
520 
521 A concrete RDataSource implementation (i.e. a class that inherits from RDataSource and implements all of its pure
522 methods) provides an adaptor that RDataFrame can leverage to read any kind of tabular data formats.
523 RDataFrame calls into RDataSource to retrieve information about the data, retrieve (thread-local) readers or "cursors" for selected columns
524 and to advance the readers to the desired data entry.
525 Some predefined RDataSources are natively provided by ROOT such as the `RCsvDS` which allows to read comma separated files:
526 ~~~{.cpp}
527 auto tdf = ROOT::RDF::MakeCsvDataFrame("MuRun2010B.csv");
528 auto filteredEvents =
529  tdf.Filter("Q1 * Q2 == -1")
530  .Define("m", "sqrt(pow(E1 + E2, 2) - (pow(px1 + px2, 2) + pow(py1 + py2, 2) + pow(pz1 + pz2, 2)))");
531 auto h = filteredEvents.Histo1D("m");
532 h->Draw();
533 ~~~
534 
535 ### <a name="callgraphs"></a>Call graphs (storing and reusing sets of transformations)
536 **Sets of transformations can be stored as variables** and reused multiple times to create **call graphs** in which
537 several paths of filtering/creation of columns are executed simultaneously; we often refer to this as "storing the
538 state of the chain".
539 
540 This feature can be used, for example, to create a temporary column once and use it in several subsequent filters or
541 actions, or to apply a strict filter to the data-set *before* executing several other transformations and actions,
542 effectively reducing the amount of events processed.
543 
544 Let's try to make this clearer with a commented example:
545 ~~~{.cpp}
546 // build the data-frame and specify a default column list
547 RDataFrame d(treeName, filePtr, {"var1", "var2", "var3"});
548 
549 // apply a cut and save the state of the chain
550 auto filtered = d.Filter(myBigCut);
551 
552 // plot branch "var1" at this point of the chain
553 auto h1 = filtered.Histo1D("var1");
554 
555 // create a new branch "vec" with a vector extracted from a complex object (only for filtered entries)
556 // and save the state of the chain
557 auto newBranchFiltered = filtered.Define("vec", [](const Obj& o) { return o.getVector(); }, {"obj"});
558 
559 // apply a cut and fill a histogram with "vec"
560 auto h2 = newBranchFiltered.Filter(cut1).Histo1D("vec");
561 
562 // apply a different cut and fill a new histogram
563 auto h3 = newBranchFiltered.Filter(cut2).Histo1D("vec");
564 
565 // Inspect results
566 h2->Draw(); // first access to an action result: run event-loop!
567 h3->Draw("SAME"); // event loop does not need to be run again here..
568 std::cout << "Entries in h1: " << h1->GetEntries() << std::endl; // ..or here
569 ~~~
570 `RDataFrame` detects when several actions use the same filter or the same temporary column, and **only evaluates each
571 filter or temporary column once per event**, regardless of how many times that result is used down the call graph.
572 Objects read from each column are **built once and never copied**, for maximum efficiency.
573 When "upstream" filters are not passed, subsequent filters, temporary column expressions and actions are not evaluated,
574 so it might be advisable to put the strictest filters first in the chain.
575 
576 ## <a name="transformations"></a>Transformations
577 ### <a name="Filters"></a> Filters
578 A filter is defined through a call to `Filter(f, columnList)`. `f` can be a function, a lambda expression, a functor
579 class, or any other callable object. It must return a `bool` signalling whether the event has passed the selection
580 (`true`) or not (`false`). It must perform "read-only" actions on the columns, and should not have side-effects (e.g.
581 modification of an external or static variable) to ensure correct results when implicit multi-threading is active.
582 
583 `RDataFrame` only evaluates filters when necessary: if multiple filters are chained one after another, they are executed
584 in order and the first one returning `false` causes the event to be discarded and triggers the processing of the next
585 entry. If multiple actions or transformations depend on the same filter, that filter is not executed multiple times for
586 each entry: after the first access it simply serves a cached result.
587 
588 #### <a name="named-filters-and-cutflow-reports"></a>Named filters and cutflow reports
589 An optional string parameter `name` can be passed to the `Filter` method to create a **named filter**. Named filters
590 work as usual, but also keep track of how many entries they accept and reject.
591 
592 Statistics are retrieved through a call to the `Report` method:
593 
594 - when `Report` is called on the main `RDataFrame` object, it returns a RResultPtr<RCutFlowReport> relative to all
595 named filters declared up to that point
596 - when called on a specific node (e.g. the result of a `Define` or `Filter`), it returns a RResultPtr<RCutFlowReport>
597 relative all named filters in the section of the chain between the main `RDataFrame` and that node (included).
598 
599 Stats are stored in the same order as named filters have been added to the graph, and *refer to the latest event-loop*
600 that has been run using the relevant `RDataFrame`.
601 
602 ### <a name="ranges"></a>Ranges
603 When `RDataFrame` is not being used in a multi-thread environment (i.e. no call to `EnableImplicitMT` was made),
604 `Range` transformations are available. These act very much like filters but instead of basing their decision on
605 a filter expression, they rely on `begin`,`end` and `stride` parameters.
606 
607 - `begin`: initial entry number considered for this range.
608 - `end`: final entry number (excluded) considered for this range. 0 means that the range goes until the end of the dataset.
609 - `stride`: process one entry of the [begin, end) range every `stride` entries. Must be strictly greater than 0.
610 
611 The actual number of entries processed downstream of a `Range` node will be `(end - begin)/stride` (or less if less
612 entries than that are available).
613 
614 Note that ranges act "locally", not based on the global entry count: `Range(10,50)` means "skip the first 10 entries
615 *that reach this node*, let the next 40 entries pass, then stop processing". If a range node hangs from a filter node,
616 and the range has a `begin` parameter of 10, that means the range will skip the first 10 entries *that pass the
617 preceding filter*.
618 
619 Ranges allow "early quitting": if all branches of execution of a functional graph reached their `end` value of
620 processed entries, the event-loop is immediately interrupted. This is useful for debugging and quick data explorations.
621 
622 ### <a name="custom-columns"></a> Custom columns
623 Custom columns are created by invoking `Define(name, f, columnList)`. As usual, `f` can be any callable object
624 (function, lambda expression, functor class...); it takes the values of the columns listed in `columnList` (a list of
625 strings) as parameters, in the same order as they are listed in `columnList`. `f` must return the value that will be
626 assigned to the temporary column.
627 
628 A new variable is created called `name`, accessible as if it was contained in the dataset from subsequent
629 transformations/actions.
630 
631 Use cases include:
632 - caching the results of complex calculations for easy and efficient multiple access
633 - extraction of quantities of interest from complex objects
634 - branch aliasing, i.e. changing the name of a branch
635 
636 An exception is thrown if the `name` of the new column/branch is already in use for another branch in the `TTree`.
637 
638 It is also possible to specify the quantity to be stored in the new temporary column as a C++ expression with the method
639 `Define(name, expression)`. For example this invocation
640 
641 ~~~{.cpp}
642 tdf.Define("pt", "sqrt(px*px + py*py)");
643 ~~~
644 
645 will create a new column called "pt" the value of which is calculated starting from the columns px and py. The system
646 builds a just-in-time compiled function starting from the expression after having deduced the list of necessary branches
647 from the names of the variables specified by the user.
648 
649 #### Custom columns as function of slot and entry number
650 
651 It is possible to create custom columns also as a function of the processing slot and entry numbers. The methods that can
652 be invoked are:
653 - `DefineSlot(name, f, columnList)`. In this case the callable f has this signature `R(unsigned int, T1, T2, ...)`: the
654 first parameter is the slot number which ranges from 0 to ROOT::GetImplicitMTPoolSize() - 1.
655 - `DefineSlotEntry(name, f, columnList)`. In this case the callable f has this signature `R(unsigned int, ULong64_t,
656 T1, T2, ...)`: the first parameter is the slot number while the second one the number of the entry being processed.
657 
658 ## <a name="actions"></a>Actions
659 ### Instant and lazy actions
660 Actions can be **instant** or **lazy**. Instant actions are executed as soon as they are called, while lazy actions are
661 executed whenever the object they return is accessed for the first time. As a rule of thumb, actions with a return value
662 are lazy, the others are instant.
663 
664 ## <a name="parallel-execution"></a>Parallel execution
665 As pointed out before in this document, `RDataFrame` can transparently perform multi-threaded event loops to speed up
666 the execution of its actions. Users only have to call `ROOT::EnableImplicitMT()` *before* constructing the `RDataFrame`
667 object to indicate that it should take advantage of a pool of worker threads. **Each worker thread processes a distinct
668 subset of entries**, and their partial results are merged before returning the final values to the user.
669 
670 ### Thread safety
671 `Filter` and `Define` transformations should be inherently thread-safe: they have no side-effects and are not
672 dependent on global state.
673 Most `Filter`/`Define` functions will in fact be pure in the functional programming sense.
674 All actions are built to be thread-safe with the exception of `Foreach`, in which case users are responsible of
675 thread-safety, see [here](#generic-actions).
676 
677 <a name="reference"></a>
678 */
679 // clang-format on
680 namespace ROOT {
681 
682 ////////////////////////////////////////////////////////////////////////////
683 /// \brief Build the dataframe
684 /// \param[in] treeName Name of the tree contained in the directory
685 /// \param[in] dirPtr TDirectory where the tree is stored, e.g. a TFile.
686 /// \param[in] defaultBranches Collection of default branches.
687 ///
688 /// The default branches are looked at in case no branch is specified in the
689 /// booking of actions or transformations.
690 /// See RInterface for the documentation of the
691 /// methods available.
692 RDataFrame::RDataFrame(std::string_view treeName, TDirectory *dirPtr, const ColumnNames_t &defaultBranches)
693  : RInterface<RDFDetail::RLoopManager>(std::make_shared<RDFDetail::RLoopManager>(nullptr, defaultBranches))
694 {
695  if (!dirPtr) {
696  auto msg = "Invalid TDirectory!";
697  throw std::runtime_error(msg);
698  }
699  const std::string treeNameInt(treeName);
700  auto tree = static_cast<TTree *>(dirPtr->Get(treeNameInt.c_str()));
701  if (!tree) {
702  auto msg = "Tree \"" + treeNameInt + "\" cannot be found!";
703  throw std::runtime_error(msg);
704  }
705  GetProxiedPtr()->SetTree(std::shared_ptr<TTree>(tree, [](TTree *) {}));
706 }
707 
708 ////////////////////////////////////////////////////////////////////////////
709 /// \brief Build the dataframe
710 /// \param[in] treeName Name of the tree contained in the directory
711 /// \param[in] filenameglob TDirectory where the tree is stored, e.g. a TFile.
712 /// \param[in] defaultBranches Collection of default branches.
713 ///
714 /// The default branches are looked at in case no branch is specified in the
715 /// booking of actions or transformations.
716 /// See RInterface for the documentation of the
717 /// methods available.
718 RDataFrame::RDataFrame(std::string_view treeName, std::string_view filenameglob, const ColumnNames_t &defaultBranches)
719  : RInterface<RDFDetail::RLoopManager>(std::make_shared<RDFDetail::RLoopManager>(nullptr, defaultBranches))
720 {
721  const std::string treeNameInt(treeName);
722  const std::string filenameglobInt(filenameglob);
723  auto chain = std::make_shared<TChain>(treeNameInt.c_str());
724  chain->Add(filenameglobInt.c_str());
725  GetProxiedPtr()->SetTree(chain);
726 }
727 
728 ////////////////////////////////////////////////////////////////////////////
729 /// \brief Build the dataframe
730 /// \param[in] treeName Name of the tree contained in the directory
731 /// \param[in] filenames Collection of file names
732 /// \param[in] defaultBranches Collection of default branches.
733 ///
734 /// The default branches are looked at in case no branch is specified in the booking of actions or transformations.
735 /// See RInterface for the documentation of the methods available.
736 RDataFrame::RDataFrame(std::string_view treeName, const std::vector<std::string> &filenames,
737  const ColumnNames_t &defaultBranches)
738  : RDF::RInterface<RDFDetail::RLoopManager>(std::make_shared<RDFDetail::RLoopManager>(nullptr, defaultBranches))
739 {
740  std::string treeNameInt(treeName);
741  auto chain = std::make_shared<TChain>(treeNameInt.c_str());
742  for (auto &fileName : filenames)
743  chain->Add(fileName.c_str());
744  GetProxiedPtr()->SetTree(chain);
745 }
746 
747 ////////////////////////////////////////////////////////////////////////////
748 /// \brief Build the dataframe
749 /// \param[in] tree The tree or chain to be studied.
750 /// \param[in] defaultBranches Collection of default column names to fall back to when none is specified.
751 ///
752 /// The default branches are looked at in case no branch is specified in the
753 /// booking of actions or transformations.
754 /// See RInterface for the documentation of the
755 /// methods available.
757  : RInterface<RDFDetail::RLoopManager>(std::make_shared<RDFDetail::RLoopManager>(&tree, defaultBranches))
758 {
759 }
760 
761 //////////////////////////////////////////////////////////////////////////
762 /// \brief Build a dataframe that generates numEntries entries.
763 /// \param[in] numEntries The number of entries to generate.
764 ///
765 /// An empty-source dataframe constructed with a number of entries will
766 /// generate those entries on the fly when some action is triggered,
767 /// and it will do so for all the previously-defined temporary branches.
769  : RInterface<RDFDetail::RLoopManager>(std::make_shared<RDFDetail::RLoopManager>(numEntries))
770 {
771 }
772 
773 //////////////////////////////////////////////////////////////////////////
774 /// \brief Build dataframe associated to datasource.
775 /// \param[in] ds The data-source object.
776 /// \param[in] defaultBranches Collection of default column names to fall back to when none is specified.
777 ///
778 /// A dataframe associated to a datasource will query it to access column values.
779 RDataFrame::RDataFrame(std::unique_ptr<RDataSource> ds, const ColumnNames_t &defaultBranches)
780  : RInterface<RDFDetail::RLoopManager>(std::make_shared<RDFDetail::RLoopManager>(std::move(ds), defaultBranches))
781 {
782 }
783 
784 } // namespace ROOT
785 
786 namespace cling {
787 //////////////////////////////////////////////////////////////////////////
788 /// Print a RDataFrame at the prompt
789 std::string printValue(ROOT::RDataFrame *tdf)
790 {
791  auto df = tdf->GetLoopManager();
792  auto *tree = df->GetTree();
793  auto defBranches = df->GetDefaultColumnNames();
794 
795  std::ostringstream ret;
796  if (tree) {
797  ret << "A data frame built on top of the " << tree->GetName() << " dataset.";
798  if (!defBranches.empty()) {
799  if (defBranches.size() == 1)
800  ret << "\nDefault branch: " << defBranches[0];
801  else {
802  ret << "\nDefault branches:\n";
803  for (auto &&branch : defBranches) {
804  ret << " - " << branch << "\n";
805  }
806  }
807  }
808  } else {
809  ret << "A data frame that will create " << df->GetNEmptyEntries() << " entries\n";
810  }
811 
812  return ret.str();
813 }
814 } // namespace cling
815 
Namespace for new ROOT classes and functions.
Definition: StringConv.hxx:21
virtual TObject * Get(const char *namecycle)
Return pointer to object identified by namecycle.
Definition: TDirectory.cxx:805
std::shared_ptr< RLoopManager > GetLoopManager()
Get the RLoopManager if reachable. If not, throw.
STL namespace.
std::string printValue(const TDatime *val)
Print a TDatime at the prompt.
Definition: TDatime.cxx:514
ROOT&#39;s RDataFrame offers a high level interface for analyses of data stored in TTrees, CSV&#39;s and other data formats.
Definition: RDataFrame.hxx:42
const std::shared_ptr< RDFDetail::RLoopManager > & GetProxiedPtr() const
Describe directory structure in memory.
Definition: TDirectory.h:34
unsigned long long ULong64_t
Definition: RtypesCore.h:70
Print a TSeq at the prompt:
Definition: TDatime.h:115
basic_string_view< char > string_view
Definition: RStringView.hxx:35
RDataFrame(std::string_view treeName, std::string_view filenameglob, const ColumnNames_t &defaultBranches={})
Build the dataframe.
Definition: RDataFrame.cxx:718
Definition: tree.py:1
A TTree object has a header with a name and a title.
Definition: TTree.h:70
RDFDetail::ColumnNames_t ColumnNames_t
Definition: RDataFrame.hxx:45