Mini-Analysis on CMS OpenData with RDataFrame. This tutorial illustrates that analyzing data with RDataFrame works the same for both TTree data and RNTuple data. The RNTuple data are converted from the Events tree in http://root.cern/files/NanoAOD_DoubleMuon_CMS2011OpenData.root Based on RDataFrame's df102_NanoAODDimuonAnalysis.C
Author: The ROOT Team
This notebook tutorial was automatically generated with ROOTBOOK-izer from the macro found in the ROOT repository on Wednesday, April 17, 2024 at 11:24 AM.
%%cpp -d
#include <ROOT/RDataFrame.hxx>
#include <ROOT/RNTupleDS.hxx>
#include <ROOT/RVec.hxx>
#include <TCanvas.h>
#include <TH1D.h>
#include <TLatex.h>
#include <TStyle.h>
#include <cassert>
#include <cmath>
#include <iostream>
#include <memory>
#include <string>
#include <vector>
#include <utility>
using RNTupleDS = ROOT::Experimental::RNTupleDS;
constexpr char const* kNTupleFileName = "http://root.cern/files/tutorials/ntpl004_dimuon_v1rc2.root";
using namespace ROOT::VecOps;
Use all available CPU cores
ROOT::EnableImplicitMT();
ROOT::RDataFrame df("Events", kNTupleFileName);
Warning in <[ROOT.NTuple] Warning /home/sftnight/build/workspace/root-makedoc-master/rootspi/rdoc/src/master/tree/ntuple/v7/src/RPageStorageFile.cxx:262 in ROOT::Experimental::Internal::RPageSourceFile::InitDescriptor(const ROOT::Experimental::RNTuple&)::<lambda()>>: Pre-release format version: RC 2
The tutorial is identical to df102_NanoAODDimuonAnalysis except the use of RNTuple.
For simplicity, select only events with exactly two muons and require opposite charge
auto df_2mu = df.Filter("nMuon == 2", "Events with exactly two muons");
auto df_os = df_2mu.Filter("Muon_charge[0] != Muon_charge[1]", "Muons with opposite charge");
Compute invariant mass of the dimuon system
auto df_mass = df_os.Define("Dimuon_mass", InvariantMass<float>, {"Muon_pt", "Muon_eta", "Muon_phi", "Muon_mass"});
Make histogram of dimuon mass spectrum
auto h = df_mass.Histo1D({"Dimuon_mass", "Dimuon_mass", 30000, 0.25, 300}, "Dimuon_mass");
Request cut-flow report
auto report = df_mass.Report();
Produce plot
gStyle->SetOptStat(0); gStyle->SetTextFont(42);
auto c = new TCanvas("c", "", 800, 700);
c->SetLogx(); c->SetLogy();
h->SetTitle("");
h->GetXaxis()->SetTitle("m_{#mu#mu} (GeV)"); h->GetXaxis()->SetTitleSize(0.04);
h->GetYaxis()->SetTitle("N_{Events}"); h->GetYaxis()->SetTitleSize(0.04);
h->DrawCopy();
TLatex label; label.SetNDC(true);
label.DrawLatex(0.175, 0.740, "#eta");
label.DrawLatex(0.205, 0.775, "#rho,#omega");
label.DrawLatex(0.270, 0.740, "#phi");
label.DrawLatex(0.400, 0.800, "J/#psi");
label.DrawLatex(0.415, 0.670, "#psi'");
label.DrawLatex(0.485, 0.700, "Y(1,2,3S)");
label.DrawLatex(0.755, 0.680, "Z");
label.SetTextSize(0.040); label.DrawLatex(0.100, 0.920, "#bf{CMS Open Data}");
label.SetTextSize(0.030); label.DrawLatex(0.630, 0.920, "#sqrt{s} = 8 TeV, L_{int} = 11.6 fb^{-1}");
Print cut-flow report
report->Print();
Events with exactly two muons: pass=1973491 all=4000000 -- eff=49.34 % cumulative eff=49.34 % Muons with opposite charge: pass=1501030 all=1973491 -- eff=76.06 % cumulative eff=37.53 %
Draw all canvases
gROOT->GetListOfCanvases()->Draw()