Use RVecs to plot the transverse momentum of selected particles.
This tutorial shows how VecOps can be used to slim down the programming model typically adopted in HEP for analysis. In this case we have a dataset containing the kinematic properties of particles stored in individual arrays. We want to plot the transverse momentum of these particles if the energy is greater than 100 MeV.
Author: Danilo Piparo (CERN), Andre Vieira Silva
This notebook tutorial was automatically generated with ROOTBOOK-izer from the macro found in the ROOT repository on Thursday, March 23, 2023 at 10:44 AM.
auto filename = gROOT->GetTutorialDir() + "/dataframe/df017_vecOpsHEP.root";
auto treename = "myDataset";
using namespace ROOT;
Definition of a helper function:
%%cpp -d
void WithTTreeReader()
{
TFile f(filename);
TTreeReader tr(treename, &f);
TTreeReaderArray<double> px(tr, "px");
TTreeReaderArray<double> py(tr, "py");
TTreeReaderArray<double> E(tr, "E");
TH1F h("pt", "pt", 16, 0, 4);
while (tr.Next()) {
for (auto i=0U;i < px.GetSize(); ++i) {
if (E[i] > 100) h.Fill(sqrt(px[i]*px[i] + py[i]*py[i]));
}
}
h.DrawCopy();
}
Definition of a helper function:
%%cpp -d
void WithRDataFrame()
{
RDataFrame f(treename, filename.Data());
auto CalcPt = [](RVecD &px, RVecD &py, RVecD &E) {
RVecD v;
for (auto i=0U;i < px.size(); ++i) {
if (E[i] > 100) {
v.emplace_back(sqrt(px[i]*px[i] + py[i]*py[i]));
}
}
return v;
};
f.Define("pt", CalcPt, {"px", "py", "E"})
.Histo1D<RVecD>({"pt", "pt", 16, 0, 4}, "pt")->DrawCopy();
}
input_line_53:3:3: error: unknown type name 'RDataFrame' RDataFrame f(treename, filename.Data()); ^ input_line_53:4:21: error: unknown type name 'RVecD' auto CalcPt = [](RVecD &px, RVecD &py, RVecD &E) { ^ input_line_53:4:32: error: unknown type name 'RVecD' auto CalcPt = [](RVecD &px, RVecD &py, RVecD &E) { ^ input_line_53:4:43: error: unknown type name 'RVecD' auto CalcPt = [](RVecD &px, RVecD &py, RVecD &E) { ^ input_line_53:5:7: error: unknown type name 'RVecD' RVecD v; ^ input_line_53:14:14: error: use of undeclared identifier 'RVecD' .Histo1D<RVecD>({"pt", "pt", 16, 0, 4}, "pt")->DrawCopy(); ^
Definition of a helper function:
%%cpp -d
void WithRDataFrameVecOps()
{
RDataFrame f(treename, filename.Data());
auto CalcPt = [](RVecD &px, RVecD &py, RVecD &E) {
auto pt = sqrt(px*px + py*py);
return pt[E>100];
};
f.Define("good_pt", CalcPt, {"px", "py", "E"})
.Histo1D<RVecD>({"pt", "pt", 16, 0, 4}, "good_pt")->DrawCopy();
}
input_line_54:3:4: error: unknown type name 'RDataFrame' RDataFrame f(treename, filename.Data()); ^ input_line_54:4:21: error: unknown type name 'RVecD' auto CalcPt = [](RVecD &px, RVecD &py, RVecD &E) { ^ input_line_54:4:32: error: unknown type name 'RVecD' auto CalcPt = [](RVecD &px, RVecD &py, RVecD &E) { ^ input_line_54:4:43: error: unknown type name 'RVecD' auto CalcPt = [](RVecD &px, RVecD &py, RVecD &E) { ^ input_line_54:9:14: error: use of undeclared identifier 'RVecD' .Histo1D<RVecD>({"pt", "pt", 16, 0, 4}, "good_pt")->DrawCopy(); ^
Definition of a helper function:
%%cpp -d
void WithRDataFrameVecOpsJit()
{
RDataFrame f(treename, filename.Data());
f.Define("good_pt", "sqrt(px*px + py*py)[E>100]")
.Histo1D({"pt", "pt", 16, 0, 4}, "good_pt")->DrawCopy();
}
input_line_55:3:4: error: unknown type name 'RDataFrame' RDataFrame f(treename, filename.Data()); ^
We plot four times the same quantity, the key is to look into the implementation of the functions above.
auto c = new TCanvas();
c->Divide(2,2);
c->cd(1);
WithTTreeReader();
c->cd(2);
WithRDataFrame();
c->cd(3);
WithRDataFrameVecOps();
c->cd(4);
WithRDataFrameVecOpsJit();
input_line_59:2:3: error: use of undeclared identifier 'WithRDataFrame' (WithRDataFrame()) ^ Error in <HandleInterpreterException>: Error evaluating expression (WithRDataFrame()) Execution of your code was aborted.
Draw all canvases
gROOT->GetListOfCanvases()->Draw()