Logo ROOT   6.14/05
Reference Guide
df014_CSVDataSource.C
Go to the documentation of this file.
1 /// \file
2 /// \ingroup tutorial_dataframe
3 /// \notebook
4 /// This tutorial illustrates how use the RDataFrame in combination with a
5 /// RDataSource. In this case we use a TCsvDS. This data source allows to read
6 /// a CSV file from a RDataFrame.
7 /// As a result of running this tutorial, we will produce plots of the dimuon
8 /// spectrum starting from a subset of the CMS collision events of Run2010B.
9 /// Dataset Reference:
10 /// McCauley, T. (2014). Dimuon event information derived from the Run2010B
11 /// public Mu dataset. CERN Open Data Portal.
12 /// DOI: [10.7483/OPENDATA.CMS.CB8H.MFFA](http://opendata.cern.ch/record/700).
13 ///
14 /// \macro_code
15 ///
16 /// \date October 2017
17 /// \author Enric Tejedor
18 
20 {
21  // Let's first create a RDF that will read from the CSV file.
22  // The types of the columns will be automatically inferred.
23  auto fileName = "df014_CsvDataSource_MuRun2010B.csv";
24  auto tdf = ROOT::RDF::MakeCsvDataFrame(fileName);
25 
26  // Now we will apply a first filter based on two columns of the CSV,
27  // and we will define a new column that will contain the invariant mass.
28  // Note how the new invariant mass column is defined from several other
29  // columns that already existed in the CSV file.
30  auto filteredEvents =
31  tdf.Filter("Q1 * Q2 == -1")
32  .Define("m", "sqrt(pow(E1 + E2, 2) - (pow(px1 + px2, 2) + pow(py1 + py2, 2) + pow(pz1 + pz2, 2)))");
33 
34  // Next we create a histogram to hold the invariant mass values and we draw it.
35  auto invMass =
36  filteredEvents.Histo1D({"invMass", "CMS Opendata: #mu#mu mass;#mu#mu mass [GeV];Events", 512, 2, 110}, "m");
37 
38  auto c = new TCanvas();
39  c->SetLogx();
40  c->SetLogy();
41  invMass->DrawClone();
42 
43  // We will now produce a plot also for the J/Psi particle. We will plot
44  // on the same canvas the full spectrum and the zoom in the J/psi particle.
45  // First we will create the full spectrum histogram from the invariant mass
46  // column, using a different histogram model than before.
47  auto fullSpectrum =
48  filteredEvents.Histo1D({"Spectrum", "Subset of CMS Run 2010B;#mu#mu mass [GeV];Events", 1024, 2, 110}, "m");
49 
50  // Next we will create the histogram for the J/psi particle, applying first
51  // the corresponding cut.
52  double jpsiLow = 2.95;
53  double jpsiHigh = 3.25;
54  auto jpsiCut = [jpsiLow, jpsiHigh](double m) { return m < jpsiHigh && m > jpsiLow; };
55  auto jpsi =
56  filteredEvents.Filter(jpsiCut, {"m"})
57  .Histo1D({"jpsi", "Subset of CMS Run 2010B: J/#psi window;#mu#mu mass [GeV];Events", 128, jpsiLow, jpsiHigh},
58  "m");
59 
60  // Finally we draw the two histograms side by side.
61  auto dualCanvas = new TCanvas("DualCanvas", "DualCanvas", 800, 512);
62  dualCanvas->Divide(2, 1);
63  auto leftPad = dualCanvas->cd(1);
64  leftPad->SetLogx();
65  leftPad->SetLogy();
66  fullSpectrum->DrawClone("Hist");
67  dualCanvas->cd(2);
68  jpsi->DrawClone("HistP");
69 
70  return 0;
71 }
RDataFrame MakeCsvDataFrame(std::string_view fileName, bool readHeaders=true, char delimiter=',', Long64_t linesChunkSize=-1LL)
Factory method to create a CSV RDataFrame.
Definition: RCsvDS.cxx:439
auto * m
Definition: textangle.C:8
The Canvas class.
Definition: TCanvas.h:31
#define c(i)
Definition: RSha256.hxx:101