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";
21
22using namespace ROOT;
23
24
25void WithTTreeReader()
26{
27 TFile f(filename);
28 TTreeReader tr(treename, &f);
29 TTreeReaderArray<double> px(tr, "px");
30 TTreeReaderArray<double> py(tr, "py");
32
33 TH1F h("pt", "pt", 16, 0, 4);
34
35 while (tr.Next()) {
36 for (auto i=0U;i < px.GetSize(); ++i) {
37 if (E[i] > 100) h.Fill(sqrt(px[i]*px[i] + py[i]*py[i]));
38 }
39 }
40 h.DrawCopy();
41}
42
43void WithRDataFrame()
44{
45 RDataFrame f(treename, filename.Data());
46 auto CalcPt = [](RVecD &px, RVecD &py, RVecD &E) {
47 RVecD v;
48 for (auto i=0U;i < px.size(); ++i) {
49 if (E[i] > 100) {
50 v.emplace_back(sqrt(px[i]*px[i] + py[i]*py[i]));
51 }
52 }
53 return v;
54 };
55 f.Define("pt", CalcPt, {"px", "py", "E"})
56 .Histo1D<RVecD>({"pt", "pt", 16, 0, 4}, "pt")->DrawCopy();
57}
58
59void WithRDataFrameVecOps()
60{
61 RDataFrame f(treename, filename.Data());
62 auto CalcPt = [](RVecD &px, RVecD &py, RVecD &E) {
63 auto pt = sqrt(px*px + py*py);
64 return pt[E>100];
65 };
66 f.Define("good_pt", CalcPt, {"px", "py", "E"})
67 .Histo1D<RVecD>({"pt", "pt", 16, 0, 4}, "good_pt")->DrawCopy();
68}
69
70void WithRDataFrameVecOpsJit()
71{
72 RDataFrame f(treename, filename.Data());
73 f.Define("good_pt", "sqrt(px*px + py*py)[E>100]")
74 .Histo1D({"pt", "pt", 16, 0, 4}, "good_pt")->DrawCopy();
75}
76
77void df017_vecOpsHEP()
78{
79 // We plot four times the same quantity, the key is to look into the implementation
80 // of the functions above
81 auto c = new TCanvas();
82 c->Divide(2,2);
83 c->cd(1);
84 WithTTreeReader();
85 c->cd(2);
86 WithRDataFrame();
87 c->cd(3);
88 WithRDataFrameVecOps();
89 c->cd(4);
90 WithRDataFrameVecOpsJit();
91}
#define f(i)
Definition RSha256.hxx:104
#define c(i)
Definition RSha256.hxx:101
#define h(i)
Definition RSha256.hxx:106
#define gROOT
Definition TROOT.h:404
reference emplace_back(ArgTypes &&...Args)
Definition RVec.hxx:892
ROOT's RDataFrame offers a high level interface for analyses of data stored in TTree,...
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
VecExpr< UnaryOp< Sqrt< T >, VecExpr< A, T, D >, T >, T, D > sqrt(const VecExpr< A, T, D > &rhs)
tbb::task_arena is an alias of tbb::interface7::task_arena, which doesn't allow to forward declare tb...
constexpr Double_t E()
Base of natural log:
Definition TMath.h:96