Logo ROOT  
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-2018, 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 "ROOT/RDataSource.hxx"
16 #include "TChain.h"
17 #include "TDirectory.h"
18 
19 // clang-format off
20 /**
21 * \class ROOT::RDataFrame
22 * \ingroup dataframe
23 * \brief ROOT's RDataFrame offers a high level interface for analyses of data stored in `TTree`s, CSV's and other data formats.
24 
25 In addition, multi-threading and other low-level optimisations allow users to exploit all the resources available
26 on their machines completely transparently.<br>
27 Skip to the [class reference](#reference) or keep reading for the user guide.
28 
29 In a nutshell:
30 ~~~{.cpp}
31 ROOT::EnableImplicitMT(); // Tell ROOT you want to go parallel
32 ROOT::RDataFrame d("myTree", "file_*.root"); // Interface to TTree and TChain
33 auto myHisto = d.Histo1D("Branch_A"); // This happens in parallel!
34 myHisto->Draw();
35 ~~~
36 
37 Calculations are expressed in terms of a type-safe *functional chain of actions and transformations*, `RDataFrame` takes
38 care of their execution. The implementation automatically puts in place several low level optimisations such as
39 multi-thread parallelisation and caching.
40 
41 \htmlonly
42 <a href="https://doi.org/10.5281/zenodo.260230"><img src="https://zenodo.org/badge/DOI/10.5281/zenodo.260230.svg"
43 alt="DOI"></a>
44 \endhtmlonly
45 
46 ## For the impatient user
47 You can directly see RDataFrame in action through its [code examples](https://root.cern.ch/doc/master/group__tutorial__dataframe.html), both in C++ and Python.
48 
49 ## Table of Contents
50 - [Cheat sheet](#cheatsheet)
51 - [Introduction](#introduction)
52 - [Crash course](#crash-course)
53 - [Efficient analysis in Python](#python)
54 - [More features](#more-features)
55 - [Transformations](#transformations) -- manipulating data
56 - [Actions](#actions) -- getting results
57 - [Parallel execution](#parallel-execution) -- how to use it and common pitfalls
58 - [Class reference](#reference) -- most methods are implemented in the [RInterface](https://root.cern/doc/master/classROOT_1_1RDF_1_1RInterface.html) base class
59 
60 ## <a name="cheatsheet"></a>Cheat sheet
61 These are the operations which can be performed with RDataFrame
62 
63 ### Transformations
64 Transformations are a way to manipulate the data.
65 
66 | **Transformation** | **Description** |
67 |------------------|--------------------|
68 | [Define](classROOT_1_1RDF_1_1RInterface.html#a7d48eb23b4378e99ebccb35e94ad025a) | Creates a new column in the dataset. |
69 | [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`. |
70 | [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. |
71 | [Filter](classROOT_1_1RDF_1_1RInterface.html#a70284a3bedc72b19610aaa91b5007ebd) | Filter the rows of the dataset. |
72 | [Range](classROOT_1_1RDF_1_1RInterface.html#a1b36b7868831de2375e061bb06cfc225) | Creates a node that filters entries based on range of entries |
73 
74 ### Actions
75 Actions are a way to produce a result out of the data. Each one is described in more detail in the reference guide.
76 
77 In the following, whenever we say an action "returns" something, we always mean it returns a smart pointer to it. Also
78 note that all actions are only executed for events that pass all preceding filters.
79 
80 Lazy actions only trigger the event loop when one of the results is accessed for the first time, making it easy to
81 produce several different results in one event loop. Instant actions trigger the event loop instantly.
82 
83 
84 | **Lazy action** | **Description** |
85 |------------------|-----------------|
86 | [Aggregate](classROOT_1_1RDF_1_1RInterface.html#ae540b00addc441f9b504cbae0ef0a24d) | Execute a user-defined accumulation operation on the processed column values. |
87 | [Book](classROOT_1_1RDF_1_1RInterface.html#a9b2f61f3333d1669e57055b9ae8be9d9) | Book execution of a custom action using a user-defined helper object. |
88 | [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). |
89 | [Count](classROOT_1_1RDF_1_1RInterface.html#a37f9e00c2ece7f53fae50b740adc1456) | Return the number of events processed. |
90 | [Display](classROOT_1_1RDF_1_1RInterface.html#aee68f4411f16f00a1d46eccb6d296f01) | Obtains the events in the dataset for the requested columns. The method returns a [RDisplay](classROOT_1_1RDF_1_1RDisplay.html) instance which can be queried to get a compressed tabular representation on the standard output or a complete representation as a string. |
91 | [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, ...). |
92 | [Graph](classROOT_1_1RDF_1_1RInterface.html#a804b466ebdbddef5c7e3400cc6b89301) | Fills a TGraph with the two columns provided. If Multithread is enabled, the order of the points may not be the one expected, it is therefore suggested to sort if before drawing. |
93 | [Histo{1D,2D,3D}](classROOT_1_1RDF_1_1RInterface.html#a247ca3aeb7ce5b95015b7fae72983055) | Fill a {one,two,three}-dimensional histogram with the processed branch values. |
94 | [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.|
95 | [Mean](classROOT_1_1RDF_1_1RInterface.html#ade6b020284f2f4fe9d3b09246b5f376a) | Return the mean of processed branch values.|
96 | [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.|
97 | [Profile{1D,2D}](classROOT_1_1RDF_1_1RInterface.html#a8ef7dc16b0e9f7bc9cfbe2d9e5de0cef) | Fill a {one,two}-dimensional profile with the branch values that passed all filters. |
98 | [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. |
99 | [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. |
100 | [StdDev](classROOT_1_1RDF_1_1RInterface.html#a482c4e4f81fe1e421c016f89cd281572) | Return the unbiased standard deviation of the processed branch values. |
101 | [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. |
102 | [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. |
103 
104 | **Instant action** | **Description** |
105 |---------------------|-----------------|
106 | [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. |
107 | [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`. |
108 | [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.|
109 
110 
111 ### Other Operations
112 
113 | **Operation** | **Description** |
114 |---------------------|-----------------|
115 | [Alias](classROOT_1_1RDF_1_1RInterface.html#a31ca327e4a192dcc05a4aac240e1a725) | Introduce an alias for a particular column name. |
116 | [GetColumnNames](classROOT_1_1RDF_1_1RInterface.html#a951fe60b74d3a9fda37df59fd1dac186) | Get the names of all the available columns of the dataset. |
117 | [GetDefinedColumnNames](classROOT_1_1RDF_1_1RInterface.html#ad5c3fab8155aae8f614735df68430c58) | Get the names of all the defined columns |
118 | [GetColumnType](classROOT_1_1RDF_1_1RInterface.html#ad3ccd813d9fed014ae6a080411c5b5a8) | Return the type of a given column as a string. |
119 | [GetColumnTypeNamesList](classROOT_1_1RDF_1_1RInterface.html#a951fe60b74d3a9fda37df59fd1dac186) | Return the list of type names of columns in the dataset. |
120 | [GetFilterNames](classROOT_1_1RDF_1_1RInterface.html#a25026681111897058299161a70ad9bb2) | Get all the filters defined. If called on a root node, all filters will be returned. For any other node, only the filters upstream of that node. |
121 | [Display](classROOT_1_1RDF_1_1RInterface.html#a652f9ab3e8d2da9335b347b540a9a941) | Provides an ASCII representation of the columns types and contents of the dataset printable by the user. |
122 | [SaveGraph](namespaceROOT_1_1RDF.html#adc17882b283c3d3ba85b1a236197c533) | Store the computation graph of an RDataFrame in graphviz format for easy inspection. |
123 | [GetNRuns](classROOT_1_1RDF_1_1RInterface.html#adfb0562a9f7732c3afb123aefa07e0df) | Get the number of event loops run by this RDataFrame instance. |
124 
125 
126 ## <a name="introduction"></a>Introduction
127 Users define their analysis as a sequence of operations to be performed on the data-frame object; the framework
128 takes care of the management of the loop over entries as well as low-level details such as I/O and parallelisation.
129 `RDataFrame` provides methods to perform most common operations required by ROOT analyses;
130 at the same time, users can just as easily specify custom code that will be executed in the event loop.
131 
132 `RDataFrame` is built with a *modular* and *flexible* workflow in mind, summarised as follows:
133 
134 1. **build a data-frame** object by specifying your data-set
135 2. **apply a series of transformations** to your data
136  1. **filter** (e.g. apply some cuts) or
137  2. **define** a new column (e.g. the result of an expensive computation on branches)
138 3. **apply actions** to the transformed data to produce results (e.g. fill a histogram)
139 
140 The following table shows how analyses based on `TTreeReader` and `TTree::Draw` translate to `RDataFrame`. Follow the
141 [crash course](#crash-course) to discover more idiomatic and flexible ways to express analyses with `RDataFrame`.
142 <table>
143 <tr>
144  <td>
145  <b>TTreeReader</b>
146  </td>
147  <td>
148  <b>ROOT::RDataFrame</b>
149  </td>
150 </tr>
151 <tr>
152  <td>
153 ~~~{.cpp}
154 TTreeReader reader("myTree", file);
155 TTreeReaderValue<A_t> a(reader, "A");
156 TTreeReaderValue<B_t> b(reader, "B");
157 TTreeReaderValue<C_t> c(reader, "C");
158 while(reader.Next()) {
159  if(IsGoodEvent(*a, *b, *c))
160  DoStuff(*a, *b, *c);
161 }
162 ~~~
163  </td>
164  <td>
165 ~~~{.cpp}
166 ROOT::RDataFrame d("myTree", file, {"A", "B", "C"});
167 d.Filter(IsGoodEvent).Foreach(DoStuff);
168 ~~~
169  </td>
170 </tr>
171 <tr>
172  <td>
173  <b>TTree::Draw</b>
174  </td>
175  <td>
176  <b>ROOT::RDataFrame</b>
177  </td>
178 </tr>
179 <tr>
180  <td>
181 ~~~{.cpp}
182 auto t = file->Get<TTree>("myTree");
183 t->Draw("x", "y > 2");
184 ~~~
185  </td>
186  <td>
187 ~~~{.cpp}
188 ROOT::RDataFrame d("myTree", file);
189 auto h = d.Filter("y > 2").Histo1D("x");
190 ~~~
191  </td>
192 </tr>
193 </table>
194 
195 ## <a name="crash-course"></a> Crash course
196 All snippets of code presented in the crash course can be executed in the ROOT interpreter. Simply precede them with
197 ~~~{.cpp}
198 using namespace ROOT; // RDataFrame's namespace
199 ~~~
200 which is omitted for brevity. The terms "column" and "branch" are used interchangeably.
201 
202 ### Creating a RDataFrame
203 RDataFrame's constructor is where the user specifies the dataset and, optionally, a default set of columns that
204 operations should work with. Here are the most common methods to construct a RDataFrame object:
205 ~~~{.cpp}
206 // single file -- all ctors are equivalent
207 TFile *f = TFile::Open("file.root");
208 auto t = f.Get<TTree>("treeName");
209 
210 RDataFrame d1("treeName", "file.root");
211 RDataFrame d2("treeName", f); // same as TTreeReader
212 RDataFrame d3(*t); // TTreeReader takes a pointer, RDF takes a reference
213 
214 // multiple files -- all ctors are equivalent
215 std::vector<std::string> files = {"file1.root", "file2.root"};
216 TChain chain("myTree");
217 chain.Add("file1.root");
218 chain.Add("file2.root");
219 
220 RDataFrame d4("myTree", {"file1.root", "file2.root"});
221 RDataFrame d5("myTree", files);
222 RDataFrame d6("myTree", "file*.root"); // see TRegexp's documentation for a list of valid regexes
223 RDataFrame d7(chain);
224 ~~~
225 Additionally, users can construct a RDataFrame specifying just an integer number. This is the number of "events" that
226 will be generated by this RDataFrame.
227 ~~~{.cpp}
228 RDataFrame d(10); // a RDF with 10 entries (and no columns/branches, for now)
229 d.Foreach([] { static int i = 0; std::cout << i++ << std::endl; }); // silly example usage: count to ten
230 ~~~
231 This is useful to generate simple data-sets on the fly: the contents of each event can be specified via the `Define`
232 transformation (explained below). For example, we have used this method to generate Pythia events (with a `Define`
233 transformation) and write them to disk in parallel (with the `Snapshot` action).
234 
235 ### Filling a histogram
236 Let's now tackle a very common task, filling a histogram:
237 ~~~{.cpp}
238 // Fill a TH1D with the "MET" branch
239 RDataFrame d("myTree", "file.root");
240 auto h = d.Histo1D("MET");
241 h->Draw();
242 ~~~
243 The first line creates a `RDataFrame` associated to the `TTree` "myTree". This tree has a branch named "MET".
244 
245 `Histo1D` is an *action*; it returns a smart pointer (a `RResultPtr` to be precise) to a `TH1D` histogram filled
246 with the `MET` of all events. If the quantity stored in the branch is a collection (e.g. a vector or an array), the
247 histogram is filled with its elements.
248 
249 You can use the objects returned by actions as if they were pointers to the desired results. There are many other
250 possible [actions](#overview), and all their results are wrapped in smart pointers; we'll see why in a minute.
251 
252 ### Applying a filter
253 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:
254 ~~~{.cpp}
255 RDataFrame d("myTree", "file.root");
256 auto c = d.Filter("MET > 4.").Count();
257 std::cout << *c << std::endl;
258 ~~~
259 The filter string (which must contain a valid c++ expression) is applied to the specified branches for each event;
260 the name and types of the columns are inferred automatically. The string expression is required to return a `bool`
261 which signals whether the event passes the filter (`true`) or not (`false`).
262 
263 You can think of your data as "flowing" through the chain of calls, being transformed, filtered and finally used to
264 perform actions. Multiple `Filter` calls can be chained one after another.
265 
266 Using string filters is nice for simple things, but they are limited to specifying the equivalent of a single return
267 statement or the body of a lambda, so it's cumbersome to use strings with more complex filters. They also add a small
268 runtime overhead, as ROOT needs to just-in-time compile the string into C++ code. When more freedom is required or
269 runtime performance is very important, a C++ callable can be specified instead (a lambda in the following snippet,
270 but it can be any kind of function or even a functor class), together with a list of branch names.
271 This snippet is analogous to the one above:
272 ~~~{.cpp}
273 RDataFrame d("myTree", "file.root");
274 auto metCut = [](double x) { return x > 4.; }; // a c++11 lambda function checking "x > 4"
275 auto c = d.Filter(metCut, {"MET"}).Count();
276 std::cout << *c << std::endl;
277 ~~~
278 
279 An example of a more complex filter expressed as a string containing C++ code is shown below
280 
281 ~~~{.cpp}
282 RDataFrame d("myTree", "file.root");
283 auto df = d.Define("p", "std::array<double, 4> p{px, py, pz}; return p;")
284  .Filter("double p2 = 0.0; for (auto&& x : p) p2 += x*x; return sqrt(p2) < 10.0;");
285 ~~~
286 
287 The code snippet above defines a column `p` that is a fixed-size array using the component column names and then
288 filters on its magnitude by looping over its elements. It must be noted that the usage of strings to define columns
289 like the one above is a major advantage when using PyROOT. However, only constants and data coming from other columns
290 in the dataset can be involved in the code passed as a string. Local variables and functions cannot be used, since
291 the interpreter will not know how to find them. When capturing local state is necessary, a C++ callable can be used.
292 
293 More information on filters and how to use them to automatically generate cutflow reports can be found [below](#Filters).
294 
295 ### Defining custom columns
296 Let's now consider the case in which "myTree" contains two quantities "x" and "y", but our analysis relies on a derived
297 quantity `z = sqrt(x*x + y*y)`. Using the `Define` transformation, we can create a new column in the data-set containing
298 the variable "z":
299 ~~~{.cpp}
300 RDataFrame d("myTree", "file.root");
301 auto sqrtSum = [](double x, double y) { return sqrt(x*x + y*y); };
302 auto zMean = d.Define("z", sqrtSum, {"x","y"}).Mean("z");
303 std::cout << *zMean << std::endl;
304 ~~~
305 `Define` creates the variable "z" by applying `sqrtSum` to "x" and "y". Later in the chain of calls we refer to
306 variables created with `Define` as if they were actual tree branches/columns, but they are evaluated on demand, at most
307 once per event. As with filters, `Define` calls can be chained with other transformations to create multiple custom
308 columns. `Define` and `Filter` transformations can be concatenated and intermixed at will.
309 
310 As with filters, it is possible to specify new columns as string expressions. This snippet is analogous to the one above:
311 ~~~{.cpp}
312 RDataFrame d("myTree", "file.root");
313 auto zMean = d.Define("z", "sqrt(x*x + y*y)").Mean("z");
314 std::cout << *zMean << std::endl;
315 ~~~
316 Again the names of the branches used in the expression and their types are inferred automatically. The string must be
317 valid c++ and is just-in-time compiled by the ROOT interpreter, cling -- the process has a small runtime overhead.
318 
319 Previously, when showing the different ways a RDataFrame can be created, we showed a constructor that only takes a
320 number of entries a parameter. In the following example we show how to combine such an "empty" `RDataFrame` with `Define`
321 transformations to create a data-set on the fly. We then save the generated data on disk using the `Snapshot` action.
322 ~~~{.cpp}
323 RDataFrame d(100); // a RDF that will generate 100 entries (currently empty)
324 int x = -1;
325 auto d_with_columns = d.Define("x", [&x] { return ++x; })
326  .Define("xx", [&x] { return x*x; });
327 d_with_columns.Snapshot("myNewTree", "newfile.root");
328 ~~~
329 This example is slightly more advanced than what we have seen so far: for starters, it makes use of lambda captures (a
330 simple way to make external variables available inside the body of c++ lambdas) to act on the same variable `x` from
331 both `Define` transformations. Secondly we have *stored* the transformed data-frame in a variable. This is always
332 possible: at each point of the transformation chain, users can store the status of the data-frame for further use (more
333 on this [below](#callgraphs)).
334 
335 You can read more about defining new columns [here](#custom-columns).
336 
337 \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."
338 
339 ### Running on a range of entries
340 It is sometimes necessary to limit the processing of the dataset to a range of entries. For this reason, the RDataFrame
341 offers the concept of ranges as a node of the RDataFrame chain of transformations; this means that filters, columns and
342 actions can be concatenated to and intermixed with `Range`s. If a range is specified after a filter, the range will act
343 exclusively on the entries passing the filter -- it will not even count the other entries! The same goes for a `Range`
344 hanging from another `Range`. Here are some commented examples:
345 ~~~{.cpp}
346 RDataFrame d("myTree", "file.root");
347 // Here we store a data-frame that loops over only the first 30 entries in a variable
348 auto d30 = d.Range(30);
349 // This is how you pick all entries from 15 onwards
350 auto d15on = d.Range(15, 0);
351 // We can specify a stride too, in this case we pick an event every 3
352 auto d15each3 = d.Range(0, 15, 3);
353 ~~~
354 Note that ranges are not available when multi-threading is enabled. More information on ranges is available
355 [here](#ranges).
356 
357 ### Executing multiple actions in the same event loop
358 As a final example let us apply two different cuts on branch "MET" and fill two different histograms with the "pt\_v" of
359 the filtered events.
360 By now, you should be able to easily understand what's happening:
361 ~~~{.cpp}
362 RDataFrame d("treeName", "file.root");
363 auto h1 = d.Filter("MET > 10").Histo1D("pt_v");
364 auto h2 = d.Histo1D("pt_v");
365 h1->Draw(); // event loop is run once here
366 h2->Draw("SAME"); // no need to run the event loop again
367 ~~~
368 `RDataFrame` executes all above actions by **running the event-loop only once**. The trick is that actions are not
369 executed at the moment they are called, but they are **lazy**, i.e. delayed until the moment one of their results is
370 accessed through the smart pointer. At that time, the event loop is triggered and *all* results are produced
371 simultaneously.
372 
373 It is therefore good practice to declare all your transformations and actions *before* accessing their results, allowing
374 `RDataFrame` to run the loop once and produce all results in one go.
375 
376 ### Going parallel
377 Let's say we would like to run the previous examples in parallel on several cores, dividing events fairly between cores.
378 The only modification required to the snippets would be the addition of this line *before* constructing the main
379 data-frame object:
380 ~~~{.cpp}
381 ROOT::EnableImplicitMT();
382 ~~~
383 Simple as that. More details are given [below](#parallel-execution).
384 
385 
386 ## <a name="python"></a>Efficient analysis in Python
387 
388 You can use `RDataFrame` in Python due to the dynamic C++/Python translation of PyROOT. In general, the interface
389 is the same as for C++, a simple example follows.
390 
391 ~~~{.python}
392 df = ROOT.RDataFrame("myTree", "myFile.root")
393 sum = df.Filter("x > 10").Sum("y")
394 print(sum.GetValue())
395 ~~~
396 
397 ### Simple usage of efficient C++ code in Python
398 
399 To perform more complex operations in the `RDataFrame` graph, e.g., in `Filter` and `Define` nodes, which don't
400 fit into a simple expression string, you can just-in-time compile such functions directly in the Python script
401 via the C++ interpreter cling. This approach has the advantage that you get the efficiency of compiled C++ code
402 combined with the convenient workflow of a Python script. See the following snippet for an example of how to
403 use a jitted C++ function from Python.
404 
405 ~~~{.python}
406 ROOT.gInterpreter.Declare("""
407 bool myFilter(float x) {
408  return x > 10;
409 }
410 """)
411 
412 df = ROOT.RDataFrame("myTree", "myFile.root")
413 sum = df.Filter("myFilter(x)").Sum("y")
414 print(sum.GetValue())
415 ~~~
416 
417 To increase the performance even further, you can also precompile a C++ library with full code optimizations
418 and load the function into the `RDataFrame` computation as follows.
419 
420 ~~~{.python}
421 ROOT.gSystem.Load("path/to/myLibrary.so") # Library with the myFilter function
422 ROOT.gInterpreter.Declare('#include "myLibrary.h"') # Header with the definition of the myFilter function
423 df = ROOT.RDataFrame("myTree", "myFile.root")
424 sum = df.Filter("myFilter(x)").Sum("y")
425 print(sum.GetValue())
426 ~~~
427 
428 ### Just-in-time compilation of Python callables with numba
429 
430 ROOT also offers the option to compile Python callables with fundamental types and arrays thereof using numba and then
431 using the function in `RDataFrame` from C++. The workflow requires the Python packages `numba` and `cffi`
432 to be installed. See the following snippet for a simple example or the full tutorial [here](pyroot004__NumbaDeclare_8py.html).
433 
434 ~~~{.python}
435 @ROOT.Numba.Declare(["float"], "bool")
436 def myFilter(x):
437  return x > 10
438 
439 df = ROOT.RDataFrame("myTree", "myFile.root")
440 sum = df.Filter("Numba::myFilter(x)").Sum("y")
441 print(sum.GetValue())
442 ~~~
443 
444 ### Conversion to numpy arrays
445 
446 Eventually, you probably would like to inspect the content of the `RDataFrame` or process the data further
447 with functionality from Python libraries. For this purpose, we provide the `AsNumpy` function, which is able
448 to provide you the columns of your `RDataFrame` as numpy arrays in Python. See a brief introduction below or
449 a full tutorial [here](df026__AsNumpyArrays_8py.html).
450 
451 ~~~{.python}
452 df = ROOT.RDataFrame("myTree", "myFile.root")
453 cols = df.Filter("x > 10").Sum("y").AsNumpy(["x", "y"])
454 print(cols["x"], cols["y"])
455 ~~~
456 
457 ## <a name="more-features"></a>More features
458 Here is a list of the most important features that have been omitted in the "Crash course" for brevity.
459 You don't need to read all these to start using `RDataFrame`, but they are useful to save typing time and runtime.
460 
461 ### Programmatically get the list of column names
462 The `GetColumnsNames()` method returns the list of valid column names for the dataset:
463 ~~~{.cpp}
464 RDataFrame d("myTree", "file.root");
465 std::vector<std::string> colNames = d.GetColumnNames();
466 ~~~
467 
468 ### Reading and manipulating collections
469 When using RDataFrame to read data from a ROOT file, users can specify that the type of a branch is `RVec<T>`
470 to indicate the branch is a c-style array, a `std::vector` or any other collection type associated to a
471 contiguous storage in memory.
472 
473 Column values of type `RVec<T>` perform no copy of the underlying array data
474 and offer a rich interface to operate on the array elements in a vectorised fashion.
475 
476 The `RVec<T>` type signals to RDataFrame that a special behaviour needs to be adopted when snapshotting
477 a dataset on disk. Indeed, if columns which are variable size C arrays are treated via the `RVec<T>`,
478 RDataFrame will correctly persistify them - if anything else is adopted, for example `std::span`, only
479 the first element of the array will be written.
480 
481 Learn more on [RVec](https://root.cern/doc/master/classROOT_1_1VecOps_1_1RVec.html).
482 
483 ### Callbacks
484 It's possible to schedule execution of arbitrary functions (callbacks) during the event loop.
485 Callbacks can be used e.g. to inspect partial results of the analysis while the event loop is running,
486 drawing a partially-filled histogram every time a certain number of new entries is processed, or event
487 displaying a progress bar while the event loop runs.
488 
489 For example one can draw an up-to-date version of a result histogram every 100 entries like this:
490 ~~~{.cpp}
491 auto h = tdf.Histo1D("x");
492 TCanvas c("c","x hist");
493 h.OnPartialResult(100, [&c](TH1D &h_) { c.cd(); h_.Draw(); c.Update(); });
494 h->Draw(); // event loop runs here, this `Draw` is executed after the event loop is finished
495 ~~~
496 
497 Callbacks are registered to a RResultPtr and must be callables that takes a reference to the result type as argument
498 and return nothing. RDataFrame will invoke registered callbacks passing partial action results as arguments to them
499 (e.g. a histogram filled with a part of the selected events).
500 
501 Read more on RResultPtr::OnPartialResult().
502 
503 ### Default branch lists
504 When constructing a `RDataFrame` object, it is possible to specify a **default column list** for your analysis, in the
505 usual form of a list of strings representing branch/column names. The default column list will be used as a fallback
506 whenever a list specific to the transformation/action is not present. RDataFrame will take as many of these columns as
507 needed, ignoring trailing extra names if present.
508 ~~~{.cpp}
509 // use "b1" and "b2" as default branches
510 RDataFrame d1("myTree", "file.root", {"b1","b2"});
511 auto h = d1.Filter([](int b1, int b2) { return b1 > b2; }) // will act on "b1" and "b2"
512  .Histo1D(); // will act on "b1"
513 
514 // just one default branch this time
515 RDataFrame d2("myTree", "file.root", {"b1"});
516 auto min = d2.Filter([](double b2) { return b2 > 0; }, {"b2"}) // we can still specify non-default branch lists
517  .Min(); // returns the minimum value of "b1" for the filtered entries
518 ~~~
519 
520 ### <a name="ImplicitColumns"></a> Implicit Columns
521 Every instance of `RDataFrame` is created with two special columns called `rdfentry_` and `rdfslot_`. The `rdfentry_`
522 column is an unsigned 64-bit integer holding the current entry number while `rdfslot_` is an unsigned 32-bit integer
523 holding the index of the current data processing slot.
524 For backwards compatibility reasons, the names `tdfentry_` and `tdfslot_` are also accepted.
525 These columns are not considered by operations such as [Cache](classROOT_1_1RDF_1_1RInterface.html#aaaa0a7bb8eb21315d8daa08c3e25f6c9)
526 or [Snapshot](classROOT_1_1RDF_1_1RInterface.html#a233b7723e498967f4340705d2c4db7f8). The _cached_ or _snapshot_ data frame
527 provides "its own" values for these columns which do not necessarily correspond to the ones of the mother data frame. This is
528 most notably the case where filters are used before deriving a cached/persistified dataframe.
529 
530 Note that in multi-thread event loops the values of `rdfentry_` _do not_ correspond to what would be the entry numbers
531 of a TChain constructed over the same set of ROOT files, as the entries are processed in an unspecified order.
532 
533 ### Branch type guessing and explicit declaration of branch types
534 C++ is a statically typed language: all types must be known at compile-time. This includes the types of the `TTree`
535 branches we want to work on. For filters, temporary columns and some of the actions, **branch types are deduced from the
536 signature** of the relevant filter function/temporary column expression/action function:
537 ~~~{.cpp}
538 // here b1 is deduced to be `int` and b2 to be `double`
539 dataFrame.Filter([](int x, double y) { return x > 0 && y < 0.; }, {"b1", "b2"});
540 ~~~
541 If we specify an incorrect type for one of the branches, an exception with an informative message will be thrown at
542 runtime, when the branch value is actually read from the `TTree`: `RDataFrame` detects type mismatches. The same would
543 happen if we swapped the order of "b1" and "b2" in the branch list passed to `Filter`.
544 
545 Certain actions, on the other hand, do not take a function as argument (e.g. `Histo1D`), so we cannot deduce the type of
546 the branch at compile-time. In this case **`RDataFrame` infers the type of the branch** from the `TTree` itself. This
547 is why we never needed to specify the branch types for all actions in the above snippets.
548 
549 When the branch type is not a common one such as `int`, `double`, `char` or `float` it is nonetheless good practice to
550 specify it as a template parameter to the action itself, like this:
551 ~~~{.cpp}
552 dataFrame.Histo1D("b1"); // OK, the type of "b1" is deduced at runtime
553 dataFrame.Min<MyNumber_t>("myObject"); // OK, "myObject" is deduced to be of type `MyNumber_t`
554 ~~~
555 
556 Deducing types at runtime requires the just-in-time compilation of the relevant actions, which has a small runtime
557 overhead, so specifying the type of the columns as template parameters to the action is good practice when performance is a goal.
558 
559 ### Generic actions
560 `RDataFrame` strives to offer a comprehensive set of standard actions that can be performed on each event. At the same
561 time, it **allows users to execute arbitrary code (i.e. a generic action) inside the event loop** through the `Foreach`
562 and `ForeachSlot` actions.
563 
564 `Foreach(f, columnList)` takes a function `f` (lambda expression, free function, functor...) and a list of columns, and
565 executes `f` on those columns for each event. The function passed must return nothing (i.e. `void`). It can be used to
566 perform actions that are not already available in the interface. For example, the following snippet evaluates the root
567 mean square of column "b":
568 ~~~{.cpp}
569 // Single-thread evaluation of RMS of column "b" using Foreach
570 double sumSq = 0.;
571 unsigned int n = 0;
572 RDataFrame d("bTree", bFilePtr);
573 d.Foreach([&sumSq, &n](double b) { ++n; sumSq += b*b; }, {"b"});
574 std::cout << "rms of b: " << std::sqrt(sumSq / n) << std::endl;
575 ~~~
576 When executing on multiple threads, users are responsible for the thread-safety of the expression passed to `Foreach`:
577 each thread will execute the expression multiple times (once per entry) in an unspecified order.
578 The code above would need to employ some resource protection mechanism to ensure non-concurrent writing of `rms`; but
579 this is probably too much head-scratch for such a simple operation.
580 
581 `ForeachSlot` can help in this situation. It is an alternative version of `Foreach` for which the function takes an
582 additional parameter besides the columns it should be applied to: an `unsigned int slot` parameter, where `slot` is a
583 number indicating which thread (0, 1, 2 , ..., poolSize - 1) the function is being run in. More specifically, RDataFrame
584 guarantees that `ForeachSlot` will invoke the user expression with different `slot` parameters for different concurrent
585 executions (there is no guarantee that a certain slot number will always correspond to a given thread id, though).
586 We can take advantage of `ForeachSlot` to evaluate a thread-safe root mean square of branch "b":
587 ~~~{.cpp}
588 // Thread-safe evaluation of RMS of branch "b" using ForeachSlot
589 ROOT::EnableImplicitMT();
590 const unsigned int nSlots = ROOT::GetThreadPoolSize();
591 std::vector<double> sumSqs(nSlots, 0.);
592 std::vector<unsigned int> ns(nSlots, 0);
593 
594 RDataFrame d("bTree", bFilePtr);
595 d.ForeachSlot([&sumSqs, &ns](unsigned int slot, double b) { sumSqs[slot] += b*b; ns[slot] += 1; }, {"b"});
596 double sumSq = std::accumulate(sumSqs.begin(), sumSqs.end(), 0.); // sum all squares
597 unsigned int n = std::accumulate(ns.begin(), ns.end(), 0); // sum all counts
598 std::cout << "rms of b: " << std::sqrt(sumSq / n) << std::endl;
599 ~~~
600 You see how we created one `double` variable for each thread in the pool, and later merged their results via
601 `std::accumulate`.
602 
603 ### Friend trees
604 Friend trees are supported by RDataFrame.
605 In order to deal with friend trees with RDataFrame, the user is required to build
606 the tree and its friends and instantiate a RDataFrame with it.
607 ~~~{.cpp}
608 TTree t([...]);
609 TTree ft([...]);
610 t.AddFriend(ft, "myFriend");
611 
612 RDataFrame d(t);
613 auto f = d.Filter("myFriend.MyCol == 42");
614 ~~~
615 
616 ### Reading file formats different from ROOT's
617 RDataFrame can be interfaced with RDataSources. The RDataSource interface defines an API that RDataFrame can use to read arbitrary data formats.
618 
619 A concrete RDataSource implementation (i.e. a class that inherits from RDataSource and implements all of its pure
620 methods) provides an adaptor that RDataFrame can leverage to read any kind of tabular data formats.
621 RDataFrame calls into RDataSource to retrieve information about the data, retrieve (thread-local) readers or "cursors" for selected columns
622 and to advance the readers to the desired data entry.
623 Some predefined RDataSources are natively provided by ROOT such as the `RCsvDS` which allows to read comma separated files:
624 ~~~{.cpp}
625 auto tdf = ROOT::RDF::MakeCsvDataFrame("MuRun2010B.csv");
626 auto filteredEvents =
627  tdf.Filter("Q1 * Q2 == -1")
628  .Define("m", "sqrt(pow(E1 + E2, 2) - (pow(px1 + px2, 2) + pow(py1 + py2, 2) + pow(pz1 + pz2, 2)))");
629 auto h = filteredEvents.Histo1D("m");
630 h->Draw();
631 ~~~
632 
633 ### <a name="callgraphs"></a>Call graphs (storing and reusing sets of transformations)
634 **Sets of transformations can be stored as variables** and reused multiple times to create **call graphs** in which
635 several paths of filtering/creation of columns are executed simultaneously; we often refer to this as "storing the
636 state of the chain".
637 
638 This feature can be used, for example, to create a temporary column once and use it in several subsequent filters or
639 actions, or to apply a strict filter to the data-set *before* executing several other transformations and actions,
640 effectively reducing the amount of events processed.
641 
642 Let's try to make this clearer with a commented example:
643 ~~~{.cpp}
644 // build the data-frame and specify a default column list
645 RDataFrame d(treeName, filePtr, {"var1", "var2", "var3"});
646 
647 // apply a cut and save the state of the chain
648 auto filtered = d.Filter(myBigCut);
649 
650 // plot branch "var1" at this point of the chain
651 auto h1 = filtered.Histo1D("var1");
652 
653 // create a new branch "vec" with a vector extracted from a complex object (only for filtered entries)
654 // and save the state of the chain
655 auto newBranchFiltered = filtered.Define("vec", [](const Obj& o) { return o.getVector(); }, {"obj"});
656 
657 // apply a cut and fill a histogram with "vec"
658 auto h2 = newBranchFiltered.Filter(cut1).Histo1D("vec");
659 
660 // apply a different cut and fill a new histogram
661 auto h3 = newBranchFiltered.Filter(cut2).Histo1D("vec");
662 
663 // Inspect results
664 h2->Draw(); // first access to an action result: run event-loop!
665 h3->Draw("SAME"); // event loop does not need to be run again here..
666 std::cout << "Entries in h1: " << h1->GetEntries() << std::endl; // ..or here
667 ~~~
668 `RDataFrame` detects when several actions use the same filter or the same temporary column, and **only evaluates each
669 filter or temporary column once per event**, regardless of how many times that result is used down the call graph.
670 Objects read from each column are **built once and never copied**, for maximum efficiency.
671 When "upstream" filters are not passed, subsequent filters, temporary column expressions and actions are not evaluated,
672 so it might be advisable to put the strictest filters first in the chain.
673 
674 ### <a name="representgraph"></a>Printing the computation graph
675 It is possible to print the computation graph from any node to obtain a dot representation either on the standard output
676 or in a file.
677 
678 Invoking the function ROOT::RDF::SaveGraph() on any node that is not the head node, the computation graph of the branch
679 the node belongs to is printed. By using the head node, the entire computation graph is printed.
680 
681 Following there is an example of usage:
682 ~~~{.cpp}
683 // First, a sample computational graph is built
684 ROOT::RDataFrame df("tree", "f.root");
685 
686 auto df2 = df.Define("x", []() { return 1; })
687  .Filter("col0 % 1 == col0")
688  .Filter([](int b1) { return b1 <2; }, {"cut1"})
689  .Define("y", []() { return 1; });
690 
691 auto count = df2.Count();
692 
693 // Prints the graph to the rd1.dot file in the current directory
694 ROOT::RDF::SaveGraph(rd1, "./mydot.dot");
695 // Prints the graph to standard output
696 ROOT::RDF::SaveGraph(rd1);
697 ~~~
698 
699 ### RDataFrame variables as function arguments and return values
700 RDataFrame variables/nodes are relatively cheap to copy and it's possible to both pass them to (or move them into)
701 functions and to return them from functions. However, in general each dataframe node will have a different C++ type,
702 which includes all available compile-time information about what that node does. One way to cope with this complication
703 is to use template functions and/or C++14 auto return types:
704 ~~~{.cpp}
705 template <typename RDF>
706 auto ApplySomeFilters(RDF df)
707 {
708  return df.Filter("x > 0").Filter([](int y) { return y < 0; }, {"y"});
709 }
710 ~~~
711 
712 A possibly simpler, C++11-compatible alternative is to take advantage of the fact that any dataframe node can be
713 converted to the common type ROOT::RDF::RNode:
714 ~~~{.cpp}
715 // a function that conditionally adds a Range to a RDataFrame node.
716 RNode MaybeAddRange(RNode df, bool mustAddRange)
717 {
718  return mustAddRange ? df.Range(1) : df;
719 }
720 // use as :
721 ROOT::RDataFrame df(10);
722 auto maybeRangedDF = MaybeAddRange(df, true);
723 ~~~
724 
725 The conversion to ROOT::RDF::RNode is cheap, but it will introduce an extra virtual call during the RDataFrame event
726 loop (in most cases, the resulting performance impact should be negligible).
727 
728 As a final note, remember that RDataFrame actions do not return another dataframe, but a RResultPtr<T>, where T is the
729 type of the result of the action.
730 
731 Read more on this topic [here](https://root.cern.ch/doc/master/classROOT_1_1RDF_1_1RInterface.html#a6909f04c05723de79f97a14b092318b1).
732 
733 ## <a name="transformations"></a>Transformations
734 ### <a name="Filters"></a> Filters
735 A filter is defined through a call to `Filter(f, columnList)`. `f` can be a function, a lambda expression, a functor
736 class, or any other callable object. It must return a `bool` signalling whether the event has passed the selection
737 (`true`) or not (`false`). It must perform "read-only" actions on the columns, and should not have side-effects (e.g.
738 modification of an external or static variable) to ensure correct results when implicit multi-threading is active.
739 
740 `RDataFrame` only evaluates filters when necessary: if multiple filters are chained one after another, they are executed
741 in order and the first one returning `false` causes the event to be discarded and triggers the processing of the next
742 entry. If multiple actions or transformations depend on the same filter, that filter is not executed multiple times for
743 each entry: after the first access it simply serves a cached result.
744 
745 #### <a name="named-filters-and-cutflow-reports"></a>Named filters and cutflow reports
746 An optional string parameter `name` can be passed to the `Filter` method to create a **named filter**. Named filters
747 work as usual, but also keep track of how many entries they accept and reject.
748 
749 Statistics are retrieved through a call to the `Report` method:
750 
751 - when `Report` is called on the main `RDataFrame` object, it returns a RResultPtr<RCutFlowReport> relative to all
752 named filters declared up to that point
753 - when called on a specific node (e.g. the result of a `Define` or `Filter`), it returns a RResultPtr<RCutFlowReport>
754 relative all named filters in the section of the chain between the main `RDataFrame` and that node (included).
755 
756 Stats are stored in the same order as named filters have been added to the graph, and *refer to the latest event-loop*
757 that has been run using the relevant `RDataFrame`.
758 
759 ### <a name="ranges"></a>Ranges
760 When `RDataFrame` is not being used in a multi-thread environment (i.e. no call to `EnableImplicitMT` was made),
761 `Range` transformations are available. These act very much like filters but instead of basing their decision on
762 a filter expression, they rely on `begin`,`end` and `stride` parameters.
763 
764 - `begin`: initial entry number considered for this range.
765 - `end`: final entry number (excluded) considered for this range. 0 means that the range goes until the end of the dataset.
766 - `stride`: process one entry of the [begin, end) range every `stride` entries. Must be strictly greater than 0.
767 
768 The actual number of entries processed downstream of a `Range` node will be `(end - begin)/stride` (or less if less
769 entries than that are available).
770 
771 Note that ranges act "locally", not based on the global entry count: `Range(10,50)` means "skip the first 10 entries
772 *that reach this node*, let the next 40 entries pass, then stop processing". If a range node hangs from a filter node,
773 and the range has a `begin` parameter of 10, that means the range will skip the first 10 entries *that pass the
774 preceding filter*.
775 
776 Ranges allow "early quitting": if all branches of execution of a functional graph reached their `end` value of
777 processed entries, the event-loop is immediately interrupted. This is useful for debugging and quick data explorations.
778 
779 ### <a name="custom-columns"></a> Custom columns
780 Custom columns are created by invoking `Define(name, f, columnList)`. As usual, `f` can be any callable object
781 (function, lambda expression, functor class...); it takes the values of the columns listed in `columnList` (a list of
782 strings) as parameters, in the same order as they are listed in `columnList`. `f` must return the value that will be
783 assigned to the temporary column.
784 
785 A new variable is created called `name`, accessible as if it was contained in the dataset from subsequent
786 transformations/actions.
787 
788 Use cases include:
789 - caching the results of complex calculations for easy and efficient multiple access
790 - extraction of quantities of interest from complex objects
791 - branch aliasing, i.e. changing the name of a branch
792 
793 An exception is thrown if the `name` of the new column/branch is already in use for another branch in the `TTree`.
794 
795 It is also possible to specify the quantity to be stored in the new temporary column as a C++ expression with the method
796 `Define(name, expression)`. For example this invocation
797 
798 ~~~{.cpp}
799 tdf.Define("pt", "sqrt(px*px + py*py)");
800 ~~~
801 
802 will create a new column called "pt" the value of which is calculated starting from the columns px and py. The system
803 builds a just-in-time compiled function starting from the expression after having deduced the list of necessary branches
804 from the names of the variables specified by the user.
805 
806 #### Custom columns as function of slot and entry number
807 
808 It is possible to create custom columns also as a function of the processing slot and entry numbers. The methods that can
809 be invoked are:
810 - `DefineSlot(name, f, columnList)`. In this case the callable f has this signature `R(unsigned int, T1, T2, ...)`: the
811 first parameter is the slot number which ranges from 0 to ROOT::GetThreadPoolSize() - 1.
812 - `DefineSlotEntry(name, f, columnList)`. In this case the callable f has this signature `R(unsigned int, ULong64_t,
813 T1, T2, ...)`: the first parameter is the slot number while the second one the number of the entry being processed.
814 
815 ## <a name="actions"></a>Actions
816 ### Instant and lazy actions
817 Actions can be **instant** or **lazy**. Instant actions are executed as soon as they are called, while lazy actions are
818 executed whenever the object they return is accessed for the first time. As a rule of thumb, actions with a return value
819 are lazy, the others are instant.
820 
821 ## <a name="parallel-execution"></a>Parallel execution
822 As pointed out before in this document, `RDataFrame` can transparently perform multi-threaded event loops to speed up
823 the execution of its actions. Users have to call `ROOT::EnableImplicitMT()` *before* constructing the `RDataFrame`
824 object to indicate that it should take advantage of a pool of worker threads. **Each worker thread processes a distinct
825 subset of entries**, and their partial results are merged before returning the final values to the user.
826 More specifically, the dataset will be divided in batches of entries, and threads will divide among themselves the
827 processing of these batches. There are no guarantees on the order the batches are processed, i.e. no guarantees in the
828 order entries of the dataset are processed. Note that this in turn means that, for multi-thread event loops, there is no
829 guarantee on the order in which `Snapshot` will _write_ entries: they could be scrambled with respect to the input dataset.
830 
831 \warning RDataFrame will by default start as many threads as the hardware supports, using up **all** the resources on
832 a machine. On a worker node of *e.g.* a batch cluster, this might not be desired if the machine is shared with other
833 users. Therefore, **when running on shared computing resources**, use
834 ```
835 ROOT::EnableImplicitMT(i)
836 ```
837 replacing `i` with the number of CPUs/slots that were allocated for this job.
838 
839 ### Thread-safety of user-defined expressions
840 RDataFrame operations such as `Histo1D` or `Snapshot` are guaranteed to work correctly in multi-thread event loops.
841 User-defined expressions, such as strings or lambdas passed to `Filter`, `Define`, `Foreach`, `Reduce` or `Aggregate`
842 will have to be thread-safe, i.e. it should be possible to call them concurrently from different threads.
843 
844 Note that simple `Filter` and `Define` transformations will inherently satisfy this requirement: `Filter`/`Define`
845 expressions will often be *pure* in the functional programming sense (no side-effects, no dependency on external state),
846 which eliminates all risks of race conditions.
847 
848 In order to facilitate writing of thread-safe operations, some RDataFrame features such as `Foreach`, `Define` or `OnPartialResult`
849 offer thread-aware counterparts (`ForeachSlot`, `DefineSlot`, `OnPartialResultSlot`): their only difference is that they
850 will pass an extra `slot` argument (an unsigned integer) to the user-defined expression. When calling user-defined code
851 concurrently, `RDataFrame` guarantees that different threads will employ different values of the `slot` parameter,
852 where `slot` will be a number between 0 and `ROOT::GetThreadPoolSize() - 1`.
853 In other words, within a slot, computation runs sequentially and events are processed sequentially.
854 Note that the same slot might be associated to different threads over the course of a single event loop, but two threads
855 will never receive the same slot at the same time.
856 This extra parameter might facilitate writing safe parallel code by having each thread write/modify a different
857 *processing slot*, e.g. a different element of a list. See [here](#generic-actions) for an example usage of `ForeachSlot`.
858 
859 <a name="reference"></a>
860 */
861 // clang-format on
862 
863 namespace ROOT {
864 
866 using ColumnNamesPtr_t = std::shared_ptr<const ColumnNames_t>;
867 
868 namespace RDFInternal = ROOT::Internal::RDF;
869 
870 ////////////////////////////////////////////////////////////////////////////
871 /// \brief Build the dataframe
872 /// \param[in] treeName Name of the tree contained in the directory
873 /// \param[in] dirPtr TDirectory where the tree is stored, e.g. a TFile.
874 /// \param[in] defaultBranches Collection of default branches.
875 ///
876 /// The default branches are looked at in case no branch is specified in the
877 /// booking of actions or transformations.
878 /// See RInterface for the documentation of the methods available.
879 RDataFrame::RDataFrame(std::string_view treeName, TDirectory *dirPtr, const ColumnNames_t &defaultBranches)
880  : RInterface(std::make_shared<RDFDetail::RLoopManager>(nullptr, defaultBranches))
881 {
882  if (!dirPtr) {
883  auto msg = "Invalid TDirectory!";
884  throw std::runtime_error(msg);
885  }
886  const std::string treeNameInt(treeName);
887  auto tree = static_cast<TTree *>(dirPtr->Get(treeNameInt.c_str()));
888  if (!tree) {
889  auto msg = "Tree \"" + treeNameInt + "\" cannot be found!";
890  throw std::runtime_error(msg);
891  }
892  GetProxiedPtr()->SetTree(std::shared_ptr<TTree>(tree, [](TTree *) {}));
893 }
894 
895 ////////////////////////////////////////////////////////////////////////////
896 /// \brief Build the dataframe
897 /// \param[in] treeName Name of the tree contained in the directory
898 /// \param[in] filenameglob TDirectory where the tree is stored, e.g. a TFile.
899 /// \param[in] defaultBranches Collection of default branches.
900 ///
901 /// The filename globbing supports the same type of expressions as TChain::Add().
902 /// The default branches are looked at in case no branch is specified in the
903 /// booking of actions or transformations.
904 /// See RInterface for the documentation of the methods available.
905 RDataFrame::RDataFrame(std::string_view treeName, std::string_view filenameglob, const ColumnNames_t &defaultBranches)
906  : RInterface(std::make_shared<RDFDetail::RLoopManager>(nullptr, defaultBranches))
907 {
908  const std::string treeNameInt(treeName);
909  const std::string filenameglobInt(filenameglob);
910  auto chain = std::make_shared<TChain>(treeNameInt.c_str());
911  chain->Add(filenameglobInt.c_str());
912  GetProxiedPtr()->SetTree(chain);
913 }
914 
915 ////////////////////////////////////////////////////////////////////////////
916 /// \brief Build the dataframe
917 /// \param[in] treeName Name of the tree contained in the directory
918 /// \param[in] fileglobs Collection of file names of filename globs
919 /// \param[in] defaultBranches Collection of default branches.
920 ///
921 /// The filename globbing supports the same type of expressions as TChain::Add().
922 /// The default branches are looked at in case no branch is specified in the booking of actions or transformations.
923 /// See RInterface for the documentation of the methods available.
924 RDataFrame::RDataFrame(std::string_view treeName, const std::vector<std::string> &fileglobs,
925  const ColumnNames_t &defaultBranches)
926  : RInterface(std::make_shared<RDFDetail::RLoopManager>(nullptr, defaultBranches))
927 {
928  std::string treeNameInt(treeName);
929  auto chain = std::make_shared<TChain>(treeNameInt.c_str());
930  for (auto &f : fileglobs)
931  chain->Add(f.c_str());
932  GetProxiedPtr()->SetTree(chain);
933 }
934 
935 ////////////////////////////////////////////////////////////////////////////
936 /// \brief Build the dataframe
937 /// \param[in] tree The tree or chain to be studied.
938 /// \param[in] defaultBranches Collection of default column names to fall back to when none is specified.
939 ///
940 /// The default branches are looked at in case no branch is specified in the
941 /// booking of actions or transformations.
942 /// See RInterface for the documentation of the methods available.
944  : RInterface(std::make_shared<RDFDetail::RLoopManager>(&tree, defaultBranches))
945 {
946 }
947 
948 //////////////////////////////////////////////////////////////////////////
949 /// \brief Build a dataframe that generates numEntries entries.
950 /// \param[in] numEntries The number of entries to generate.
951 ///
952 /// An empty-source dataframe constructed with a number of entries will
953 /// generate those entries on the fly when some action is triggered,
954 /// and it will do so for all the previously-defined temporary branches.
955 /// See RInterface for the documentation of the methods available.
957  : RInterface(std::make_shared<RDFDetail::RLoopManager>(numEntries))
958 
959 {
960 }
961 
962 //////////////////////////////////////////////////////////////////////////
963 /// \brief Build dataframe associated to datasource.
964 /// \param[in] ds The data-source object.
965 /// \param[in] defaultBranches Collection of default column names to fall back to when none is specified.
966 ///
967 /// A dataframe associated to a datasource will query it to access column values.
968 /// See RInterface for the documentation of the methods available.
969 RDataFrame::RDataFrame(std::unique_ptr<ROOT::RDF::RDataSource> ds, const ColumnNames_t &defaultBranches)
970  : RInterface(std::make_shared<RDFDetail::RLoopManager>(std::move(ds), defaultBranches))
971 {
972 }
973 
974 } // namespace ROOT
975 
976 namespace cling {
977 //////////////////////////////////////////////////////////////////////////
978 /// Print a RDataFrame at the prompt
979 std::string printValue(ROOT::RDataFrame *tdf)
980 {
981  auto &df = *tdf->GetLoopManager();
982  auto *tree = df.GetTree();
983  auto defBranches = df.GetDefaultColumnNames();
984 
985  std::ostringstream ret;
986  if (tree) {
987  ret << "A data frame built on top of the " << tree->GetName() << " dataset.";
988  if (!defBranches.empty()) {
989  if (defBranches.size() == 1)
990  ret << "\nDefault branch: " << defBranches[0];
991  else {
992  ret << "\nDefault branches:\n";
993  for (auto &&branch : defBranches) {
994  ret << " - " << branch << "\n";
995  }
996  }
997  }
998  } else if (auto ds = tdf->fDataSource) {
999  ret << "A data frame associated to the data source \"" << cling::printValue(ds) << "\"";
1000  } else {
1001  ret << "An empty data frame that will create " << df.GetNEmptyEntries() << " entries\n";
1002  }
1003 
1004  return ret.str();
1005 }
1006 } // namespace cling
ROOT::RDataFrame::ColumnNames_t
RDFDetail::ColumnNames_t ColumnNames_t
Definition: RDataFrame.hxx:44
f
#define f(i)
Definition: RSha256.hxx:104
TDirectory.h
tree
Definition: tree.py:1
ROOT::RDF::RInterface::fDataSource
RDataSource * fDataSource
Non-owning pointer to a data-source object. Null if no data-source. RLoopManager has ownership of the...
Definition: RInterface.hxx:105
string_view
basic_string_view< char > string_view
Definition: libcpp_string_view.h:785
TTree
A TTree represents a columnar dataset.
Definition: TTree.h:79
extract_docstrings.ds
ds
Definition: extract_docstrings.py:40
ROOT::RDataFrame::RDataFrame
RDataFrame(std::string_view treeName, std::string_view filenameglob, const ColumnNames_t &defaultBranches={})
Build the dataframe.
Definition: RDataFrame.cxx:905
ROOT::ColumnNamesPtr_t
std::shared_ptr< const ColumnNames_t > ColumnNamesPtr_t
Definition: RDataFrame.cxx:866
ROOT::RDataFrame
ROOT's RDataFrame offers a high level interface for analyses of data stored in TTrees,...
Definition: RDataFrame.hxx:42
ROOT::Detail::RDF
Definition: GraphUtils.hxx:28
TDirectory::Get
virtual TObject * Get(const char *namecycle)
Return pointer to object identified by namecycle.
Definition: TDirectory.cxx:806
ROOT::Detail::RDF::ColumnNames_t
std::vector< std::string > ColumnNames_t
Definition: RLoopManager.hxx:53
RDataFrame.hxx
RDataSource.hxx
TChain.h
ULong64_t
unsigned long long ULong64_t
Definition: RtypesCore.h:74
ROOT::Internal::RDF
Definition: RArrowDS.hxx:15
TDirectory
Describe directory structure in memory.
Definition: TDirectory.h:40
ROOT::RDF::RInterface< RDFDetail::RLoopManager >::GetProxiedPtr
const std::shared_ptr< RDFDetail::RLoopManager > & GetProxiedPtr() const
Definition: RInterface.hxx:2459
ROOT::RDF::RInterface::GetLoopManager
RLoopManager * GetLoopManager() const
Definition: RInterface.hxx:2457
ROOT
VSD Structures.
Definition: StringConv.hxx:21
ROOT::Detail::RDF::RLoopManager
The head node of a RDF computation graph.
Definition: RLoopManager.hxx:57