Logo ROOT  
Reference Guide
df102_NanoAODDimuonAnalysis.C
Go to the documentation of this file.
1/// \file
2/// \ingroup tutorial_dataframe
3/// \notebook -js
4/// This tutorial illustrates how NanoAOD files can be processed with ROOT
5/// dataframes. The NanoAOD-like input files are filled with 66 mio. events
6/// from CMS OpenData containing muon candidates part of 2012 dataset
7/// ([DOI: 10.7483/OPENDATA.CMS.YLIC.86ZZ](http://opendata.cern.ch/record/6004)
8/// and [DOI: 10.7483/OPENDATA.CMS.M5AD.Y3V3](http://opendata.cern.ch/record/6030)).
9/// The macro matches muon pairs and produces an histogram of the dimuon mass
10/// spectrum showing resonances up to the Z mass.
11/// Note that the bump at 30 GeV is not a resonance but a trigger effect.
12///
13/// Some more details about the dataset:
14/// - It contains about 66 millions events (muon and electron collections, plus some other information, e.g. about primary vertices)
15/// - It spans two compressed ROOT files located on EOS for about a total size of 7.5 GB.
16///
17/// \macro_image
18/// \macro_code
19/// \macro_output
20///
21/// \date August 2018
22/// \author Stefan Wunsch (KIT, CERN)
23
24#include "ROOT/RDataFrame.hxx"
25#include "ROOT/RVec.hxx"
26#include "TCanvas.h"
27#include "TH1D.h"
28#include "TLatex.h"
29#include "Math/Vector4D.h"
30#include "TStyle.h"
31
32using namespace ROOT::VecOps;
33
35{
36 // Enable multi-threading
38
39 // Create dataframe from NanoAOD files
40 ROOT::RDataFrame df("Events",
41 {"root://eospublic.cern.ch//eos/root-eos/cms_opendata_2012_nanoaod/Run2012B_DoubleMuParked.root",
42 "root://eospublic.cern.ch//eos/root-eos/cms_opendata_2012_nanoaod/Run2012C_DoubleMuParked.root"});
43
44 // For simplicity, select only events with exactly two muons and require opposite charge
45 auto df_2mu = df.Filter("nMuon == 2", "Events with exactly two muons");
46 auto df_os = df_2mu.Filter("Muon_charge[0] != Muon_charge[1]", "Muons with opposite charge");
47
48 // Compute invariant mass of the dimuon system
49 auto df_mass = df_os.Define("Dimuon_mass", InvariantMass<float>, {"Muon_pt", "Muon_eta", "Muon_phi", "Muon_mass"});
50
51 // Make histogram of dimuon mass spectrum
52 auto h = df_mass.Histo1D({"Dimuon_mass", "Dimuon_mass", 30000, 0.25, 300}, "Dimuon_mass");
53
54 // Request cut-flow report
55 auto report = df_mass.Report();
56
57 // Produce plot
59 auto c = new TCanvas("c", "", 800, 700);
60 c->SetLogx(); c->SetLogy();
61
62 h->SetTitle("");
63 h->GetXaxis()->SetTitle("m_{#mu#mu} (GeV)"); h->GetXaxis()->SetTitleSize(0.04);
64 h->GetYaxis()->SetTitle("N_{Events}"); h->GetYaxis()->SetTitleSize(0.04);
65 h->DrawClone();
66
67 TLatex label; label.SetNDC(true);
68 label.DrawLatex(0.175, 0.740, "#eta");
69 label.DrawLatex(0.205, 0.775, "#rho,#omega");
70 label.DrawLatex(0.270, 0.740, "#phi");
71 label.DrawLatex(0.400, 0.800, "J/#psi");
72 label.DrawLatex(0.415, 0.670, "#psi'");
73 label.DrawLatex(0.485, 0.700, "Y(1,2,3S)");
74 label.DrawLatex(0.755, 0.680, "Z");
75 label.SetTextSize(0.040); label.DrawLatex(0.100, 0.920, "#bf{CMS Open Data}");
76 label.SetTextSize(0.030); label.DrawLatex(0.630, 0.920, "#sqrt{s} = 8 TeV, L_{int} = 11.6 fb^{-1}");
77
78 c->SaveAs("dimuon_spectrum.pdf");
79
80 // Print cut-flow report
81 report->Print();
82}
83
84int main()
85{
87}
#define c(i)
Definition: RSha256.hxx:101
#define h(i)
Definition: RSha256.hxx:106
R__EXTERN TStyle * gStyle
Definition: TStyle.h:410
RInterface< RDFDetail::RFilter< F, Proxied >, DS_t > Filter(F f, const ColumnNames_t &columns={}, std::string_view name="")
Append a filter to the call graph.
Definition: RInterface.hxx:187
ROOT's RDataFrame offers a high level interface for analyses of data stored in TTrees,...
Definition: RDataFrame.hxx:42
virtual void SetTextFont(Font_t tfont=62)
Set the text font.
Definition: TAttText.h:45
virtual void SetTextSize(Float_t tsize=1)
Set the text size.
Definition: TAttText.h:46
The Canvas class.
Definition: TCanvas.h:27
To draw Mathematical Formula.
Definition: TLatex.h:18
TLatex * DrawLatex(Double_t x, Double_t y, const char *text)
Make a copy of this object with the new parameters And copy object attributes.
Definition: TLatex.cxx:1926
void SetOptStat(Int_t stat=1)
The type of information printed in the histogram statistics box can be selected via the parameter mod...
Definition: TStyle.cxx:1590
virtual void SetNDC(Bool_t isNDC=kTRUE)
Set NDC mode on if isNDC = kTRUE, off otherwise.
Definition: TText.cxx:813
int main(int argc, char **argv)
void EnableImplicitMT(UInt_t numthreads=0)
Enable ROOT's implicit multi-threading for all objects and methods that provide an internal paralleli...
Definition: TROOT.cxx:526