Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
df017_vecOpsHEP.C
Go to the documentation of this file.
1/// \file
2/// \ingroup tutorial_dataframe
3/// \notebook -draw
4/// Use RVecs to plot the transverse momentum of selected particles.
5///
6/// This tutorial shows how VecOps can be used to slim down the programming
7/// model typically adopted in HEP for analysis.
8/// In this case we have a dataset containing the kinematic properties of
9/// particles stored in individual arrays.
10/// We want to plot the transverse momentum of these particles if the energy is
11/// greater than 100.
12///
13/// \macro_code
14/// \macro_image
15///
16/// \date March 2018
17/// \authors Danilo Piparo (CERN), Andre Vieira Silva
18
19auto filename = gROOT->GetTutorialDir() + "/dataframe/df017_vecOpsHEP.root";
20auto treename = "myDataset";
21using doubles = ROOT::VecOps::RVec<double>;
22using RDF = ROOT::RDataFrame;
23
24void WithTTreeReader()
25{
26 TFile f(filename);
27 TTreeReader tr(treename, &f);
28 TTreeReaderArray<double> px(tr, "px");
29 TTreeReaderArray<double> py(tr, "py");
31
32 TH1F h("pt", "pt", 16, 0, 4);
33
34 while (tr.Next()) {
35 for (auto i=0U;i < px.GetSize(); ++i) {
36 if (E[i] > 100) h.Fill(sqrt(px[i]*px[i] + py[i]*py[i]));
37 }
38 }
39 h.DrawCopy();
40}
41
42void WithRDataFrame()
43{
44 RDF f(treename, filename.Data());
45 auto CalcPt = [](doubles &px, doubles &py, doubles &E) {
46 doubles v;
47 for (auto i=0U;i < px.size(); ++i) {
48 if (E[i] > 100) {
49 v.emplace_back(sqrt(px[i]*px[i] + py[i]*py[i]));
50 }
51 }
52 return v;
53 };
54 f.Define("pt", CalcPt, {"px", "py", "E"})
55 .Histo1D<doubles>({"pt", "pt", 16, 0, 4}, "pt")->DrawCopy();
56}
57
58void WithRDataFrameVecOps()
59{
60 RDF f(treename, filename.Data());
61 auto CalcPt = [](doubles &px, doubles &py, doubles &E) {
62 auto pt = sqrt(px*px + py*py);
63 return pt[E>100];
64 };
65 f.Define("good_pt", CalcPt, {"px", "py", "E"})
66 .Histo1D<doubles>({"pt", "pt", 16, 0, 4}, "good_pt")->DrawCopy();
67}
68
69void WithRDataFrameVecOpsJit()
70{
71 RDF f(treename, filename.Data());
72 f.Define("good_pt", "sqrt(px*px + py*py)[E>100]")
73 .Histo1D({"pt", "pt", 16, 0, 4}, "good_pt")->DrawCopy();
74}
75
76void df017_vecOpsHEP()
77{
78 // We plot four times the same quantity, the key is to look into the implementation
79 // of the functions above
80 auto c = new TCanvas();
81 c->Divide(2,2);
82 c->cd(1);
83 WithTTreeReader();
84 c->cd(2);
85 WithRDataFrame();
86 c->cd(3);
87 WithRDataFrameVecOps();
88 c->cd(4);
89 WithRDataFrameVecOpsJit();
90}
#define f(i)
Definition RSha256.hxx:104
#define c(i)
Definition RSha256.hxx:101
#define h(i)
Definition RSha256.hxx:106
double sqrt(double)
#define gROOT
Definition TROOT.h:406
ROOT's RDataFrame offers a high level interface for analyses of data stored in TTrees,...
A "std::vector"-like collection of values implementing handy operation to analyse them.
Definition RVec.hxx:296
The Canvas class.
Definition TCanvas.h:23
A ROOT file is a suite of consecutive data records (TKey instances) with a well defined format.
Definition TFile.h:54
1-D histogram with a float per channel (see TH1 documentation)}
Definition TH1.h:575
An interface for reading collections stored in ROOT columnar datasets.
A simple, robust and fast interface to read values from ROOT columnar datasets such as TTree,...
Definition TTreeReader.h:44
TPaveText * pt
constexpr Double_t E()
Base of natural log:
Definition TMath.h:96