Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
df016_vecOps.C
Go to the documentation of this file.
1/// \file
2/// \ingroup tutorial_dataframe
3/// \notebook -draw
4/// Process collections in RDataFrame with the help of RVec.
5///
6/// This tutorial shows the potential of the VecOps approach for treating collections
7/// stored in datasets, a situation very common in HEP data analysis.
8///
9/// \macro_code
10/// \macro_image
11///
12/// \date February 2018
13/// \author Danilo Piparo (CERN)
14
15using namespace ROOT;
16
17int df016_vecOps()
18{
19 // We re-create a set of points in a square.
20 // This is a technical detail, just to create a dataset to play with!
21 auto unifGen = [](double) { return gRandom->Uniform(-1.0, 1.0); };
22 auto vGen = [&](int len) {
23 RVecD v(len);
24 std::transform(v.begin(), v.end(), v.begin(), unifGen);
25 return v;
26 };
27 RDataFrame d(1024);
28 auto d0 = d.Define("len", []() { return (int)gRandom->Uniform(0, 16); })
29 .Define("x", vGen, {"len"})
30 .Define("y", vGen, {"len"});
31
32 // Now we have in our hands d, a RDataFrame with two columns, x and y, which
33 // hold collections of coordinates. The sizes of these collections vary.
34 // Let's now define radii from the x and y coordinates. We'll do it treating
35 // the collections stored in the columns without looping on the individual elements.
36 auto d1 = d0.Define("r", "sqrt(x*x + y*y)");
37
38 // Now we want to plot 2 quarters of a ring with radii .5 and 1.
39 // Note how the cuts are performed on RVecs, comparing them with integers and
40 // among themselves.
41 auto ring_h = d1.Define("rInFig", "r > .5 && r < 1 && x*y < 0")
42 .Define("yFig", "y[rInFig]")
43 .Define("xFig", "x[rInFig]")
44 .Histo2D({"fig", "Two quarters of a ring", 64, -1.1, 1.1, 64, -1.1, 1.1}, "xFig", "yFig");
45
46 auto cring = new TCanvas();
47 ring_h->DrawCopy("Colz");
48
49 return 0;
50}
#define d(i)
Definition RSha256.hxx:102
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char Pixmap_t Pixmap_t PictureAttributes_t attr const char char ret_data h unsigned char height h Atom_t Int_t ULong_t ULong_t unsigned char prop_list Atom_t Atom_t Atom_t Time_t UChar_t len
R__EXTERN TRandom * gRandom
Definition TRandom.h:62
ROOT's RDataFrame offers a modern, high-level interface for analysis of data stored in TTree ,...
The Canvas class.
Definition TCanvas.h:23
virtual Double_t Uniform(Double_t x1=1)
Returns a uniform deviate on the interval (0, x1).
Definition TRandom.cxx:682
tbb::task_arena is an alias of tbb::interface7::task_arena, which doesn't allow to forward declare tb...