Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
STrackMatcherComparator.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file STrackMatcherComparator.cc
1 // ----------------------------------------------------------------------------
2 // 'STrackMatcherComparator.cc'
3 // Derek Anderson
4 // 01.30.2024
5 //
6 // Small module to produce plots from the output of the
7 // SvtxEvaluator and FillMatchedClusTree modules.
8 // ----------------------------------------------------------------------------
9 
10 // c++ utilities
11 #include <cmath>
12 #include <array>
13 #include <string>
14 #include <vector>
15 #include <cassert>
16 #include <utility>
17 #include <iostream>
18 // root libraries
19 #include <TH1.h>
20 #include <TH2.h>
21 #include <TPad.h>
22 #include <TFile.h>
23 #include <TLine.h>
24 #include <TError.h>
25 #include <TString.h>
26 #include <TLegend.h>
27 #include <TCanvas.h>
28 #include <TPaveText.h>
29 #include <TDirectory.h>
30 // analysis utilities
31 #include "STrackMatcherComparator.h"
32 #include "STrackMatcherComparatorConfig.h"
33 #include "STrackMatcherComparatorHistDef.h"
34 #include "STrackMatcherComparatorHistContent.h"
35 
36 // make common namespaces implicit
37 using namespace std;
38 
39 
40 
41 // ctor/dtor ------------------------------------------------------------------
42 
43 STrackMatcherComparator::STrackMatcherComparator(optional<STrackMatcherComparatorConfig> config) {
44 
45  if (config.has_value()) {
46  m_config = config.value();
47  }
48 
49  // make sure vectors are empty
50  m_vecHistDirs.clear();
51  m_vecRatioDirs.clear();
52  m_vecPlotDirs.clear();
53  m_vecTreeHists1D.clear();
54  m_vecTupleHists1D.clear();
55  m_vecOldHists1D.clear();
56  m_vecTreeHists2D.clear();
57  m_vecTupleHists2D.clear();
58  m_vecOldHists2D.clear();
59 
60 } // end ctor
61 
62 
63 
64 STrackMatcherComparator::~STrackMatcherComparator() {
65 
66  /* nothing to do */
67 
68 } // end dtor
69 
70 
71 
72 // public methods -------------------------------------------------------------
73 
74 void STrackMatcherComparator::Init() {
75 
76  // announce start
77  cout << "\n Beginning track matcher comparison...\n"
78  << " Initialization:"
79  << endl;
80 
81  // run initialization routines
82  OpenOutput();
83  OpenInput();
84  InitHists();
85  return;
86 
87 } // end 'Init()'
88 
89 
90 
91 void STrackMatcherComparator::Analyze() {
92 
93  // announce analyzing
94  cout << " Analysis:" << endl;
95 
96  // run analysis routines
97  GetNewTreeHists();
98  GetNewTupleHists();
99  GetOldTupleHists();
100  return;
101 
102 } // end 'Analysis()'
103 
104 
105 
106 void STrackMatcherComparator::End() {
107 
108  // announce end
109  cout << " End:" << endl;
110 
111  // run ending routines
112  MakeRatiosAndPlots(m_vecTreeHists1D, m_vecTreeHists2D, Src::NewTree, "VsNewTree");
113  MakeRatiosAndPlots(m_vecTupleHists1D, m_vecTupleHists2D, Src::NewTuple, "VsNewTuple");
114  SaveHistograms();
115  CloseInput();
116  CloseOutput();
117 
118  // exit
119  cout << " Finished track matcher comparison!\n" << endl;
120  return;
121 }
122 
123 
124 
125 // internal methods -----------------------------------------------------------
126 
127 void STrackMatcherComparator::OpenOutput() {
128 
129  // open output files
130  m_outFile = new TFile(m_config.outFileName.data(), "recreate");
131  if (!m_outFile) {
132  cerr << "PANIC: couldn't open output file!" << endl;
133  assert(m_outFile);
134  }
135  cout << " Opened output file." << endl;
136 
137  // create output directories
138  m_vecHistDirs.resize(m_const.nDirHist);
139  m_vecRatioDirs.resize(m_const.nDirRatio);
140  m_vecPlotDirs.resize(m_const.nDirPlot);
141 
142  m_vecHistDirs[Src::NewTree] = (TDirectory*) m_outFile -> mkdir(m_config.histDirNames.at(Src::NewTree).data());
143  m_vecHistDirs[Src::NewTuple] = (TDirectory*) m_outFile -> mkdir(m_config.histDirNames.at(Src::NewTuple).data());
144  m_vecHistDirs[Src::OldTuple] = (TDirectory*) m_outFile -> mkdir(m_config.histDirNames.at(Src::OldTuple).data());
145  m_vecRatioDirs[Src::NewTree] = (TDirectory*) m_outFile -> mkdir(m_config.ratioDirNames.at(Src::NewTree).data());
146  m_vecRatioDirs[Src::NewTuple] = (TDirectory*) m_outFile -> mkdir(m_config.ratioDirNames.at(Src::NewTuple).data());
147  m_vecPlotDirs[Src::NewTree] = (TDirectory*) m_outFile -> mkdir(m_config.plotDirNames.at(Src::NewTree).data());
148  m_vecPlotDirs[Src::NewTuple] = (TDirectory*) m_outFile -> mkdir(m_config.plotDirNames.at(Src::NewTuple).data());
149 
150  // announce end and return
151  cout << " Created output directories." << endl;
152  return;
153 
154 } // end 'OpenOutput()'
155 
156 
157 
158 void STrackMatcherComparator::OpenInput() {
159 
160  // open true files
161  m_treeInFileTrue = new TFile(m_config.newInFileName.data(), "read");
162  m_tupleInFileTrue = new TFile(m_config.newInFileName.data(), "read");
163  m_oldInFileTrue = new TFile(m_config.oldInFileName.data(), "read");
164  if (!m_treeInFileTrue || !m_tupleInFileTrue || !m_oldInFileTrue) {
165  cerr << "PANIC: couldn't open an input file!\n"
166  << " m_treeInFileTrue = " << m_treeInFileTrue << "\n"
167  << " m_tupleInFileTrue = " << m_tupleInFileTrue << "\n"
168  << " m_oldInFileTrue = " << m_oldInFileTrue
169  << endl;
170  assert(m_treeInFileTrue && m_tupleInFileTrue && m_oldInFileTrue);
171  }
172  cout << " Opened truth input files." << endl;
173 
174  // open reco files
175  m_treeInFileReco = new TFile(m_config.newInFileName.data(), "read");
176  m_tupleInFileReco = new TFile(m_config.newInFileName.data(), "read");
177  m_oldInFileReco = new TFile(m_config.oldInFileName.data(), "read");
178  if (!m_treeInFileReco || !m_tupleInFileReco || !m_oldInFileReco) {
179  cerr << "PANIC: couldn't open an input file!\n"
180  << " m_treeInFileReco = " << m_treeInFileReco << "\n"
181  << " m_tupleInFileReco = " << m_tupleInFileReco << "\n"
182  << " m_oldInFileReco = " << m_oldInFileReco
183  << endl;
184  assert(m_treeInFileReco && m_tupleInFileReco && m_oldInFileReco);
185  }
186  cout << " Opened reco input files." << endl;
187 
188  // grab true input trees/tuples
189  m_tTreeTrue = (TTree*) m_treeInFileTrue -> Get(m_config.newTreeTrueName.data());
190  m_ntTupleTrue = (TNtuple*) m_tupleInFileTrue -> Get(m_config.newTupleTrueName.data());
191  m_ntOldTrue = (TNtuple*) m_oldInFileTrue -> Get(m_config.oldTupleTrueName.data());
192  if (!m_tTreeTrue || !m_ntTupleTrue || !m_ntOldTrue) {
193  cerr << "PANIC: couldn't grab an input truth tree/tuple!\n"
194  << " m_tTreeTrue = " << m_tTreeTrue << "\n"
195  << " m_ntTupleTrue = " << m_ntTupleTrue << "\n"
196  << " m_ntOldTrue = " << m_ntOldTrue
197  << endl;
198  assert(m_tTreeTrue && m_ntTupleTrue && m_ntOldTrue);
199  }
200  cout << " Grabbed input truth trees/tuples." << endl;
201 
202  // grab reco input trees/tuples
203  m_tTreeReco = (TTree*) m_treeInFileReco -> Get(m_config.newTreeRecoName.data());
204  m_ntTupleReco = (TNtuple*) m_tupleInFileReco -> Get(m_config.newTupleRecoName.data());
205  m_ntOldReco = (TNtuple*) m_oldInFileReco -> Get(m_config.oldTupleRecoName.data());
206  if (!m_tTreeReco || !m_ntTupleReco || !m_ntOldReco) {
207  cerr << "PANIC: couldn't grab an input truth tree/tuple!\n"
208  << " m_tTreeReco = " << m_tTreeReco << "\n"
209  << " m_ntTupleReco = " << m_ntTupleReco << "\n"
210  << " m_ntOldReco = " << m_ntOldReco
211  << endl;
212  assert(m_tTreeReco && m_ntTupleReco && m_ntOldReco);
213  }
214 
215  // announce end and return
216  cout << " Grabbed input truth trees/tuples." << endl;
217  return;
218 
219 } // end 'OpenInput()'
220 
221 
222 
223 void STrackMatcherComparator::InitHists() {
224 
225  // announce start of routine
226  const size_t nHistRows = m_hist.vecNameBase.size();
227  const size_t nHistTypes = m_hist.vecNameBase[0].size();
228  const size_t nVsHistMods = m_hist.vecVsModifiers.size();
229  cout << " Initializing histograms." << endl;
230 
231  // output histogram base binning definitions
232  vector<tuple<uint32_t, pair<float, float>>> vecBaseHistBins = m_hist.GetVecHistBins();
233  vector<tuple<uint32_t, pair<float, float>>> vecVsHistBins = m_hist.GetVecVsHistBins();
234 
235  // turn on errors
236  TH1::SetDefaultSumw2(true);
237  TH2::SetDefaultSumw2(true);
238 
239  // set up 1d output histograms
240  m_vecTreeHists1D.resize(nHistRows);
241  m_vecTupleHists1D.resize(nHistRows);
242  m_vecOldHists1D.resize(nHistRows);
243  for (size_t iHistRow = 0; iHistRow < nHistRows; iHistRow++) {
244  for (size_t iHistType = 0; iHistType < nHistTypes; iHistType++) {
245 
246  // create names
247  const string sTreeHist1D = m_hist.vecNameBase[iHistRow][iHistType] + "_" + m_config.newTreeLabel;
248  const string sTupleHist1D = m_hist.vecNameBase[iHistRow][iHistType] + "_" + m_config.newTupleLabel;
249  const string sOldHist1D = m_hist.vecNameBase[iHistRow][iHistType] + "_" + m_config.oldTupleLabel;
250 
251  // create output histograms
252  m_vecTreeHists1D[iHistRow].push_back(
253  new TH1D(
254  sTreeHist1D.data(),
255  "",
256  get<0>(vecBaseHistBins[iHistRow]),
257  get<1>(vecBaseHistBins[iHistRow]).first,
258  get<1>(vecBaseHistBins[iHistRow]).second
259  )
260  );
261  m_vecTupleHists1D[iHistRow].push_back(
262  new TH1D(
263  sTupleHist1D.data(),
264  "",
265  get<0>(vecBaseHistBins[iHistRow]),
266  get<1>(vecBaseHistBins[iHistRow]).first,
267  get<1>(vecBaseHistBins[iHistRow]).second
268  )
269  );
270  m_vecOldHists1D[iHistRow].push_back(
271  new TH1D(
272  sOldHist1D.data(),
273  "",
274  get<0>(vecBaseHistBins[iHistRow]),
275  get<1>(vecBaseHistBins[iHistRow]).first,
276  get<1>(vecBaseHistBins[iHistRow]).second
277  )
278  );
279 
280  // set axis titles
281  m_vecTreeHists1D[iHistRow].back() -> GetXaxis() -> SetTitle(m_hist.vecBaseAxisVars.at(iHistRow).data());
282  m_vecTupleHists1D[iHistRow].back() -> GetXaxis() -> SetTitle(m_hist.vecBaseAxisVars.at(iHistRow).data());
283  m_vecOldHists1D[iHistRow].back() -> GetXaxis() -> SetTitle(m_hist.vecBaseAxisVars.at(iHistRow).data());
284  m_vecTreeHists1D[iHistRow].back() -> GetYaxis() -> SetTitle(m_config.count.data());
285  m_vecTupleHists1D[iHistRow].back() -> GetYaxis() -> SetTitle(m_config.count.data());
286  m_vecOldHists1D[iHistRow].back() -> GetYaxis() -> SetTitle(m_config.count.data());
287  } // end type loop
288  } // end row loop
289  cout << " Initialized 1D histograms." << endl;
290 
291  // set up 2d output histograms
292  m_vecTreeHists2D.resize(nHistRows);
293  m_vecTupleHists2D.resize(nHistRows);
294  m_vecOldHists2D.resize(nHistRows);
295  for (size_t iHistRow = 0; iHistRow < nHistRows; iHistRow++) {
296 
297  // reserve space for all histograms in row
298  m_vecTreeHists2D[iHistRow].resize(nHistTypes);
299  m_vecTupleHists2D[iHistRow].resize(nHistTypes);
300  m_vecOldHists2D[iHistRow].resize(nHistTypes);
301  for (size_t iHistType = 0; iHistType < nHistTypes; iHistType++) {
302  for (size_t iVsHistMod = 0; iVsHistMod < nVsHistMods; iVsHistMod++) {
303 
304  // create names
305  const string sTreeHist2D = m_hist.vecNameBase[iHistRow][iHistType] + m_hist.vecVsModifiers[iVsHistMod] + "_" + m_config.newTreeLabel;
306  const string sTupleHist2D = m_hist.vecNameBase[iHistRow][iHistType] + m_hist.vecVsModifiers[iVsHistMod] + "_" + m_config.newTupleLabel;
307  const string sOldHist2D = m_hist.vecNameBase[iHistRow][iHistType] + m_hist.vecVsModifiers[iVsHistMod] + "_" + m_config.oldTupleLabel;
308 
309  // create output histograms
310  m_vecTreeHists2D[iHistRow][iHistType].push_back(
311  new TH2D(
312  sTreeHist2D.data(),
313  "",
314  get<0>(vecVsHistBins[iVsHistMod]),
315  get<1>(vecVsHistBins[iVsHistMod]).first,
316  get<1>(vecVsHistBins[iVsHistMod]).second,
317  get<0>(vecBaseHistBins[iHistRow]),
318  get<1>(vecBaseHistBins[iHistRow]).first,
319  get<1>(vecBaseHistBins[iHistRow]).second
320  )
321  );
322  m_vecTupleHists2D[iHistRow][iHistType].push_back(
323  new TH2D(
324  sTupleHist2D.data(),
325  "",
326  get<0>(vecVsHistBins[iVsHistMod]),
327  get<1>(vecVsHistBins[iVsHistMod]).first,
328  get<1>(vecVsHistBins[iVsHistMod]).second,
329  get<0>(vecBaseHistBins[iHistRow]),
330  get<1>(vecBaseHistBins[iHistRow]).first,
331  get<1>(vecBaseHistBins[iHistRow]).second
332  )
333  );
334  m_vecOldHists2D[iHistRow][iHistType].push_back(
335  new TH2D(
336  sOldHist2D.data(),
337  "",
338  get<0>(vecVsHistBins[iVsHistMod]),
339  get<1>(vecVsHistBins[iVsHistMod]).first,
340  get<1>(vecVsHistBins[iVsHistMod]).second,
341  get<0>(vecBaseHistBins[iHistRow]),
342  get<1>(vecBaseHistBins[iHistRow]).first,
343  get<1>(vecBaseHistBins[iHistRow]).second
344  )
345  );
346 
347  // set axis titles
348  m_vecTreeHists2D[iHistRow][iHistType].back() -> GetXaxis() -> SetTitle(m_hist.vecVsAxisVars.at(iVsHistMod).data());
349  m_vecTupleHists2D[iHistRow][iHistType].back() -> GetXaxis() -> SetTitle(m_hist.vecVsAxisVars.at(iVsHistMod).data());
350  m_vecOldHists2D[iHistRow][iHistType].back() -> GetXaxis() -> SetTitle(m_hist.vecVsAxisVars.at(iVsHistMod).data());
351  m_vecTreeHists2D[iHistRow][iHistType].back() -> GetYaxis() -> SetTitle(m_hist.vecBaseAxisVars.at(iHistRow).data());
352  m_vecTupleHists2D[iHistRow][iHistType].back() -> GetYaxis() -> SetTitle(m_hist.vecBaseAxisVars.at(iHistRow).data());
353  m_vecOldHists2D[iHistRow][iHistType].back() -> GetYaxis() -> SetTitle(m_hist.vecBaseAxisVars.at(iHistRow).data());
354  m_vecTreeHists2D[iHistRow][iHistType].back() -> GetZaxis() -> SetTitle(m_config.count.data());
355  m_vecTupleHists2D[iHistRow][iHistType].back() -> GetZaxis() -> SetTitle(m_config.count.data());
356  m_vecOldHists2D[iHistRow][iHistType].back() -> GetZaxis() -> SetTitle(m_config.count.data());
357  } // end modifier loop
358  } // end type loop
359  } // end row loop
360 
361  // announce end and exit
362  cout << " Initialized 2D output histograms." << endl;
363  return;
364 
365 } // end 'InitHists()'
366 
367 
368 
369 void STrackMatcherComparator::GetNewTreeHists() {
370 
371  // announce start of routine
372  cout << " Grabbing new matcher tree histograms:" << endl;
373 
374  // declare input leaves (ignore cluster branches for now)
375  int tru_event;
376  int tru_nphg4_part;
377  float tru_centrality;
378  int tru_ntrackmatches;
379  int tru_nphg4;
380  int tru_nsvtx;
381  int tru_trackid;
382  bool tru_is_G4track;
383  bool tru_is_Svtrack;
384  bool tru_is_matched;
385  float tru_trkpt;
386  float tru_trketa;
387  float tru_trkphi;
388  int tru_nclus;
389  int tru_nclustpc;
390  int tru_nclusmvtx;
391  int tru_nclusintt;
392  float tru_matchrat;
393  float tru_matchrat_intt;
394  float tru_matchrat_mvtx;
395  float tru_matchrat_tpc;
396 
397  int rec_event;
398  int rec_nphg4_part;
399  float rec_centrality;
400  int rec_ntrackmatches;
401  int rec_nphg4;
402  int rec_nsvtx;
403  int rec_trackid;
404  bool rec_is_G4track;
405  bool rec_is_Svtrack;
406  bool rec_is_matched;
407  float rec_trkpt;
408  float rec_trketa;
409  float rec_trkphi;
410  int rec_nclus;
411  int rec_nclustpc;
412  int rec_nclusmvtx;
413  int rec_nclusintt;
414  float rec_matchrat;
415  float rec_matchrat_intt;
416  float rec_matchrat_mvtx;
417  float rec_matchrat_tpc;
418 
419  // set branch addresses (ignore cluster branches for now)
420  m_tTreeTrue -> SetBranchAddress("event", &tru_event);
421  m_tTreeTrue -> SetBranchAddress("nphg4_part", &tru_nphg4_part);
422  m_tTreeTrue -> SetBranchAddress("centrality", &tru_centrality);
423  m_tTreeTrue -> SetBranchAddress("ntrackmatches", &tru_ntrackmatches);
424  m_tTreeTrue -> SetBranchAddress("nphg4", &tru_nphg4);
425  m_tTreeTrue -> SetBranchAddress("nsvtx", &tru_nsvtx);
426  m_tTreeTrue -> SetBranchAddress("trackid", &tru_trackid);
427  m_tTreeTrue -> SetBranchAddress("is_G4track", &tru_is_G4track);
428  m_tTreeTrue -> SetBranchAddress("is_Svtrack", &tru_is_Svtrack);
429  m_tTreeTrue -> SetBranchAddress("is_matched", &tru_is_matched);
430  m_tTreeTrue -> SetBranchAddress("trkpt", &tru_trkpt);
431  m_tTreeTrue -> SetBranchAddress("trketa", &tru_trketa);
432  m_tTreeTrue -> SetBranchAddress("trkphi", &tru_trkphi);
433  m_tTreeTrue -> SetBranchAddress("nclus", &tru_nclus);
434  m_tTreeTrue -> SetBranchAddress("nclustpc", &tru_nclustpc);
435  m_tTreeTrue -> SetBranchAddress("nclusmvtx", &tru_nclusmvtx);
436  m_tTreeTrue -> SetBranchAddress("nclusintt", &tru_nclusintt);
437  m_tTreeTrue -> SetBranchAddress("matchrat", &tru_matchrat);
438  m_tTreeTrue -> SetBranchAddress("matchrat_intt", &tru_matchrat_intt);
439  m_tTreeTrue -> SetBranchAddress("matchrat_mvtx", &tru_matchrat_mvtx);
440  m_tTreeTrue -> SetBranchAddress("matchrat_tpc", &tru_matchrat_tpc);
441 
442  m_tTreeReco -> SetBranchAddress("event", &rec_event);
443  m_tTreeReco -> SetBranchAddress("nphg4_part", &rec_nphg4_part);
444  m_tTreeReco -> SetBranchAddress("centrality", &rec_centrality);
445  m_tTreeReco -> SetBranchAddress("ntrackmatches", &rec_ntrackmatches);
446  m_tTreeReco -> SetBranchAddress("nphg4", &rec_nphg4);
447  m_tTreeReco -> SetBranchAddress("nsvtx", &rec_nsvtx);
448  m_tTreeReco -> SetBranchAddress("trackid", &rec_trackid);
449  m_tTreeReco -> SetBranchAddress("is_G4track", &rec_is_G4track);
450  m_tTreeReco -> SetBranchAddress("is_Svtrack", &rec_is_Svtrack);
451  m_tTreeReco -> SetBranchAddress("is_matched", &rec_is_matched);
452  m_tTreeReco -> SetBranchAddress("trkpt", &rec_trkpt);
453  m_tTreeReco -> SetBranchAddress("trketa", &rec_trketa);
454  m_tTreeReco -> SetBranchAddress("trkphi", &rec_trkphi);
455  m_tTreeReco -> SetBranchAddress("nclus", &rec_nclus);
456  m_tTreeReco -> SetBranchAddress("nclustpc", &rec_nclustpc);
457  m_tTreeReco -> SetBranchAddress("nclusmvtx", &rec_nclusmvtx);
458  m_tTreeReco -> SetBranchAddress("nclusintt", &rec_nclusintt);
459  m_tTreeReco -> SetBranchAddress("matchrat", &rec_matchrat);
460  m_tTreeReco -> SetBranchAddress("matchrat_intt", &rec_matchrat_intt);
461  m_tTreeReco -> SetBranchAddress("matchrat_mvtx", &rec_matchrat_mvtx);
462  m_tTreeReco -> SetBranchAddress("matchrat_tpc", &rec_matchrat_tpc);
463  cout << " Set input branches." << endl;
464 
465  // grab no. of entries
466  const int64_t nTrueEntries = m_tTreeTrue -> GetEntries();
467  const int64_t nRecoEntries = m_tTreeReco -> GetEntries();
468  cout << " Beginning truth particle loop: " << nTrueEntries << " to process" << endl;
469 
470  // loop over truth particles
471  int64_t nTrueBytes = 0;
472  for (int64_t iTrueEntry = 0; iTrueEntry < nTrueEntries; iTrueEntry++) {
473 
474  // grab truth particle entry
475  const int64_t trueBytes = m_tTreeTrue -> GetEntry(iTrueEntry);
476  if (trueBytes < 0) {
477  cerr << "PANIC: issue with entry " << iTrueEntry << "! Aborting loop!\n" << endl;
478  break;
479  } else {
480  nTrueBytes += trueBytes;
481  }
482 
483  const int64_t iTrueProg = iTrueEntry + 1;
484  if (iTrueProg == nTrueEntries) {
485  cout << " Processing entry " << iTrueProg << "/" << nTrueEntries << "..." << endl;
486  } else {
487  cout << " Processing entry " << iTrueProg << "/" << nTrueEntries << "...\r" << flush;
488  }
489 
490  // select truth particles
491  // - FIXME add cuts
492  if (!tru_is_G4track) continue;
493 
494  // bundle histogram values to fill
495  const STrackMatcherComparatorHistContent content {
496  (double) tru_nclus,
497  (double) tru_nclusintt,
498  (double) tru_nclusmvtx,
499  (double) tru_nclustpc,
500  1.,
501  1.,
502  1.,
503  1.,
504  tru_trkphi,
505  tru_trketa,
506  tru_trkpt,
507  1.,
508  1.,
509  1.,
510  1.,
511  1.,
512  1.,
513  1.,
514  1.
515  };
516 
517  // fill truth histograms
518  FillHistogram1D(content, Type::Truth, m_vecTreeHists1D);
519  FillHistogram2D(content, Type::Truth, Comp::VsTruthPt, tru_trkpt, m_vecTreeHists2D);
520  FillHistogram2D(content, Type::Truth, Comp::VsNumTpc, tru_nclustpc, m_vecTreeHists2D);
521  } // end truth particle loop
522 
523  // announce next entry loop
524  cout << " Finished truth particle loop.\n"
525  << " Beginning reconstructed track loop: " << nRecoEntries << " to process"
526  << endl;
527 
528  // loop over reco tracks
529  // - TODO identify matched truth particles
530  int64_t nRecoBytes = 0;
531  for (int64_t iRecoEntry = 0; iRecoEntry < nRecoEntries; iRecoEntry++) {
532 
533  // grab reco track entry
534  const int64_t recoBytes = m_tTreeReco -> GetEntry(iRecoEntry);
535  if (recoBytes < 0) {
536  cerr << "PANIC: issue with entry " << iRecoEntry << "! Aborting loop!\n" << endl;
537  break;
538  } else {
539  nRecoBytes += recoBytes;
540  }
541 
542  const int64_t iRecoProg = iRecoEntry + 1;
543  if (iRecoProg == nRecoEntries) {
544  cout << " Processing entry " << iRecoProg << "/" << nRecoEntries << "..." << endl;
545  } else {
546  cout << " Processing entry " << iRecoProg << "/" << nRecoEntries << "...\r" << flush;
547  }
548 
549  // select only tracks matched to truth particle
550  // - FIXME add cuts
551  if (!rec_is_matched || !rec_is_Svtrack) continue;
552 
553  // bundle histogram values to fill
554  // - FIXME actually calculate pt fraction
555  // - FIXME include quality, errors, and resolutions
556  const STrackMatcherComparatorHistContent content {
557  (double) rec_nclus,
558  (double) rec_nclusintt,
559  (double) rec_nclusmvtx,
560  (double) rec_nclustpc,
561  rec_matchrat,
562  rec_matchrat_intt,
563  rec_matchrat_mvtx,
564  rec_matchrat_tpc,
565  rec_trkphi,
566  rec_trketa,
567  rec_trkpt,
568  1.,
569  1.,
570  1.,
571  1.,
572  1.,
573  1.,
574  1.,
575  1.
576  };
577 
578  // fill all matched reco histograms
579  // - FIXME use actual truth pt & ntpc of matched particle
580  FillHistogram1D(content, Type::Track, m_vecTreeHists1D);
581  FillHistogram2D(content, Type::Track, Comp::VsTruthPt, rec_trkpt, m_vecTreeHists2D);
582  FillHistogram2D(content, Type::Track, Comp::VsTruthPt, rec_nclustpc, m_vecTreeHists2D);
583 
584  // fill weird and normal matched reco 1D histograms
585  // - FIXME actually check if is a weird track
586  const bool isWeirdTrack = true;
587  if (isWeirdTrack) {
588  FillHistogram1D(content, Type::Weird, m_vecTreeHists1D);
589  FillHistogram2D(content, Type::Weird, Comp::VsTruthPt, rec_trkpt, m_vecTreeHists2D);
590  FillHistogram2D(content, Type::Weird, Comp::VsNumTpc, rec_nclustpc, m_vecTreeHists2D);
591  } else {
592  FillHistogram1D(content, Type::Normal, m_vecTreeHists1D);
593  FillHistogram2D(content, Type::Normal, Comp::VsTruthPt, rec_trkpt, m_vecTreeHists2D);
594  FillHistogram2D(content, Type::Normal, Comp::VsNumTpc, rec_nclustpc, m_vecTreeHists2D);
595  }
596  } // end reco track loop
597 
598  // announce method end and return
599  cout << " Finished reconstructed track loop.\n"
600  << " Finished getting histograms from new matcher cluster tree."
601  << endl;
602  return;
603 
604 } // end 'GetNewTreeHists()'
605 
606 
607 
608 void STrackMatcherComparator::GetNewTupleHists() {
609 
610  // announce start of method
611  cout << " Grabbing new matcher tuple histograms:" << endl;
612 
613  // declare input leaves
614  float tru_evtid;
615  float tru_trkid;
616  float tru_pt;
617  float tru_eta;
618  float tru_phi;
619  float tru_nmvtxclust_trkmatcher;
620  float tru_ninttclust_trkmatcher;
621  float tru_ntpclust_trkmatcher;
622  float tru_nmvtxclust_manual;
623  float tru_ninttclust_manual;
624  float tru_ntpcclust_manual;
625  float tru_nmvtxlayer_;
626  float tru_ninttlayer;
627  float tru_ntpclayer;
628  float tru_deltapt;
629  float tru_quality;
630  float tru_dcaxy;
631  float tru_dcaz;
632  float tru_sigmadcaxy;
633  float tru_sigmadacz;
634  float tru_vx;
635  float tru_vy;
636  float tru_vz;
637  float tru_gtrkid;
638  float tru_gpt;
639  float tru_geta;
640  float tru_gphi;
641  float tru_gnmvtxclust_trkmatcher;
642  float tru_gninttclust_trkmatcher;
643  float tru_gntpclust_trkmatcher;
644  float tru_gnmvtxclust_manual;
645  float tru_gninttclust_manual;
646  float tru_gntpcclust_manual;
647  float tru_gnmvtxlayer;
648  float tru_gninttlayer;
649  float tru_gntpclayer;
650  float tru_gdeltapt;
651  float tru_gquality;
652  float tru_gdcaxy;
653  float tru_gdcaz;
654  float tru_gsigmadcaxy;
655  float tru_gsigmadacz;
656  float tru_gvx;
657  float tru_gvy;
658  float tru_gvz;
659  float tru_fracnmvtxmatched;
660  float tru_fracninttmatched;
661  float tru_fracntpcmatched;
662 
663  float rec_evtid;
664  float rec_trkid;
665  float rec_pt;
666  float rec_eta;
667  float rec_phi;
668  float rec_nmvtxclust_trkmatcher;
669  float rec_ninttclust_trkmatcher;
670  float rec_ntpclust_trkmatcher;
671  float rec_nmvtxclust_manual;
672  float rec_ninttclust_manual;
673  float rec_ntpcclust_manual;
674  float rec_nmvtxlayer_;
675  float rec_ninttlayer;
676  float rec_ntpclayer;
677  float rec_deltapt;
678  float rec_quality;
679  float rec_dcaxy;
680  float rec_dcaz;
681  float rec_sigmadcaxy;
682  float rec_sigmadacz;
683  float rec_vx;
684  float rec_vy;
685  float rec_vz;
686  float rec_gtrkid;
687  float rec_gpt;
688  float rec_geta;
689  float rec_gphi;
690  float rec_gnmvtxclust_trkmatcher;
691  float rec_gninttclust_trkmatcher;
692  float rec_gntpclust_trkmatcher;
693  float rec_gnmvtxclust_manual;
694  float rec_gninttclust_manual;
695  float rec_gntpcclust_manual;
696  float rec_gnmvtxlayer;
697  float rec_gninttlayer;
698  float rec_gntpclayer;
699  float rec_gdeltapt;
700  float rec_gquality;
701  float rec_gdcaxy;
702  float rec_gdcaz;
703  float rec_gsigmadcaxy;
704  float rec_gsigmadacz;
705  float rec_gvx;
706  float rec_gvy;
707  float rec_gvz;
708  float rec_fracnmvtxmatched;
709  float rec_fracninttmatched;
710  float rec_fracntpcmatched;
711 
712  // set branch addresses
713  m_ntTupleTrue -> SetBranchAddress("evtid", &tru_evtid);
714  m_ntTupleTrue -> SetBranchAddress("trkid", &tru_trkid);
715  m_ntTupleTrue -> SetBranchAddress("pt", &tru_pt);
716  m_ntTupleTrue -> SetBranchAddress("eta", &tru_eta);
717  m_ntTupleTrue -> SetBranchAddress("phi", &tru_phi);
718  m_ntTupleTrue -> SetBranchAddress("nmvtxclust_trkmatcher", &tru_nmvtxclust_trkmatcher);
719  m_ntTupleTrue -> SetBranchAddress("ninttclust_trkmatcher", &tru_ninttclust_trkmatcher);
720  m_ntTupleTrue -> SetBranchAddress("ntpclust_trkmatcher", &tru_ntpclust_trkmatcher);
721  m_ntTupleTrue -> SetBranchAddress("nmvtxclust_manual", &tru_nmvtxclust_manual);
722  m_ntTupleTrue -> SetBranchAddress("ninttclust_manual", &tru_ninttclust_manual);
723  m_ntTupleTrue -> SetBranchAddress("ntpcclust_manual", &tru_ntpcclust_manual);
724  m_ntTupleTrue -> SetBranchAddress("nmvtxlayer_", &tru_nmvtxlayer_);
725  m_ntTupleTrue -> SetBranchAddress("ninttlayer", &tru_ninttlayer);
726  m_ntTupleTrue -> SetBranchAddress("ntpclayer", &tru_ntpclayer);
727  m_ntTupleTrue -> SetBranchAddress("deltapt", &tru_deltapt);
728  m_ntTupleTrue -> SetBranchAddress("quality", &tru_quality);
729  m_ntTupleTrue -> SetBranchAddress("dcaxy", &tru_dcaxy);
730  m_ntTupleTrue -> SetBranchAddress("dcaz", &tru_dcaz);
731  m_ntTupleTrue -> SetBranchAddress("sigmadcaxy", &tru_sigmadcaxy);
732  m_ntTupleTrue -> SetBranchAddress("sigmadacz", &tru_sigmadacz);
733  m_ntTupleTrue -> SetBranchAddress("vx", &tru_vx);
734  m_ntTupleTrue -> SetBranchAddress("vy", &tru_vy);
735  m_ntTupleTrue -> SetBranchAddress("vz", &tru_vz);
736  m_ntTupleTrue -> SetBranchAddress("gtrkid", &tru_gtrkid);
737  m_ntTupleTrue -> SetBranchAddress("gpt", &tru_gpt);
738  m_ntTupleTrue -> SetBranchAddress("geta", &tru_geta);
739  m_ntTupleTrue -> SetBranchAddress("gphi", &tru_gphi);
740  m_ntTupleTrue -> SetBranchAddress("gnmvtxclust_trkmatcher",&tru_gnmvtxclust_trkmatcher);
741  m_ntTupleTrue -> SetBranchAddress("gninttclust_trkmatcher",&tru_gninttclust_trkmatcher);
742  m_ntTupleTrue -> SetBranchAddress("gntpclust_trkmatcher", &tru_gntpclust_trkmatcher);
743  m_ntTupleTrue -> SetBranchAddress("gnmvtxclust_manual", &tru_gnmvtxclust_manual);
744  m_ntTupleTrue -> SetBranchAddress("gninttclust_manual", &tru_gninttclust_manual);
745  m_ntTupleTrue -> SetBranchAddress("gntpcclust_manual", &tru_gntpcclust_manual);
746  m_ntTupleTrue -> SetBranchAddress("gnmvtxlayer", &tru_gnmvtxlayer);
747  m_ntTupleTrue -> SetBranchAddress("gninttlayer", &tru_gninttlayer);
748  m_ntTupleTrue -> SetBranchAddress("gntpclayer", &tru_gntpclayer);
749  m_ntTupleTrue -> SetBranchAddress("gdeltapt", &tru_gdeltapt);
750  m_ntTupleTrue -> SetBranchAddress("gquality", &tru_gquality);
751  m_ntTupleTrue -> SetBranchAddress("gdcaxy", &tru_gdcaxy);
752  m_ntTupleTrue -> SetBranchAddress("gdcaz", &tru_gdcaz);
753  m_ntTupleTrue -> SetBranchAddress("gsigmadcaxy", &tru_gsigmadcaxy);
754  m_ntTupleTrue -> SetBranchAddress("gsigmadacz", &tru_gsigmadacz);
755  m_ntTupleTrue -> SetBranchAddress("gvx", &tru_gvx);
756  m_ntTupleTrue -> SetBranchAddress("gvy", &tru_gvy);
757  m_ntTupleTrue -> SetBranchAddress("gvz", &tru_gvz);
758  m_ntTupleTrue -> SetBranchAddress("fracnmvtxmatched", &tru_fracnmvtxmatched);
759  m_ntTupleTrue -> SetBranchAddress("fracninttmatched", &tru_fracninttmatched);
760  m_ntTupleTrue -> SetBranchAddress("fracntpcmatched", &tru_fracntpcmatched);
761 
762  m_ntTupleReco -> SetBranchAddress("evtid", &rec_evtid);
763  m_ntTupleReco -> SetBranchAddress("trkid", &rec_trkid);
764  m_ntTupleReco -> SetBranchAddress("pt", &rec_pt);
765  m_ntTupleReco -> SetBranchAddress("eta", &rec_eta);
766  m_ntTupleReco -> SetBranchAddress("phi", &rec_phi);
767  m_ntTupleReco -> SetBranchAddress("nmvtxclust_trkmatcher", &rec_nmvtxclust_trkmatcher);
768  m_ntTupleReco -> SetBranchAddress("ninttclust_trkmatcher", &rec_ninttclust_trkmatcher);
769  m_ntTupleReco -> SetBranchAddress("ntpclust_trkmatcher", &rec_ntpclust_trkmatcher);
770  m_ntTupleReco -> SetBranchAddress("nmvtxclust_manual", &rec_nmvtxclust_manual);
771  m_ntTupleReco -> SetBranchAddress("ninttclust_manual", &rec_ninttclust_manual);
772  m_ntTupleReco -> SetBranchAddress("ntpcclust_manual", &rec_ntpcclust_manual);
773  m_ntTupleReco -> SetBranchAddress("nmvtxlayer_", &rec_nmvtxlayer_);
774  m_ntTupleReco -> SetBranchAddress("ninttlayer", &rec_ninttlayer);
775  m_ntTupleReco -> SetBranchAddress("ntpclayer", &rec_ntpclayer);
776  m_ntTupleReco -> SetBranchAddress("deltapt", &rec_deltapt);
777  m_ntTupleReco -> SetBranchAddress("quality", &rec_quality);
778  m_ntTupleReco -> SetBranchAddress("dcaxy", &rec_dcaxy);
779  m_ntTupleReco -> SetBranchAddress("dcaz", &rec_dcaz);
780  m_ntTupleReco -> SetBranchAddress("sigmadcaxy", &rec_sigmadcaxy);
781  m_ntTupleReco -> SetBranchAddress("sigmadacz", &rec_sigmadacz);
782  m_ntTupleReco -> SetBranchAddress("vx", &rec_vx);
783  m_ntTupleReco -> SetBranchAddress("vy", &rec_vy);
784  m_ntTupleReco -> SetBranchAddress("vz", &rec_vz);
785  m_ntTupleReco -> SetBranchAddress("gtrkid", &rec_gtrkid);
786  m_ntTupleReco -> SetBranchAddress("gpt", &rec_gpt);
787  m_ntTupleReco -> SetBranchAddress("geta", &rec_geta);
788  m_ntTupleReco -> SetBranchAddress("gphi", &rec_gphi);
789  m_ntTupleReco -> SetBranchAddress("gnmvtxclust_trkmatcher", &rec_gnmvtxclust_trkmatcher);
790  m_ntTupleReco -> SetBranchAddress("gninttclust_trkmatcher", &rec_gninttclust_trkmatcher);
791  m_ntTupleReco -> SetBranchAddress("gntpclust_trkmatcher", &rec_gntpclust_trkmatcher);
792  m_ntTupleReco -> SetBranchAddress("gnmvtxclust_manual", &rec_gnmvtxclust_manual);
793  m_ntTupleReco -> SetBranchAddress("gninttclust_manual", &rec_gninttclust_manual);
794  m_ntTupleReco -> SetBranchAddress("gntpcclust_manual", &rec_gntpcclust_manual);
795  m_ntTupleReco -> SetBranchAddress("gnmvtxlayer", &rec_gnmvtxlayer);
796  m_ntTupleReco -> SetBranchAddress("gninttlayer", &rec_gninttlayer);
797  m_ntTupleReco -> SetBranchAddress("gntpclayer", &rec_gntpclayer);
798  m_ntTupleReco -> SetBranchAddress("gdeltapt", &rec_gdeltapt);
799  m_ntTupleReco -> SetBranchAddress("gquality", &rec_gquality);
800  m_ntTupleReco -> SetBranchAddress("gdcaxy", &rec_gdcaxy);
801  m_ntTupleReco -> SetBranchAddress("gdcaz", &rec_gdcaz);
802  m_ntTupleReco -> SetBranchAddress("gsigmadcaxy", &rec_gsigmadcaxy);
803  m_ntTupleReco -> SetBranchAddress("gsigmadacz", &rec_gsigmadacz);
804  m_ntTupleReco -> SetBranchAddress("gvx", &rec_gvx);
805  m_ntTupleReco -> SetBranchAddress("gvy", &rec_gvy);
806  m_ntTupleReco -> SetBranchAddress("gvz", &rec_gvz);
807  m_ntTupleReco -> SetBranchAddress("fracnmvtxmatched", &rec_fracnmvtxmatched);
808  m_ntTupleReco -> SetBranchAddress("fracninttmatched", &rec_fracninttmatched);
809  m_ntTupleReco -> SetBranchAddress("fracntpcmatched", &rec_fracntpcmatched);
810  cout << " Set input branches." << endl;
811 
812  // grab no. of entries
813  const int64_t nTrueEntries = m_ntTupleTrue -> GetEntries();
814  const int64_t nRecoEntries = m_ntTupleReco -> GetEntries();
815  cout << " Beginning truth particle loop: " << nTrueEntries << " to process" << endl;
816 
817  // loop over truth particles
818  int64_t nTrueBytes = 0;
819  for (int64_t iTrueEntry = 0; iTrueEntry < nTrueEntries; iTrueEntry++) {
820 
821  // grab truth particle entry
822  const int64_t trueBytes = m_ntTupleTrue -> GetEntry(iTrueEntry);
823  if (trueBytes < 0) {
824  cerr << "PANIC: issue with entry " << iTrueEntry << "! Aborting loop!\n" << endl;
825  break;
826  } else {
827  nTrueBytes += trueBytes;
828  }
829 
830  const int64_t iTrueProg = iTrueEntry + 1;
831  if (iTrueProg == nTrueEntries) {
832  cout << " Processing entry " << iTrueProg << "/" << nTrueEntries << "..." << endl;
833  } else {
834  cout << " Processing entry " << iTrueProg << "/" << nTrueEntries << "...\r" << flush;
835  }
836 
837  // check if near sector
838  bool isNearSector = false;
839  if (m_config.doPhiCut) {
840  isNearSector = IsNearSector(tru_gphi);
841  }
842 
843  // apply cuts
844  const bool isInZVtxCut = ((tru_gvz >= m_config.zVtxRange.first) && (tru_gvz <= m_config.zVtxRange.second));
845  if (m_config.doZVtxCut && !isInZVtxCut) continue;
846  if (m_config.doPhiCut && isNearSector) continue;
847 
848  // run calculations
849  const float tru_gnclust = tru_gnmvtxclust_trkmatcher + tru_gninttclust_trkmatcher + tru_gntpclust_trkmatcher;
850 
851  // bundle histogram values to fill
852  STrackMatcherComparatorHistContent content {
853  tru_gnclust,
854  tru_gninttclust_trkmatcher,
855  tru_gnmvtxclust_trkmatcher,
856  tru_gntpclust_trkmatcher,
857  1.,
858  1.,
859  1.,
860  1.,
861  tru_gphi,
862  tru_geta,
863  tru_gpt,
864  1.,
865  1.,
866  1.,
867  1.,
868  1.,
869  1.,
870  1.,
871  1.
872  };
873 
874  // fill truth histograms
875  FillHistogram1D(content, Type::Truth, m_vecTupleHists1D);
876  FillHistogram2D(content, Type::Truth, Comp::VsTruthPt, tru_gpt, m_vecTupleHists2D);
877  FillHistogram2D(content, Type::Truth, Comp::VsNumTpc, tru_ntpclust_trkmatcher, m_vecTupleHists2D);
878  } // end truth particle loop
879 
880  // announce next entry loop
881  cout << " Finished truth particle loop.\n"
882  << " Beginning reconstructed track loop: " << nRecoEntries << " to process"
883  << endl;
884 
885  // loop over reco tracks
886  int64_t nRecoBytes = 0;
887  for (int64_t iRecoEntry = 0; iRecoEntry < nRecoEntries; iRecoEntry++) {
888 
889  // grab reco track entry
890  const int64_t recoBytes = m_ntTupleReco -> GetEntry(iRecoEntry);
891  if (recoBytes < 0) {
892  cerr << "PANIC: issue with entry " << iRecoEntry << "! Aborting loop!\n" << endl;
893  break;
894  } else {
895  nRecoBytes += recoBytes;
896  }
897 
898  const int64_t iRecoProg = iRecoEntry + 1;
899  if (iRecoProg == nRecoEntries) {
900  cout << " Processing entry " << iRecoProg << "/" << nRecoEntries << "..." << endl;
901  } else {
902  cout << " Processing entry " << iRecoProg << "/" << nRecoEntries << "...\r" << flush;
903  }
904 
905  // check if near sector
906  bool isNearSector = false;
907  if (m_config.doPhiCut) {
908  isNearSector = IsNearSector(rec_phi);
909  }
910 
911  // apply cuts
912  const bool isInZVtxCut = ((rec_vz >= m_config.zVtxRange.first) && (rec_vz <= m_config.zVtxRange.second));
913  if (m_config.doZVtxCut && !isInZVtxCut) continue;
914  if (m_config.doPhiCut && isNearSector) continue;
915 
916  // run calculations
917  // - FIXME add other errors
918  const double rec_nclus = rec_ninttclust_trkmatcher + rec_nmvtxclust_trkmatcher + rec_ntpclust_trkmatcher;
919  const double rec_gnclus = rec_gninttclust_trkmatcher + rec_gnmvtxclust_trkmatcher + rec_gntpcclust_manual;
920  const double rec_rnclus = rec_nclus / rec_gnclus;
921  const double rec_rintt = rec_ninttclust_trkmatcher / rec_gninttclust_trkmatcher;
922  const double rec_rmaps = rec_nmvtxclust_trkmatcher / rec_gnmvtxclust_trkmatcher;
923  const double rec_rtpc = rec_ntpclust_trkmatcher / rec_gntpcclust_manual;
924  const double rec_ptfrac = rec_pt / rec_gpt;
925  const double rec_ptper = rec_deltapt / rec_pt;
926  const double rec_etaper = 1.;
927  const double rec_phiper = 1.;
928  const double rec_ptres = abs(rec_pt - rec_gpt) / rec_gpt;
929  const double rec_etares = abs(rec_eta - rec_geta) / rec_geta;
930  const double rec_phires = abs(rec_phi - rec_gphi) / rec_gphi;
931 
932  // bundle histogram values to fill
933  const STrackMatcherComparatorHistContent content {
934  rec_nclus,
935  rec_ninttclust_trkmatcher,
936  rec_nmvtxclust_trkmatcher,
937  rec_ntpclust_trkmatcher,
938  rec_rnclus,
939  rec_rintt,
940  rec_rmaps,
941  rec_rtpc,
942  rec_phi,
943  rec_eta,
944  rec_pt,
945  rec_ptfrac,
946  rec_quality,
947  rec_ptper,
948  rec_etaper,
949  rec_phiper,
950  rec_ptres,
951  rec_etares,
952  rec_phires
953  };
954 
955  // fill all matched reco histograms
956  FillHistogram1D(content, Type::Track, m_vecTupleHists1D);
957  FillHistogram2D(content, Type::Track, Comp::VsTruthPt, rec_gpt, m_vecTupleHists2D);
958  FillHistogram2D(content, Type::Track, Comp::VsNumTpc, rec_ntpclust_trkmatcher, m_vecTupleHists2D);
959 
960  // fill weird and normal matched reco 1D histograms
961  const bool isNormalTrack = ((rec_ptfrac >= m_config.oddPtFrac.first) && (rec_ptfrac <= m_config.oddPtFrac.second));
962  if (isNormalTrack) {
963  FillHistogram1D(content, Type::Normal, m_vecTupleHists1D);
964  FillHistogram2D(content, Type::Normal, Comp::VsTruthPt, rec_gpt, m_vecTupleHists2D);
965  FillHistogram2D(content, Type::Normal, Comp::VsNumTpc, rec_ntpclust_trkmatcher, m_vecTupleHists2D);
966  } else {
967  FillHistogram1D(content, Type::Weird, m_vecTupleHists1D);
968  FillHistogram2D(content, Type::Weird, Comp::VsTruthPt, rec_gpt, m_vecTupleHists2D);
969  FillHistogram2D(content, Type::Weird, Comp::VsNumTpc, rec_ntpclust_trkmatcher, m_vecTupleHists2D);
970  }
971  } // end reco track loop
972 
973  // announce method end and return
974  cout << " Finished reconstructed track loop.\n"
975  << " Finished getting histograms from new matcher track tuple."
976  << endl;
977  return;
978 
979 } // end 'GetNewTupleHists()'
980 
981 
982 
983 void STrackMatcherComparator::GetOldTupleHists() {
984 
985  // announce start of method
986  cout << " Grabbing old matcher tuple histograms:" << endl;
987 
988  // declare input leaves
989  float tru_event;
990  float tru_seed;
991  float tru_gntracks;
992  float tru_gnchghad;
993  float tru_gtrackID;
994  float tru_gflavor;
995  float tru_gnhits;
996  float tru_gnmaps;
997  float tru_gnintt;
998  float tru_gnmms;
999  float tru_gnintt1;
1000  float tru_gnintt2;
1001  float tru_gnintt3;
1002  float tru_gnintt4;
1003  float tru_gnintt5;
1004  float tru_gnintt6;
1005  float tru_gnintt7;
1006  float tru_gnintt8;
1007  float tru_gntpc;
1008  float tru_gnlmaps;
1009  float tru_gnlintt;
1010  float tru_gnltpc;
1011  float tru_gnlmms;
1012  float tru_gpx;
1013  float tru_gpy;
1014  float tru_gpz;
1015  float tru_gpt;
1016  float tru_geta;
1017  float tru_gphi;
1018  float tru_gvx;
1019  float tru_gvy;
1020  float tru_gvz;
1021  float tru_gvt;
1022  float tru_gfpx;
1023  float tru_gfpy;
1024  float tru_gfpz;
1025  float tru_gfx;
1026  float tru_gfy;
1027  float tru_gfz;
1028  float tru_gembed;
1029  float tru_gprimary;
1030  float tru_trackID;
1031  float tru_px;
1032  float tru_py;
1033  float tru_pz;
1034  float tru_pt;
1035  float tru_eta;
1036  float tru_phi;
1037  float tru_deltapt;
1038  float tru_deltaeta;
1039  float tru_deltaphi;
1040  float tru_siqr;
1041  float tru_siphi;
1042  float tru_sithe;
1043  float tru_six0;
1044  float tru_siy0;
1045  float tru_tpqr;
1046  float tru_tpphi;
1047  float tru_tpthe;
1048  float tru_tpx0;
1049  float tru_tpy0;
1050  float tru_charge;
1051  float tru_quality;
1052  float tru_chisq;
1053  float tru_ndf;
1054  float tru_nhits;
1055  float tru_layers;
1056  float tru_nmaps;
1057  float tru_nintt;
1058  float tru_ntpc;
1059  float tru_nmms;
1060  float tru_ntpc1;
1061  float tru_ntpc11;
1062  float tru_ntpc2;
1063  float tru_ntpc3;
1064  float tru_nlmaps;
1065  float tru_nlintt;
1066  float tru_nltpc;
1067  float tru_nlmms;
1068  float tru_vertexID;
1069  float tru_vx;
1070  float tru_vy;
1071  float tru_vz;
1072  float tru_dca2d;
1073  float tru_dca2dsigma;
1074  float tru_dca3dxy;
1075  float tru_dca3dxysigma;
1076  float tru_dca3dz;
1077  float tru_dca3dzsigma;
1078  float tru_pcax;
1079  float tru_pcay;
1080  float tru_pcaz;
1081  float tru_nfromtruth;
1082  float tru_nwrong;
1083  float tru_ntrumaps;
1084  float tru_nwrongmaps;
1085  float tru_ntruintt;
1086  float tru_nwrongintt;
1087  float tru_ntrutpc;
1088  float tru_nwrongtpc;
1089  float tru_ntrumms;
1090  float tru_nwrongmms;
1091  float tru_ntrutpc1;
1092  float tru_nwrongtpc1;
1093  float tru_ntrutpc11;
1094  float tru_nwrongtpc11;
1095  float tru_ntrutpc2;
1096  float tru_nwrongtpc2;
1097  float tru_ntrutpc3;
1098  float tru_nwrongtpc3;
1099  float tru_layersfromtruth;
1100  float tru_npedge;
1101  float tru_nredge;
1102  float tru_nbig;
1103  float tru_novlp;
1104  float tru_merr;
1105  float tru_msize;
1106  float tru_nhittpcall;
1107  float tru_nhittpcin;
1108  float tru_nhittpcmid;
1109  float tru_nhittpcout;
1110  float tru_nclusall;
1111  float tru_nclustpc;
1112  float tru_nclusintt;
1113  float tru_nclusmaps;
1114  float tru_nclusmms;
1115 
1116  float rec_event;
1117  float rec_seed;
1118  float rec_trackID;
1119  float rec_crossing;
1120  float rec_px;
1121  float rec_py;
1122  float rec_pz;
1123  float rec_pt;
1124  float rec_eta;
1125  float rec_phi;
1126  float rec_deltapt;
1127  float rec_deltaeta;
1128  float rec_deltaphi;
1129  float rec_siqr;
1130  float rec_siphi;
1131  float rec_sithe;
1132  float rec_six0;
1133  float rec_siy0;
1134  float rec_tpqr;
1135  float rec_tpphi;
1136  float rec_tpthe;
1137  float rec_tpx0;
1138  float rec_tpy0;
1139  float rec_charge;
1140  float rec_quality;
1141  float rec_chisq;
1142  float rec_ndf;
1143  float rec_nhits;
1144  float rec_nmaps;
1145  float rec_nintt;
1146  float rec_ntpc;
1147  float rec_nmms;
1148  float rec_ntpc1;
1149  float rec_ntpc11;
1150  float rec_ntpc2;
1151  float rec_ntpc3;
1152  float rec_nlmaps;
1153  float rec_nlintt;
1154  float rec_nltpc;
1155  float rec_nlmms;
1156  float rec_layers;
1157  float rec_vertexID;
1158  float rec_vx;
1159  float rec_vy;
1160  float rec_vz;
1161  float rec_dca2d;
1162  float rec_dca2dsigma;
1163  float rec_dca3dxy;
1164  float rec_dca3dxysigma;
1165  float rec_dca3dz;
1166  float rec_dca3dzsigma;
1167  float rec_pcax;
1168  float rec_pcay;
1169  float rec_pcaz;
1170  float rec_gtrackID;
1171  float rec_gflavor;
1172  float rec_gnhits;
1173  float rec_gnmaps;
1174  float rec_gnintt;
1175  float rec_gntpc;
1176  float rec_gnmms;
1177  float rec_gnlmaps;
1178  float rec_gnlintt;
1179  float rec_gnltpc;
1180  float rec_gnlmms;
1181  float rec_gpx;
1182  float rec_gpy;
1183  float rec_gpz;
1184  float rec_gpt;
1185  float rec_geta;
1186  float rec_gphi;
1187  float rec_gvx;
1188  float rec_gvy;
1189  float rec_gvz;
1190  float rec_gvt;
1191  float rec_gfpx;
1192  float rec_gfpy;
1193  float rec_gfpz;
1194  float rec_gfx;
1195  float rec_gfy;
1196  float rec_gfz;
1197  float rec_gembed;
1198  float rec_gprimary;
1199  float rec_nfromtruth;
1200  float rec_nwrong;
1201  float rec_ntrumaps;
1202  float rec_nwrongmaps;
1203  float rec_ntruintt;
1204  float rec_nwrongintt;
1205  float rec_ntrutpc;
1206  float rec_nwrongtpc;
1207  float rec_ntrumms;
1208  float rec_nwrongmms;
1209  float rec_ntrutpc1;
1210  float rec_nwrongtpc1;
1211  float rec_ntrutpc11;
1212  float rec_nwrongtpc11;
1213  float rec_ntrutpc2;
1214  float rec_nwrongtpc2;
1215  float rec_ntrutpc3;
1216  float rec_nwrongtpc3;
1217  float rec_layersfromtruth;
1218  float rec_npedge;
1219  float rec_nredge;
1220  float rec_nbig;
1221  float rec_novlp;
1222  float rec_merr;
1223  float rec_msize;
1224  float rec_nhittpcall;
1225  float rec_nhittpcin;
1226  float rec_nhittpcmid;
1227  float rec_nhittpcout;
1228  float rec_nclusall;
1229  float rec_nclustpc;
1230  float rec_nclusintt;
1231  float rec_nclusmaps;
1232  float rec_nclusmms;
1233 
1234  // Set branch addresses.
1235  m_ntOldTrue -> SetBranchAddress("event", &tru_event);
1236  m_ntOldTrue -> SetBranchAddress("seed", &tru_seed);
1237  m_ntOldTrue -> SetBranchAddress("gntracks", &tru_gntracks);
1238  m_ntOldTrue -> SetBranchAddress("gnchghad", &tru_gnchghad);
1239  m_ntOldTrue -> SetBranchAddress("gtrackID", &tru_gtrackID);
1240  m_ntOldTrue -> SetBranchAddress("gflavor", &tru_gflavor);
1241  m_ntOldTrue -> SetBranchAddress("gnhits", &tru_gnhits);
1242  m_ntOldTrue -> SetBranchAddress("gnmaps", &tru_gnmaps);
1243  m_ntOldTrue -> SetBranchAddress("gnintt", &tru_gnintt);
1244  m_ntOldTrue -> SetBranchAddress("gnmms", &tru_gnmms);
1245  m_ntOldTrue -> SetBranchAddress("gnintt1", &tru_gnintt1);
1246  m_ntOldTrue -> SetBranchAddress("gnintt2", &tru_gnintt2);
1247  m_ntOldTrue -> SetBranchAddress("gnintt3", &tru_gnintt3);
1248  m_ntOldTrue -> SetBranchAddress("gnintt4", &tru_gnintt4);
1249  m_ntOldTrue -> SetBranchAddress("gnintt5", &tru_gnintt5);
1250  m_ntOldTrue -> SetBranchAddress("gnintt6", &tru_gnintt6);
1251  m_ntOldTrue -> SetBranchAddress("gnintt7", &tru_gnintt7);
1252  m_ntOldTrue -> SetBranchAddress("gnintt8", &tru_gnintt8);
1253  m_ntOldTrue -> SetBranchAddress("gntpc", &tru_gntpc);
1254  m_ntOldTrue -> SetBranchAddress("gnlmaps", &tru_gnlmaps);
1255  m_ntOldTrue -> SetBranchAddress("gnlintt", &tru_gnlintt);
1256  m_ntOldTrue -> SetBranchAddress("gnltpc", &tru_gnltpc);
1257  m_ntOldTrue -> SetBranchAddress("gnlmms", &tru_gnlmms);
1258  m_ntOldTrue -> SetBranchAddress("gpx", &tru_gpx);
1259  m_ntOldTrue -> SetBranchAddress("gpy", &tru_gpy);
1260  m_ntOldTrue -> SetBranchAddress("gpz", &tru_gpz);
1261  m_ntOldTrue -> SetBranchAddress("gpt", &tru_gpt);
1262  m_ntOldTrue -> SetBranchAddress("geta", &tru_geta);
1263  m_ntOldTrue -> SetBranchAddress("gphi", &tru_gphi);
1264  m_ntOldTrue -> SetBranchAddress("gvx", &tru_gvx);
1265  m_ntOldTrue -> SetBranchAddress("gvy", &tru_gvy);
1266  m_ntOldTrue -> SetBranchAddress("gvz", &tru_gvz);
1267  m_ntOldTrue -> SetBranchAddress("gvt", &tru_gvt);
1268  m_ntOldTrue -> SetBranchAddress("gfpx", &tru_gfpx);
1269  m_ntOldTrue -> SetBranchAddress("gfpy", &tru_gfpy);
1270  m_ntOldTrue -> SetBranchAddress("gfpz", &tru_gfpz);
1271  m_ntOldTrue -> SetBranchAddress("gfx", &tru_gfx);
1272  m_ntOldTrue -> SetBranchAddress("gfy", &tru_gfy);
1273  m_ntOldTrue -> SetBranchAddress("gfz", &tru_gfz);
1274  m_ntOldTrue -> SetBranchAddress("gembed", &tru_gembed);
1275  m_ntOldTrue -> SetBranchAddress("gprimary", &tru_gprimary);
1276  m_ntOldTrue -> SetBranchAddress("trackID", &tru_trackID);
1277  m_ntOldTrue -> SetBranchAddress("px", &tru_px);
1278  m_ntOldTrue -> SetBranchAddress("py", &tru_py);
1279  m_ntOldTrue -> SetBranchAddress("pz", &tru_pz);
1280  m_ntOldTrue -> SetBranchAddress("pt", &tru_pt);
1281  m_ntOldTrue -> SetBranchAddress("eta", &tru_eta);
1282  m_ntOldTrue -> SetBranchAddress("phi", &tru_phi);
1283  m_ntOldTrue -> SetBranchAddress("deltapt", &tru_deltapt);
1284  m_ntOldTrue -> SetBranchAddress("deltaeta", &tru_deltaeta);
1285  m_ntOldTrue -> SetBranchAddress("deltaphi", &tru_deltaphi);
1286  m_ntOldTrue -> SetBranchAddress("siqr", &tru_siqr);
1287  m_ntOldTrue -> SetBranchAddress("siphi", &tru_siphi);
1288  m_ntOldTrue -> SetBranchAddress("sithe", &tru_sithe);
1289  m_ntOldTrue -> SetBranchAddress("six0", &tru_six0);
1290  m_ntOldTrue -> SetBranchAddress("siy0", &tru_siy0);
1291  m_ntOldTrue -> SetBranchAddress("tpqr", &tru_tpqr);
1292  m_ntOldTrue -> SetBranchAddress("tpphi", &tru_tpphi);
1293  m_ntOldTrue -> SetBranchAddress("tpthe", &tru_tpthe);
1294  m_ntOldTrue -> SetBranchAddress("tpx0", &tru_tpx0);
1295  m_ntOldTrue -> SetBranchAddress("tpy0", &tru_tpy0);
1296  m_ntOldTrue -> SetBranchAddress("charge", &tru_charge);
1297  m_ntOldTrue -> SetBranchAddress("quality", &tru_quality);
1298  m_ntOldTrue -> SetBranchAddress("chisq", &tru_chisq);
1299  m_ntOldTrue -> SetBranchAddress("ndf", &tru_ndf);
1300  m_ntOldTrue -> SetBranchAddress("nhits", &tru_nhits);
1301  m_ntOldTrue -> SetBranchAddress("layers", &tru_layers);
1302  m_ntOldTrue -> SetBranchAddress("nmaps", &tru_nmaps);
1303  m_ntOldTrue -> SetBranchAddress("nintt", &tru_nintt);
1304  m_ntOldTrue -> SetBranchAddress("ntpc", &tru_ntpc);
1305  m_ntOldTrue -> SetBranchAddress("nmms", &tru_nmms);
1306  m_ntOldTrue -> SetBranchAddress("ntpc1", &tru_ntpc1);
1307  m_ntOldTrue -> SetBranchAddress("ntpc11", &tru_ntpc11);
1308  m_ntOldTrue -> SetBranchAddress("ntpc2", &tru_ntpc2);
1309  m_ntOldTrue -> SetBranchAddress("ntpc3", &tru_ntpc3);
1310  m_ntOldTrue -> SetBranchAddress("nlmaps", &tru_nlmaps);
1311  m_ntOldTrue -> SetBranchAddress("nlintt", &tru_nlintt);
1312  m_ntOldTrue -> SetBranchAddress("nltpc", &tru_nltpc);
1313  m_ntOldTrue -> SetBranchAddress("nlmms", &tru_nlmms);
1314  m_ntOldTrue -> SetBranchAddress("vertexID", &tru_vertexID);
1315  m_ntOldTrue -> SetBranchAddress("vx", &tru_vx);
1316  m_ntOldTrue -> SetBranchAddress("vy", &tru_vy);
1317  m_ntOldTrue -> SetBranchAddress("vz", &tru_vz);
1318  m_ntOldTrue -> SetBranchAddress("dca2d", &tru_dca2d);
1319  m_ntOldTrue -> SetBranchAddress("dca2dsigma", &tru_dca2dsigma);
1320  m_ntOldTrue -> SetBranchAddress("dca3dxy", &tru_dca3dxy);
1321  m_ntOldTrue -> SetBranchAddress("dca3dxysigma", &tru_dca3dxysigma);
1322  m_ntOldTrue -> SetBranchAddress("dca3dz", &tru_dca3dz);
1323  m_ntOldTrue -> SetBranchAddress("dca3dzsigma", &tru_dca3dzsigma);
1324  m_ntOldTrue -> SetBranchAddress("pcax", &tru_pcax);
1325  m_ntOldTrue -> SetBranchAddress("pcay", &tru_pcay);
1326  m_ntOldTrue -> SetBranchAddress("pcaz", &tru_pcaz);
1327  m_ntOldTrue -> SetBranchAddress("nfromtruth", &tru_nfromtruth);
1328  m_ntOldTrue -> SetBranchAddress("nwrong", &tru_nwrong);
1329  m_ntOldTrue -> SetBranchAddress("ntrumaps", &tru_ntrumaps);
1330  m_ntOldTrue -> SetBranchAddress("nwrongmaps", &tru_nwrongmaps);
1331  m_ntOldTrue -> SetBranchAddress("ntruintt", &tru_ntruintt);
1332  m_ntOldTrue -> SetBranchAddress("nwrongintt", &tru_nwrongintt);
1333  m_ntOldTrue -> SetBranchAddress("ntrutpc", &tru_ntrutpc);
1334  m_ntOldTrue -> SetBranchAddress("nwrongtpc", &tru_nwrongtpc);
1335  m_ntOldTrue -> SetBranchAddress("ntrumms", &tru_ntrumms);
1336  m_ntOldTrue -> SetBranchAddress("nwrongmms", &tru_nwrongmms);
1337  m_ntOldTrue -> SetBranchAddress("ntrutpc1", &tru_ntrutpc1);
1338  m_ntOldTrue -> SetBranchAddress("nwrongtpc1", &tru_nwrongtpc1);
1339  m_ntOldTrue -> SetBranchAddress("ntrutpc11", &tru_ntrutpc11);
1340  m_ntOldTrue -> SetBranchAddress("nwrongtpc11", &tru_nwrongtpc11);
1341  m_ntOldTrue -> SetBranchAddress("ntrutpc2", &tru_ntrutpc2);
1342  m_ntOldTrue -> SetBranchAddress("nwrongtpc2", &tru_nwrongtpc2);
1343  m_ntOldTrue -> SetBranchAddress("ntrutpc3", &tru_ntrutpc3);
1344  m_ntOldTrue -> SetBranchAddress("nwrongtpc3", &tru_nwrongtpc3);
1345  m_ntOldTrue -> SetBranchAddress("layersfromtruth", &tru_layersfromtruth);
1346  m_ntOldTrue -> SetBranchAddress("npedge", &tru_npedge);
1347  m_ntOldTrue -> SetBranchAddress("nredge", &tru_nredge);
1348  m_ntOldTrue -> SetBranchAddress("nbig", &tru_nbig);
1349  m_ntOldTrue -> SetBranchAddress("novlp", &tru_novlp);
1350  m_ntOldTrue -> SetBranchAddress("merr", &tru_merr);
1351  m_ntOldTrue -> SetBranchAddress("msize", &tru_msize);
1352  m_ntOldTrue -> SetBranchAddress("nhittpcall", &tru_nhittpcall);
1353  m_ntOldTrue -> SetBranchAddress("nhittpcin", &tru_nhittpcin);
1354  m_ntOldTrue -> SetBranchAddress("nhittpcmid", &tru_nhittpcmid);
1355  m_ntOldTrue -> SetBranchAddress("nhittpcout", &tru_nhittpcout);
1356  m_ntOldTrue -> SetBranchAddress("nclusall", &tru_nclusall);
1357  m_ntOldTrue -> SetBranchAddress("nclustpc", &tru_nclustpc);
1358  m_ntOldTrue -> SetBranchAddress("nclusintt", &tru_nclusintt);
1359  m_ntOldTrue -> SetBranchAddress("nclusmaps", &tru_nclusmaps);
1360  m_ntOldTrue -> SetBranchAddress("nclusmms", &tru_nclusmms);
1361 
1362  m_ntOldReco -> SetBranchAddress("event", &rec_event);
1363  m_ntOldReco -> SetBranchAddress("seed", &rec_seed);
1364  m_ntOldReco -> SetBranchAddress("trackID", &rec_trackID);
1365  m_ntOldReco -> SetBranchAddress("crossing", &rec_crossing);
1366  m_ntOldReco -> SetBranchAddress("px", &rec_px);
1367  m_ntOldReco -> SetBranchAddress("py", &rec_py);
1368  m_ntOldReco -> SetBranchAddress("pz", &rec_pz);
1369  m_ntOldReco -> SetBranchAddress("pt", &rec_pt);
1370  m_ntOldReco -> SetBranchAddress("eta", &rec_eta);
1371  m_ntOldReco -> SetBranchAddress("phi", &rec_phi);
1372  m_ntOldReco -> SetBranchAddress("deltapt", &rec_deltapt);
1373  m_ntOldReco -> SetBranchAddress("deltaeta", &rec_deltaeta);
1374  m_ntOldReco -> SetBranchAddress("deltaphi", &rec_deltaphi);
1375  m_ntOldReco -> SetBranchAddress("siqr", &rec_siqr);
1376  m_ntOldReco -> SetBranchAddress("siphi", &rec_siphi);
1377  m_ntOldReco -> SetBranchAddress("sithe", &rec_sithe);
1378  m_ntOldReco -> SetBranchAddress("six0", &rec_six0);
1379  m_ntOldReco -> SetBranchAddress("siy0", &rec_siy0);
1380  m_ntOldReco -> SetBranchAddress("tpqr", &rec_tpqr);
1381  m_ntOldReco -> SetBranchAddress("tpphi", &rec_tpphi);
1382  m_ntOldReco -> SetBranchAddress("tpthe", &rec_tpthe);
1383  m_ntOldReco -> SetBranchAddress("tpx0", &rec_tpx0);
1384  m_ntOldReco -> SetBranchAddress("tpy0", &rec_tpy0);
1385  m_ntOldReco -> SetBranchAddress("charge", &rec_charge);
1386  m_ntOldReco -> SetBranchAddress("quality", &rec_quality);
1387  m_ntOldReco -> SetBranchAddress("chisq", &rec_chisq);
1388  m_ntOldReco -> SetBranchAddress("ndf", &rec_ndf);
1389  m_ntOldReco -> SetBranchAddress("nhits", &rec_nhits);
1390  m_ntOldReco -> SetBranchAddress("nmaps", &rec_nmaps);
1391  m_ntOldReco -> SetBranchAddress("nintt", &rec_nintt);
1392  m_ntOldReco -> SetBranchAddress("ntpc", &rec_ntpc);
1393  m_ntOldReco -> SetBranchAddress("nmms", &rec_nmms);
1394  m_ntOldReco -> SetBranchAddress("ntpc1", &rec_ntpc1);
1395  m_ntOldReco -> SetBranchAddress("ntpc11", &rec_ntpc11);
1396  m_ntOldReco -> SetBranchAddress("ntpc2", &rec_ntpc2);
1397  m_ntOldReco -> SetBranchAddress("ntpc3", &rec_ntpc3);
1398  m_ntOldReco -> SetBranchAddress("nlmaps", &rec_nlmaps);
1399  m_ntOldReco -> SetBranchAddress("nlintt", &rec_nlintt);
1400  m_ntOldReco -> SetBranchAddress("nltpc", &rec_nltpc);
1401  m_ntOldReco -> SetBranchAddress("nlmms", &rec_nlmms);
1402  m_ntOldReco -> SetBranchAddress("layers", &rec_layers);
1403  m_ntOldReco -> SetBranchAddress("vertexID", &rec_vertexID);
1404  m_ntOldReco -> SetBranchAddress("vx", &rec_vx);
1405  m_ntOldReco -> SetBranchAddress("vy", &rec_vy);
1406  m_ntOldReco -> SetBranchAddress("vz", &rec_vz);
1407  m_ntOldReco -> SetBranchAddress("dca2d", &rec_dca2d);
1408  m_ntOldReco -> SetBranchAddress("dca2dsigma", &rec_dca2dsigma);
1409  m_ntOldReco -> SetBranchAddress("dca3dxy", &rec_dca3dxy);
1410  m_ntOldReco -> SetBranchAddress("dca3dxysigma", &rec_dca3dxysigma);
1411  m_ntOldReco -> SetBranchAddress("dca3dz", &rec_dca3dz);
1412  m_ntOldReco -> SetBranchAddress("dca3dzsigma", &rec_dca3dzsigma);
1413  m_ntOldReco -> SetBranchAddress("pcax", &rec_pcax);
1414  m_ntOldReco -> SetBranchAddress("pcay", &rec_pcay);
1415  m_ntOldReco -> SetBranchAddress("pcaz", &rec_pcaz);
1416  m_ntOldReco -> SetBranchAddress("gtrackID", &rec_gtrackID);
1417  m_ntOldReco -> SetBranchAddress("gflavor", &rec_gflavor);
1418  m_ntOldReco -> SetBranchAddress("gnhits", &rec_gnhits);
1419  m_ntOldReco -> SetBranchAddress("gnmaps", &rec_gnmaps);
1420  m_ntOldReco -> SetBranchAddress("gnintt", &rec_gnintt);
1421  m_ntOldReco -> SetBranchAddress("gntpc", &rec_gntpc);
1422  m_ntOldReco -> SetBranchAddress("gnmms", &rec_gnmms);
1423  m_ntOldReco -> SetBranchAddress("gnlmaps", &rec_gnlmaps);
1424  m_ntOldReco -> SetBranchAddress("gnlintt", &rec_gnlintt);
1425  m_ntOldReco -> SetBranchAddress("gnltpc", &rec_gnltpc);
1426  m_ntOldReco -> SetBranchAddress("gnlmms", &rec_gnlmms);
1427  m_ntOldReco -> SetBranchAddress("gpx", &rec_gpx);
1428  m_ntOldReco -> SetBranchAddress("gpy", &rec_gpy);
1429  m_ntOldReco -> SetBranchAddress("gpz", &rec_gpz);
1430  m_ntOldReco -> SetBranchAddress("gpt", &rec_gpt);
1431  m_ntOldReco -> SetBranchAddress("geta", &rec_geta);
1432  m_ntOldReco -> SetBranchAddress("gphi", &rec_gphi);
1433  m_ntOldReco -> SetBranchAddress("gvx", &rec_gvx);
1434  m_ntOldReco -> SetBranchAddress("gvy", &rec_gvy);
1435  m_ntOldReco -> SetBranchAddress("gvz", &rec_gvz);
1436  m_ntOldReco -> SetBranchAddress("gvt", &rec_gvt);
1437  m_ntOldReco -> SetBranchAddress("gfpx", &rec_gfpx);
1438  m_ntOldReco -> SetBranchAddress("gfpy", &rec_gfpy);
1439  m_ntOldReco -> SetBranchAddress("gfpz", &rec_gfpz);
1440  m_ntOldReco -> SetBranchAddress("gfx", &rec_gfx);
1441  m_ntOldReco -> SetBranchAddress("gfy", &rec_gfy);
1442  m_ntOldReco -> SetBranchAddress("gfz", &rec_gfz);
1443  m_ntOldReco -> SetBranchAddress("gembed", &rec_gembed);
1444  m_ntOldReco -> SetBranchAddress("gprimary", &rec_gprimary);
1445  m_ntOldReco -> SetBranchAddress("nfromtruth", &rec_nfromtruth);
1446  m_ntOldReco -> SetBranchAddress("nwrong", &rec_nwrong);
1447  m_ntOldReco -> SetBranchAddress("ntrumaps", &rec_ntrumaps);
1448  m_ntOldReco -> SetBranchAddress("nwrongmaps", &rec_nwrongmaps);
1449  m_ntOldReco -> SetBranchAddress("ntruintt", &rec_ntruintt);
1450  m_ntOldReco -> SetBranchAddress("nwrongintt", &rec_nwrongintt);
1451  m_ntOldReco -> SetBranchAddress("ntrutpc", &rec_ntrutpc);
1452  m_ntOldReco -> SetBranchAddress("nwrongtpc", &rec_nwrongtpc);
1453  m_ntOldReco -> SetBranchAddress("ntrumms", &rec_ntrumms);
1454  m_ntOldReco -> SetBranchAddress("nwrongmms", &rec_nwrongmms);
1455  m_ntOldReco -> SetBranchAddress("ntrutpc1", &rec_ntrutpc1);
1456  m_ntOldReco -> SetBranchAddress("nwrongtpc1", &rec_nwrongtpc1);
1457  m_ntOldReco -> SetBranchAddress("ntrutpc11", &rec_ntrutpc11);
1458  m_ntOldReco -> SetBranchAddress("nwrongtpc11", &rec_nwrongtpc11);
1459  m_ntOldReco -> SetBranchAddress("ntrutpc2", &rec_ntrutpc2);
1460  m_ntOldReco -> SetBranchAddress("nwrongtpc2", &rec_nwrongtpc2);
1461  m_ntOldReco -> SetBranchAddress("ntrutpc3", &rec_ntrutpc3);
1462  m_ntOldReco -> SetBranchAddress("nwrongtpc3", &rec_nwrongtpc3);
1463  m_ntOldReco -> SetBranchAddress("layersfromtruth", &rec_layersfromtruth);
1464  m_ntOldReco -> SetBranchAddress("npedge", &rec_npedge);
1465  m_ntOldReco -> SetBranchAddress("nredge", &rec_nredge);
1466  m_ntOldReco -> SetBranchAddress("nbig", &rec_nbig);
1467  m_ntOldReco -> SetBranchAddress("novlp", &rec_novlp);
1468  m_ntOldReco -> SetBranchAddress("merr", &rec_merr);
1469  m_ntOldReco -> SetBranchAddress("msize", &rec_msize);
1470  m_ntOldReco -> SetBranchAddress("nhittpcall", &rec_nhittpcall);
1471  m_ntOldReco -> SetBranchAddress("nhittpcin", &rec_nhittpcin);
1472  m_ntOldReco -> SetBranchAddress("nhittpcmid", &rec_nhittpcmid);
1473  m_ntOldReco -> SetBranchAddress("nhittpcout", &rec_nhittpcout);
1474  m_ntOldReco -> SetBranchAddress("nclusall", &rec_nclusall);
1475  m_ntOldReco -> SetBranchAddress("nclustpc", &rec_nclustpc);
1476  m_ntOldReco -> SetBranchAddress("nclusintt", &rec_nclusintt);
1477  m_ntOldReco -> SetBranchAddress("nclusmaps", &rec_nclusmaps);
1478  m_ntOldReco -> SetBranchAddress("nclusmms", &rec_nclusmms);
1479  cout << " Set input branches." << endl;
1480 
1481  // grab no. of entries
1482  const int64_t nTrueEntries = m_ntOldTrue -> GetEntries();
1483  const int64_t nRecoEntries = m_ntOldReco -> GetEntries();
1484  cout << " Beginning truth particle loop: " << nTrueEntries << " to process" << endl;
1485 
1486  // loop over truth particles
1487  int64_t nTrueBytes = 0;
1488  for (int64_t iTrueEntry = 0; iTrueEntry < nTrueEntries; iTrueEntry++) {
1489 
1490  // grab truth particle entry
1491  const int64_t trueBytes = m_ntOldTrue -> GetEntry(iTrueEntry);
1492  if (trueBytes < 0) {
1493  cerr << "PANIC: issue with entry " << iTrueEntry << "! Aborting loop!\n" << endl;
1494  break;
1495  } else {
1496  nTrueBytes += trueBytes;
1497  }
1498 
1499  const int64_t iTrueProg = iTrueEntry + 1;
1500  if (iTrueProg == nTrueEntries) {
1501  cout << " Processing entry " << iTrueProg << "/" << nTrueEntries << "..." << endl;
1502  } else {
1503  cout << " Processing entry " << iTrueProg << "/" << nTrueEntries << "...\r" << flush;
1504  }
1505 
1506  // skip nan's
1507  if (isnan(tru_trackID)) continue;
1508 
1509  // check if near sector
1510  bool isNearSector = false;
1511  if (m_config.doPhiCut) {
1512  isNearSector = IsNearSector(tru_gphi);
1513  }
1514 
1515  // apply cuts
1516  const bool isPrimary = (tru_gprimary == 1);
1517  const bool isInZVtxCut = ((tru_gvz > m_config.zVtxRange.first) && (tru_gvz < m_config.zVtxRange.second));
1518  if (m_config.useOnlyPrimTrks && !isPrimary) continue;
1519  if (m_config.doZVtxCut && !isInZVtxCut) continue;
1520  if (m_config.doPhiCut && isNearSector) continue;
1521 
1522  // run calculations
1523  const double tru_gntot = tru_gnintt + tru_gnmaps + tru_gntpc;
1524 
1525  // bundle histogram values to fill
1526  STrackMatcherComparatorHistContent content {
1527  tru_gntot,
1528  tru_gnintt,
1529  tru_gnmaps,
1530  tru_gntpc,
1531  1.,
1532  1.,
1533  1.,
1534  1.,
1535  tru_gphi,
1536  tru_geta,
1537  tru_gpt,
1538  1.,
1539  1.,
1540  1.,
1541  1.,
1542  1.,
1543  1.,
1544  1.,
1545  1.
1546  };
1547 
1548  // fill truth histograms
1549  FillHistogram1D(content, Type::Truth, m_vecOldHists1D);
1550  FillHistogram2D(content, Type::Truth, Comp::VsTruthPt, tru_gpt, m_vecOldHists2D);
1551  FillHistogram2D(content, Type::Truth, Comp::VsNumTpc, tru_gntpc, m_vecOldHists2D);
1552  } // end truth particle loop
1553 
1554  // announce next entry loop
1555  cout << " Finished truth particle loop.\n"
1556  << " Beginning reconstructed track loop: " << nRecoEntries << " to process"
1557  << endl;
1558 
1559  // loop over reco tracks
1560  int64_t nRecoBytes = 0;
1561  for (int64_t iRecoEntry = 0; iRecoEntry < nRecoEntries; iRecoEntry++) {
1562 
1563  // grab reco track entry
1564  const int64_t recoBytes = m_ntOldReco -> GetEntry(iRecoEntry);
1565  if (recoBytes < 0) {
1566  cerr << "PANIC: issue with entry " << iRecoEntry << "! Aborting loop!\n" << endl;
1567  break;
1568  } else {
1569  nRecoBytes += recoBytes;
1570  }
1571 
1572  const int64_t iRecoProg = iRecoEntry + 1;
1573  if (iRecoProg == nRecoEntries) {
1574  cout << " Processing entry " << iRecoProg << "/" << nRecoEntries << "..." << endl;
1575  } else {
1576  cout << " Processing entry " << iRecoProg << "/" << nRecoEntries << "...\r" << flush;
1577  }
1578 
1579  // skip nan's
1580  if (isnan(rec_gpt)) continue;
1581 
1582  // check if near sector
1583  bool isNearSector = false;
1584  if (m_config.doPhiCut) {
1585  isNearSector = IsNearSector(rec_phi);
1586  }
1587 
1588  // apply cuts
1589  const bool isPrimary = (rec_gprimary == 1);
1590  const bool isInZVtxCut = ((rec_vz > m_config.zVtxRange.first) && (rec_vz < m_config.zVtxRange.second));
1591  if (m_config.useOnlyPrimTrks && !isPrimary) continue;
1592  if (m_config.doZVtxCut && !isInZVtxCut) continue;
1593  if (m_config.doPhiCut && isNearSector) continue;
1594 
1595  // run calculations
1596  const double rec_ntot = rec_nintt + rec_nmaps + rec_ntpc;
1597  const double rec_gntot = rec_gnintt + rec_gnmaps + rec_gntpc;
1598  const double rec_rtot = rec_ntot / rec_gntot;
1599  const double rec_rintt = rec_nintt / rec_gnintt;
1600  const double rec_rmaps = rec_nmaps / rec_gnmaps;
1601  const double rec_rtpc = rec_ntpc / rec_gntpc;
1602  const double rec_ptfrac = rec_pt / rec_gpt;
1603  const double rec_ptper = rec_deltapt / rec_pt;
1604  const double rec_etaper = rec_deltaeta / rec_eta;
1605  const double rec_phiper = rec_deltaphi / rec_phi;
1606  const double rec_ptres = abs(rec_pt - rec_gpt) / rec_gpt;
1607  const double rec_etares = abs(rec_eta - rec_geta) / rec_geta;
1608  const double rec_phires = abs(rec_phi - rec_gphi) / rec_gphi;
1609 
1610  // bundle histogram values to fill
1611  STrackMatcherComparatorHistContent content {
1612  rec_ntot,
1613  rec_nintt,
1614  rec_nmaps,
1615  rec_ntpc,
1616  rec_rtot,
1617  rec_rintt,
1618  rec_rmaps,
1619  rec_rtpc,
1620  rec_phi,
1621  rec_eta,
1622  rec_pt,
1623  rec_ptfrac,
1624  rec_quality,
1625  rec_ptper,
1626  rec_etaper,
1627  rec_phiper,
1628  rec_ptres,
1629  rec_etares,
1630  rec_phires
1631  };
1632 
1633  // fill all matched reco histograms
1634  FillHistogram1D(content, Type::Track, m_vecOldHists1D);
1635  FillHistogram2D(content, Type::Track, Comp::VsTruthPt, rec_gpt, m_vecOldHists2D);
1636  FillHistogram2D(content, Type::Track, Comp::VsNumTpc, rec_gntpc, m_vecOldHists2D);
1637 
1638  // fill weird and normal matched reco 1D histograms
1639  const bool isNormalTrack = ((rec_ptfrac >= m_config.oddPtFrac.first) && (rec_ptfrac <= m_config.oddPtFrac.second));
1640  if (isNormalTrack) {
1641  FillHistogram1D(content, Type::Normal, m_vecOldHists1D);
1642  FillHistogram2D(content, Type::Normal, Comp::VsTruthPt, rec_gpt, m_vecOldHists2D);
1643  FillHistogram2D(content, Type::Normal, Comp::VsNumTpc, rec_gntpc, m_vecOldHists2D);
1644  } else {
1645  FillHistogram1D(content, Type::Weird, m_vecOldHists1D);
1646  FillHistogram2D(content, Type::Weird, Comp::VsTruthPt, rec_gpt, m_vecOldHists2D);
1647  FillHistogram2D(content, Type::Weird, Comp::VsNumTpc, rec_gntpc, m_vecOldHists2D);
1648  }
1649  } // end reco track loop
1650 
1651  // announce method end and return
1652  cout << " Finished reconstructed track loop.\n"
1653  << " Finished getting histograms from old matcher track tuple."
1654  << endl;
1655  return;
1656 
1657 } // end 'GetOldTupleHists()'
1658 
1659 
1660 
1661 void STrackMatcherComparator::MakeRatiosAndPlots(
1662  const vector<vector<TH1D*>> vecNewHists1D,
1663  const vector<vector<vector<TH2D*>>> vecNewHists2D,
1664  const int iDir,
1665  const string sLabel
1666 ) {
1667 
1668  // announce start of routine
1669  const size_t nHistRows1D = m_vecOldHists1D.size();
1670  const size_t nHistRows2D = m_vecOldHists2D.size();
1671  const size_t nHistTypes1D = m_vecOldHists1D.front().size();
1672  const size_t nHistTypes2D = m_vecOldHists2D.front().size();
1673  const size_t nHist2D = m_vecOldHists2D.front().front().size();
1674  cout << " Making ratios and plots:" << endl;
1675 
1676  // create 1d canvas names
1677  vector<vector<string>> vecCanvasNames1D(nHistRows1D);
1678  for (size_t iHistRow1D = 0; iHistRow1D < nHistRows1D; iHistRow1D++) {
1679  vecCanvasNames1D[iHistRow1D].resize(nHistTypes1D);
1680  for (size_t iHistType1D = 0; iHistType1D < nHistTypes1D; iHistType1D++) {
1681  vecCanvasNames1D[iHistRow1D][iHistType1D] = m_hist.vecNameBase[iHistRow1D][iHistType1D];
1682  vecCanvasNames1D[iHistRow1D][iHistType1D].replace(0, 1, "c");
1683  }
1684  }
1685 
1686  // create 2d canvas names
1687  vector<vector<vector<string>>> vecCanvasNames2D(nHistRows2D);
1688  for (size_t iHistRow2D = 0; iHistRow2D < nHistRows2D; iHistRow2D++) {
1689  vecCanvasNames2D[iHistRow2D].resize(nHistRows2D);
1690  for (size_t iHistType2D = 0; iHistType2D < nHistTypes2D; iHistType2D++) {
1691  vecCanvasNames2D[iHistRow2D][iHistType2D].resize(nHist2D);
1692  for (size_t iHist2D = 0; iHist2D < nHist2D; iHist2D++) {
1693  vecCanvasNames2D[iHistRow2D][iHistType2D][iHist2D] = m_hist.vecNameBase[iHistRow2D][iHistType2D];
1694  vecCanvasNames2D[iHistRow2D][iHistType2D][iHist2D] += m_hist.vecVsModifiers[iHist2D];
1695  vecCanvasNames2D[iHistRow2D][iHistType2D][iHist2D].replace(0, 1, "c");
1696  }
1697  }
1698  }
1699  cout << " Set up canvas names." << endl;
1700 
1701  // normalize histograms and set axes' ranges as needed ----------------------
1702 
1703  // normalize by integral if needed
1704  if (m_config.doIntNorm) {
1705  for (auto oldHistRow1D : m_vecOldHists1D) {
1706  for (auto hOldHist1D : oldHistRow1D) {
1707  const double intOld1D = hOldHist1D -> Integral();
1708  if (intOld1D > 0.) {
1709  hOldHist1D -> Scale(1. / intOld1D);
1710  }
1711  }
1712  }
1713  for (auto newHistRow1D : vecNewHists1D) {
1714  for (auto hNewHist1D : newHistRow1D) {
1715  const double intNew1D = hNewHist1D -> Integral();
1716  if (intNew1D > 0.) {
1717  hNewHist1D -> Scale(1. / intNew1D);
1718  }
1719  }
1720  }
1721  for (auto oldHistRow2D : m_vecOldHists2D) {
1722  for (auto oldHistTypes2D : oldHistRow2D) {
1723  for (auto hOldHist2D : oldHistTypes2D) {
1724  const double intOld2D = hOldHist2D -> Integral();
1725  if (intOld2D > 0.) {
1726  hOldHist2D -> Scale(1. / intOld2D);
1727  }
1728  }
1729  }
1730  }
1731  for (auto newHistRow2D : vecNewHists2D) {
1732  for (auto newHistTypes2D : newHistRow2D) {
1733  for (auto hNewHist2D : newHistTypes2D) {
1734  const double intNew2D = hNewHist2D -> Integral();
1735  if (intNew2D > 0.) {
1736  hNewHist2D -> Scale(1. / intNew2D);
1737  }
1738  }
1739  }
1740  }
1741  cout << " Normalized histograms." << endl;
1742  }
1743 
1744  // set z-axis ranges if needed
1745  if (m_config.matchVertScales) {
1746  for (size_t iHistRow2D = 0; iHistRow2D < nHistRows2D; iHistRow2D++) {
1747  for (size_t iHistType2D = 0; iHistType2D < nHistTypes2D; iHistType2D++) {
1748  for (size_t iHist2D = 0; iHist2D < nHist2D; iHist2D++) {
1749  const float oldMin = m_vecOldHists2D[iHistRow2D][iHistType2D][iHist2D] -> GetMinimum(0.);
1750  const float oldMax = m_vecOldHists2D[iHistRow2D][iHistType2D][iHist2D] -> GetMaximum();
1751  const float newMin = vecNewHists2D[iHistRow2D][iHistType2D][iHist2D] -> GetMinimum(0.);
1752  const float newMax = vecNewHists2D[iHistRow2D][iHistType2D][iHist2D] -> GetMaximum();
1753  const float setMin = min(oldMin, newMin);
1754  const float setMax = max(oldMax, newMax);
1755  m_vecOldHists2D[iHistRow2D][iHistType2D][iHist2D] -> GetZaxis() -> SetRangeUser(setMin, setMax);
1756  vecNewHists2D[iHistRow2D][iHistType2D][iHist2D] -> GetZaxis() -> SetRangeUser(setMin, setMax);
1757  }
1758  }
1759  }
1760  cout << " Adjusted z-axis scales to match." << endl;
1761  }
1762 
1763  // pick relevant count
1764  string countUse("");
1765  if (m_config.doIntNorm) {
1766  countUse = m_config.norm;
1767  } else {
1768  countUse = m_config.count;
1769  }
1770 
1771  // calculate ratios ---------------------------------------------------------
1772 
1773  // calculate 1d ratios
1774  vector<vector<TH1D*>> vecRatios1D(nHistRows1D);
1775  for (size_t iHistRow1D = 0; iHistRow1D < nHistRows1D; iHistRow1D++) {
1776  vecRatios1D[iHistRow1D].resize(nHistTypes1D);
1777  for (size_t iHist1D = 0; iHist1D < nHistTypes1D; iHist1D++) {
1778 
1779  // make histogram name
1780  TString sRatioName = m_vecOldHists1D[iHistRow1D][iHist1D] -> GetName();
1781  sRatioName.Append("_");
1782  sRatioName.Append(sLabel.data());
1783 
1784  vecRatios1D[iHistRow1D][iHist1D] = (TH1D*) m_vecOldHists1D[iHistRow1D][iHist1D] -> Clone();
1785  vecRatios1D[iHistRow1D][iHist1D] -> SetName(sRatioName.Data());
1786  vecRatios1D[iHistRow1D][iHist1D] -> Reset("ICES");
1787  vecRatios1D[iHistRow1D][iHist1D] -> Divide(m_vecOldHists1D[iHistRow1D][iHist1D], vecNewHists1D[iHistRow1D][iHist1D], 1., 1.);
1788  }
1789  }
1790 
1791  // calculate 2d ratios
1792  vector<vector<vector<TH2D*>>> vecRatios2D(nHistRows2D);
1793  for (size_t iHistRow2D = 0; iHistRow2D < nHistRows2D; iHistRow2D++) {
1794  vecRatios2D[iHistRow2D].resize(nHistTypes2D);
1795  for (size_t iHistType2D = 0; iHistType2D < nHistTypes2D; iHistType2D++) {
1796  vecRatios2D[iHistRow2D][iHistType2D].resize(nHist2D);
1797  for (size_t iHist2D = 0; iHist2D < nHist2D; iHist2D++) {
1798 
1799  // make histogram name
1800  TString sRatioName = m_vecOldHists2D[iHistRow2D][iHistType2D][iHist2D] -> GetName();
1801  sRatioName.Append("_");
1802  sRatioName.Append(sLabel.data());
1803 
1804  vecRatios2D[iHistRow2D][iHistType2D][iHist2D] = (TH2D*) m_vecOldHists2D[iHistRow2D][iHistType2D][iHist2D] -> Clone();
1805  vecRatios2D[iHistRow2D][iHistType2D][iHist2D] -> SetName(sRatioName.Data());
1806  vecRatios2D[iHistRow2D][iHistType2D][iHist2D] -> Reset("ICES");
1807  vecRatios2D[iHistRow2D][iHistType2D][iHist2D] -> Divide(m_vecOldHists2D[iHistRow2D][iHistType2D][iHist2D], vecNewHists2D[iHistRow2D][iHistType2D][iHist2D], 1., 1.);
1808  }
1809  }
1810  }
1811  cout << " Calculated ratios." << endl;
1812 
1813  // set histogram styles -----------------------------------------------------
1814 
1815  // style options
1816  const uint32_t fTxt(42);
1817  const uint32_t fAln(12);
1818  const uint32_t fCnt(1);
1819  const float fLabH[m_const.nAxes] = {0.04, 0.04, 0.03};
1820  const float fLabR1[m_const.nAxes] = {0.074, 0.074, 0.056};
1821  const float fLabR2[m_const.nAxes] = {0.04, 0.04, 0.03};
1822  const float fTitH[m_const.nAxes] = {0.04, 0.04, 0.04};
1823  const float fTitR1[m_const.nAxes] = {0.074, 0.074, 0.056};
1824  const float fTitR2[m_const.nAxes] = {0.04, 0.04, 0.04};
1825  const float fOffTH[m_const.nAxes] = {1.0, 1.3, 1.2};
1826  const float fOffTR1[m_const.nAxes] = {0.8, 0.8, 1.0};
1827  const float fOffTR2[m_const.nAxes] = {1.0, 1.3, 1.2};
1828  const float fOffLH[m_const.nAxes] = {0.005, 0.005, 0.000};
1829  const float fOffLR1[m_const.nAxes] = {0.005, 0.005, 0.000};
1830  const float fOffLR2[m_const.nAxes] = {0.005, 0.005, 0.000};
1831 
1832  // set old histogram styles
1833  for (auto oldHistRow1D : m_vecOldHists1D) {
1834  for (auto hOldHist1D : oldHistRow1D) {
1835  hOldHist1D -> SetMarkerColor(m_config.fCol.first);
1836  hOldHist1D -> SetMarkerStyle(m_config.fMar.first);
1837  hOldHist1D -> SetFillColor(m_config.fCol.first);
1838  hOldHist1D -> SetFillStyle(m_config.fFil.first);
1839  hOldHist1D -> SetLineColor(m_config.fCol.first);
1840  hOldHist1D -> SetLineStyle(m_config.fLin.first);
1841  hOldHist1D -> SetLineWidth(m_config.fWid.first);
1842  hOldHist1D -> SetTitle("");
1843  hOldHist1D -> SetTitleFont(fTxt);
1844  hOldHist1D -> GetXaxis() -> SetTitleFont(fTxt);
1845  hOldHist1D -> GetXaxis() -> SetTitleSize(fTitH[0]);
1846  hOldHist1D -> GetXaxis() -> SetTitleOffset(fOffTH[0]);
1847  hOldHist1D -> GetXaxis() -> SetLabelFont(fTxt);
1848  hOldHist1D -> GetXaxis() -> SetLabelSize(fLabH[0]);
1849  hOldHist1D -> GetXaxis() -> SetLabelOffset(fOffLH[0]);
1850  hOldHist1D -> GetXaxis() -> CenterTitle(fCnt);
1851  hOldHist1D -> GetYaxis() -> SetTitle(countUse.data());
1852  hOldHist1D -> GetYaxis() -> SetTitleFont(fTxt);
1853  hOldHist1D -> GetYaxis() -> SetTitleSize(fTitH[1]);
1854  hOldHist1D -> GetYaxis() -> SetTitleOffset(fOffTH[1]);
1855  hOldHist1D -> GetYaxis() -> SetLabelFont(fTxt);
1856  hOldHist1D -> GetYaxis() -> SetLabelSize(fLabH[1]);
1857  hOldHist1D -> GetYaxis() -> SetLabelOffset(fOffLH[1]);
1858  hOldHist1D -> GetYaxis() -> CenterTitle(fCnt);
1859  }
1860  }
1861  for (auto oldHistRow2D : m_vecOldHists2D) {
1862  for (auto oldHistTypes2D : oldHistRow2D) {
1863  for (auto hOldHist2D : oldHistTypes2D) {
1864  hOldHist2D -> SetMarkerColor(m_config.fCol.first);
1865  hOldHist2D -> SetMarkerStyle(m_config.fMar.first);
1866  hOldHist2D -> SetFillColor(m_config.fCol.first);
1867  hOldHist2D -> SetFillStyle(m_config.fFil.first);
1868  hOldHist2D -> SetLineColor(m_config.fCol.first);
1869  hOldHist2D -> SetLineStyle(m_config.fLin.first);
1870  hOldHist2D -> SetLineWidth(m_config.fWid.first);
1871  hOldHist2D -> SetTitle(m_config.legOld.data());
1872  hOldHist2D -> SetTitleFont(fTxt);
1873  hOldHist2D -> GetXaxis() -> SetTitleFont(fTxt);
1874  hOldHist2D -> GetXaxis() -> SetTitleSize(fTitH[0]);
1875  hOldHist2D -> GetXaxis() -> SetTitleOffset(fOffTH[0]);
1876  hOldHist2D -> GetXaxis() -> SetLabelFont(fTxt);
1877  hOldHist2D -> GetXaxis() -> SetLabelSize(fLabH[0]);
1878  hOldHist2D -> GetXaxis() -> SetLabelOffset(fOffLH[0]);
1879  hOldHist2D -> GetXaxis() -> CenterTitle(fCnt);
1880  hOldHist2D -> GetYaxis() -> SetTitleFont(fTxt);
1881  hOldHist2D -> GetYaxis() -> SetTitleSize(fTitH[1]);
1882  hOldHist2D -> GetYaxis() -> SetTitleOffset(fOffTH[1]);
1883  hOldHist2D -> GetYaxis() -> SetLabelFont(fTxt);
1884  hOldHist2D -> GetYaxis() -> SetLabelSize(fLabH[1]);
1885  hOldHist2D -> GetYaxis() -> SetLabelOffset(fOffLH[1]);
1886  hOldHist2D -> GetYaxis() -> CenterTitle(fCnt);
1887  hOldHist2D -> GetZaxis() -> SetTitle(countUse.data());
1888  hOldHist2D -> GetZaxis() -> SetTitleFont(fTxt);
1889  hOldHist2D -> GetZaxis() -> SetTitleSize(fTitH[2]);
1890  hOldHist2D -> GetZaxis() -> SetTitleOffset(fOffTH[2]);
1891  hOldHist2D -> GetZaxis() -> SetLabelFont(fTxt);
1892  hOldHist2D -> GetZaxis() -> SetLabelSize(fLabH[2]);
1893  hOldHist2D -> GetZaxis() -> SetLabelOffset(fOffLH[2]);
1894  hOldHist2D -> GetZaxis() -> CenterTitle(fCnt);
1895  }
1896  }
1897  }
1898 
1899  // set new histogram styles
1900  for (auto newHistRow1D : vecNewHists1D) {
1901  for (auto hNewHist1D : newHistRow1D) {
1902  hNewHist1D -> SetMarkerColor(m_config.fCol.second);
1903  hNewHist1D -> SetMarkerStyle(m_config.fMar.second);
1904  hNewHist1D -> SetFillColor(m_config.fCol.second);
1905  hNewHist1D -> SetFillStyle(m_config.fFil.second);
1906  hNewHist1D -> SetLineColor(m_config.fCol.second);
1907  hNewHist1D -> SetLineStyle(m_config.fLin.second);
1908  hNewHist1D -> SetLineWidth(m_config.fWid.second);
1909  hNewHist1D -> SetTitle("");
1910  hNewHist1D -> SetTitleFont(fTxt);
1911  hNewHist1D -> GetXaxis() -> SetTitleFont(fTxt);
1912  hNewHist1D -> GetXaxis() -> SetTitleSize(fTitH[0]);
1913  hNewHist1D -> GetXaxis() -> SetTitleOffset(fOffTH[0]);
1914  hNewHist1D -> GetXaxis() -> SetLabelFont(fTxt);
1915  hNewHist1D -> GetXaxis() -> SetLabelSize(fLabH[0]);
1916  hNewHist1D -> GetXaxis() -> SetLabelOffset(fOffLH[0]);
1917  hNewHist1D -> GetXaxis() -> CenterTitle(fCnt);
1918  hNewHist1D -> GetYaxis() -> SetTitle(countUse.data());
1919  hNewHist1D -> GetYaxis() -> SetTitleFont(fTxt);
1920  hNewHist1D -> GetYaxis() -> SetTitleSize(fTitH[1]);
1921  hNewHist1D -> GetYaxis() -> SetTitleOffset(fOffTH[1]);
1922  hNewHist1D -> GetYaxis() -> SetLabelFont(fTxt);
1923  hNewHist1D -> GetYaxis() -> SetLabelSize(fLabH[1]);
1924  hNewHist1D -> GetYaxis() -> SetLabelOffset(fOffLH[1]);
1925  hNewHist1D -> GetYaxis() -> CenterTitle(fCnt);
1926  }
1927  }
1928  for (auto newHistRow2D : vecNewHists2D) {
1929  for (auto newHistTypes2D : newHistRow2D) {
1930  for (auto hNewHist2D : newHistTypes2D) {
1931  hNewHist2D -> SetMarkerColor(m_config.fCol.first);
1932  hNewHist2D -> SetMarkerStyle(m_config.fMar.first);
1933  hNewHist2D -> SetFillColor(m_config.fCol.first);
1934  hNewHist2D -> SetFillStyle(m_config.fFil.first);
1935  hNewHist2D -> SetLineColor(m_config.fCol.first);
1936  hNewHist2D -> SetLineStyle(m_config.fLin.first);
1937  hNewHist2D -> SetLineWidth(m_config.fWid.first);
1938  hNewHist2D -> SetTitle(m_config.legNew.data());
1939  hNewHist2D -> SetTitleFont(fTxt);
1940  hNewHist2D -> GetXaxis() -> SetTitleFont(fTxt);
1941  hNewHist2D -> GetXaxis() -> SetTitleSize(fTitH[0]);
1942  hNewHist2D -> GetXaxis() -> SetTitleOffset(fOffTH[0]);
1943  hNewHist2D -> GetXaxis() -> SetLabelFont(fTxt);
1944  hNewHist2D -> GetXaxis() -> SetLabelSize(fLabH[0]);
1945  hNewHist2D -> GetXaxis() -> SetLabelOffset(fOffLH[0]);
1946  hNewHist2D -> GetXaxis() -> CenterTitle(fCnt);
1947  hNewHist2D -> GetYaxis() -> SetTitleFont(fTxt);
1948  hNewHist2D -> GetYaxis() -> SetTitleSize(fTitH[1]);
1949  hNewHist2D -> GetYaxis() -> SetTitleOffset(fOffTH[1]);
1950  hNewHist2D -> GetYaxis() -> SetLabelFont(fTxt);
1951  hNewHist2D -> GetYaxis() -> SetLabelSize(fLabH[1]);
1952  hNewHist2D -> GetYaxis() -> SetLabelOffset(fOffLH[1]);
1953  hNewHist2D -> GetYaxis() -> CenterTitle(fCnt);
1954  hNewHist2D -> GetZaxis() -> SetTitle(countUse.data());
1955  hNewHist2D -> GetZaxis() -> SetTitleFont(fTxt);
1956  hNewHist2D -> GetZaxis() -> SetTitleSize(fTitH[2]);
1957  hNewHist2D -> GetZaxis() -> SetTitleOffset(fOffTH[2]);
1958  hNewHist2D -> GetZaxis() -> SetLabelFont(fTxt);
1959  hNewHist2D -> GetZaxis() -> SetLabelSize(fLabH[2]);
1960  hNewHist2D -> GetZaxis() -> SetLabelOffset(fOffLH[2]);
1961  hNewHist2D -> GetZaxis() -> CenterTitle(fCnt);
1962  }
1963  }
1964  }
1965 
1966  // set ratio styles
1967  for (auto ratioRow1D : vecRatios1D) {
1968  for (auto hRatio1D : ratioRow1D) {
1969  hRatio1D -> SetMarkerColor(m_config.fCol.first);
1970  hRatio1D -> SetMarkerStyle(m_config.fMar.first);
1971  hRatio1D -> SetFillColor(m_config.fCol.first);
1972  hRatio1D -> SetFillStyle(m_config.fFil.first);
1973  hRatio1D -> SetLineColor(m_config.fCol.first);
1974  hRatio1D -> SetLineStyle(m_config.fLin.first);
1975  hRatio1D -> SetLineWidth(m_config.fWid.first);
1976  hRatio1D -> SetTitle("");
1977  hRatio1D -> SetTitleFont(fTxt);
1978  hRatio1D -> GetXaxis() -> SetTitleFont(fTxt);
1979  hRatio1D -> GetXaxis() -> SetTitleSize(fTitR1[0]);
1980  hRatio1D -> GetXaxis() -> SetTitleOffset(fOffTR1[0]);
1981  hRatio1D -> GetXaxis() -> SetLabelFont(fTxt);
1982  hRatio1D -> GetXaxis() -> SetLabelSize(fLabR1[0]);
1983  hRatio1D -> GetXaxis() -> SetLabelOffset(fOffLR1[0]);
1984  hRatio1D -> GetXaxis() -> CenterTitle(fCnt);
1985  hRatio1D -> GetYaxis() -> SetTitle(m_config.ratio.data());
1986  hRatio1D -> GetYaxis() -> SetTitleFont(fTxt);
1987  hRatio1D -> GetYaxis() -> SetTitleSize(fTitR1[1]);
1988  hRatio1D -> GetYaxis() -> SetTitleOffset(fOffTR1[1]);
1989  hRatio1D -> GetYaxis() -> SetLabelFont(fTxt);
1990  hRatio1D -> GetYaxis() -> SetLabelSize(fLabR1[1]);
1991  hRatio1D -> GetYaxis() -> SetLabelOffset(fOffLR1[1]);
1992  hRatio1D -> GetYaxis() -> CenterTitle(fCnt);
1993  }
1994  }
1995  for (auto ratioRow2D : vecRatios2D) {
1996  for (auto ratioTypes2D : ratioRow2D) {
1997  for (auto hRatio2D : ratioTypes2D) {
1998  hRatio2D -> SetMarkerColor(m_config.fCol.first);
1999  hRatio2D -> SetMarkerStyle(m_config.fMar.first);
2000  hRatio2D -> SetFillColor(m_config.fCol.first);
2001  hRatio2D -> SetFillStyle(m_config.fFil.first);
2002  hRatio2D -> SetLineColor(m_config.fCol.first);
2003  hRatio2D -> SetLineStyle(m_config.fLin.first);
2004  hRatio2D -> SetLineWidth(m_config.fWid.first);
2005  hRatio2D -> SetTitle("");
2006  hRatio2D -> SetTitleFont(fTxt);
2007  hRatio2D -> GetXaxis() -> SetTitleFont(fTxt);
2008  hRatio2D -> GetXaxis() -> SetTitleSize(fTitR2[0]);
2009  hRatio2D -> GetXaxis() -> SetTitleOffset(fOffTR2[0]);
2010  hRatio2D -> GetXaxis() -> SetLabelFont(fTxt);
2011  hRatio2D -> GetXaxis() -> SetLabelSize(fLabR2[0]);
2012  hRatio2D -> GetXaxis() -> SetLabelOffset(fOffLR2[0]);
2013  hRatio2D -> GetXaxis() -> CenterTitle(fCnt);
2014  hRatio2D -> GetYaxis() -> SetTitleFont(fTxt);
2015  hRatio2D -> GetYaxis() -> SetTitleSize(fTitR2[1]);
2016  hRatio2D -> GetYaxis() -> SetTitleOffset(fOffTR2[1]);
2017  hRatio2D -> GetYaxis() -> SetLabelFont(fTxt);
2018  hRatio2D -> GetYaxis() -> SetLabelSize(fLabR2[1]);
2019  hRatio2D -> GetYaxis() -> SetLabelOffset(fOffLR2[1]);
2020  hRatio2D -> GetYaxis() -> CenterTitle(fCnt);
2021  hRatio2D -> GetZaxis() -> SetTitle(m_config.ratio.data());
2022  hRatio2D -> GetZaxis() -> SetTitleFont(fTxt);
2023  hRatio2D -> GetZaxis() -> SetTitleSize(fTitR2[2]);
2024  hRatio2D -> GetZaxis() -> SetTitleOffset(fOffTR2[2]);
2025  hRatio2D -> GetZaxis() -> SetLabelFont(fTxt);
2026  hRatio2D -> GetZaxis() -> SetLabelSize(fLabR2[2]);
2027  hRatio2D -> GetZaxis() -> SetLabelOffset(fOffLR2[2]);
2028  hRatio2D -> GetZaxis() -> CenterTitle(fCnt);
2029  }
2030  }
2031  }
2032  cout << " Set styles." << endl;
2033 
2034  // make legends and text boxes ----------------------------------------------
2035 
2036  // make legend
2037  const uint32_t fColLe(0);
2038  const uint32_t fFilLe(0);
2039  const uint32_t fLinLe(0);
2040  const float fLegXY[m_const.nVtx] = {0.1, 0.1, 0.3, 0.2};
2041 
2042  TLegend *leg = new TLegend(fLegXY[0], fLegXY[1], fLegXY[2], fLegXY[3], m_config.header.data());
2043  leg -> SetFillColor(fColLe);
2044  leg -> SetFillStyle(fFilLe);
2045  leg -> SetLineColor(fColLe);
2046  leg -> SetLineStyle(fLinLe);
2047  leg -> SetTextFont(fTxt);
2048  leg -> SetTextAlign(fAln);
2049  leg -> AddEntry(m_vecOldHists1D.front().front(), m_config.legOld.data(), "pf");
2050  leg -> AddEntry(vecNewHists1D.front().front(), m_config.legNew.data(), "pf");
2051  cout << " Made legend." << endl;
2052 
2053  // make text
2054  const uint32_t fColTx(0);
2055  const uint32_t fFilTx(0);
2056  const uint32_t fLinTx(0);
2057  const float fTxtXY[m_const.nVtx] = {0.3, 0.1, 0.5, 0.25};
2058 
2059  TPaveText *txt = new TPaveText(fTxtXY[0], fTxtXY[1], fTxtXY[2], fTxtXY[3], "NDC NB");
2060  txt -> SetFillColor(fColTx);
2061  txt -> SetFillStyle(fFilTx);
2062  txt -> SetLineColor(fColTx);
2063  txt -> SetLineStyle(fLinTx);
2064  txt -> SetTextFont(fTxt);
2065  txt -> SetTextAlign(fAln);
2066  for (string txtLine : m_config.info) {
2067  txt -> AddText(txtLine.data());
2068  }
2069  cout << " Made text." << endl;
2070 
2071  // make plots ---------------------------------------------------------------
2072 
2073  // canvas parameters
2074  const uint32_t width1D(750);
2075  const uint32_t width2D(1500);
2076  const uint32_t width2DR(2250);
2077  const uint32_t height(750);
2078  const uint32_t heightR1(900);
2079  const uint32_t heightR2(500);
2080  const uint32_t fMode(0);
2081  const uint32_t fBord(2);
2082  const uint32_t fGrid(0);
2083  const uint32_t fTick(1);
2084  const uint32_t fLogX(0);
2085  const uint32_t fLogY(1);
2086  const uint32_t fLogZ(1);
2087  const uint32_t fFrame(0);
2088  const string sOldPadName("pOld");
2089  const string sNewPadName("pNew");
2090  const string sHistPadName("pHist");
2091  const string sRatPadName("pRatio");
2092  const float fMargin1D[m_const.nSide] = {0.02, 0.02, 0.15, 0.15};
2093  const float fMargin1DH[m_const.nSide] = {0.02, 0.02, 0.005, 0.15};
2094  const float fMargin1DR[m_const.nSide] = {0.005, 0.02, 0.2, 0.15};
2095  const float fMargin2D[m_const.nSide] = {0.10, 0.15, 0.15, 0.15};
2096  const float xyOldPad[m_const.nVtx] = {0.0, 0.0, 0.5, 1.0};
2097  const float xyOldPadR[m_const.nVtx] = {0.0, 0.0, 0.33, 1.0};
2098  const float xyNewPad[m_const.nVtx] = {0.5, 0.0, 1.0, 1.0};
2099  const float xyNewPadR[m_const.nVtx] = {0.33, 0.0, 0.66, 1.0};
2100  const float xyHistPadR[m_const.nVtx] = {0.0, 0.33, 1.0, 1.0};
2101  const float xyRatPadR1[m_const.nVtx] = {0.0, 0.0, 1.0, 0.33};
2102  const float xyRatPadR2[m_const.nVtx] = {0.66, 0.0, 1.0, 1.0};
2103 
2104  // make 1d plot without ratio
2105  for (size_t iHistRow1D = 0; iHistRow1D < nHistRows1D; iHistRow1D++) {
2106  for (size_t iHist1D = 0; iHist1D < nHistTypes1D; iHist1D++) {
2107 
2108  // make new name
2109  const string sName = vecCanvasNames1D[iHistRow1D][iHist1D] + "_" + sLabel;
2110 
2111 
2112  // construct canvas
2113  TCanvas* cPlot1D = new TCanvas(sName.data(), "", width1D, height);
2114  cPlot1D -> SetGrid(fGrid, fGrid);
2115  cPlot1D -> SetTicks(fTick, fTick);
2116  cPlot1D -> SetBorderMode(fMode);
2117  cPlot1D -> SetBorderSize(fBord);
2118  cPlot1D -> SetFrameBorderMode(fFrame);
2119  cPlot1D -> SetTopMargin(fMargin1D[0]);
2120  cPlot1D -> SetRightMargin(fMargin1D[1]);
2121  cPlot1D -> SetBottomMargin(fMargin1D[2]);
2122  cPlot1D -> SetLeftMargin(fMargin1D[3]);
2123  cPlot1D -> SetLogx(fLogX);
2124  cPlot1D -> SetLogy(fLogY);
2125  cPlot1D -> cd();
2126 
2127  // draw old vs. new histograms
2128  m_vecOldHists1D[iHistRow1D][iHist1D] -> Draw();
2129  vecNewHists1D[iHistRow1D][iHist1D] -> Draw("same");
2130 
2131  // draw text and save
2132  leg -> Draw();
2133  txt -> Draw();
2134  m_vecPlotDirs.at(iDir) -> cd();
2135  cPlot1D -> Write();
2136  cPlot1D -> Close();
2137  }
2138  }
2139 
2140  // make 1d plot with ratio
2141  for (size_t iHistRow1D = 0; iHistRow1D < nHistRows1D; iHistRow1D++) {
2142  for (size_t iHist1D = 0; iHist1D < nHistTypes1D; iHist1D++) {
2143 
2144  // make new name
2145  const string sNameWithRatio = vecCanvasNames1D[iHistRow1D][iHist1D] + "WithRatio_" + sLabel;
2146 
2147  // construct canvas
2148  TCanvas* cPlot1D = new TCanvas(sNameWithRatio.data(), "", width1D, heightR1);
2149  TPad* pPadH1D = new TPad(sHistPadName.data(), "", xyHistPadR[0], xyHistPadR[1], xyHistPadR[2], xyHistPadR[3]);
2150  TPad* pPadR1D = new TPad(sRatPadName.data(), "", xyRatPadR1[0], xyRatPadR1[1], xyRatPadR1[2], xyRatPadR1[3]);
2151  cPlot1D -> SetGrid(fGrid, fGrid);
2152  cPlot1D -> SetTicks(fTick, fTick);
2153  cPlot1D -> SetBorderMode(fMode);
2154  cPlot1D -> SetBorderSize(fBord);
2155  pPadH1D -> SetGrid(fGrid, fGrid);
2156  pPadH1D -> SetTicks(fTick, fTick);
2157  pPadH1D -> SetBorderMode(fMode);
2158  pPadH1D -> SetBorderSize(fBord);
2159  pPadH1D -> SetFrameBorderMode(fFrame);
2160  pPadH1D -> SetTopMargin(fMargin1DH[0]);
2161  pPadH1D -> SetRightMargin(fMargin1DH[1]);
2162  pPadH1D -> SetBottomMargin(fMargin1DH[2]);
2163  pPadH1D -> SetLeftMargin(fMargin1DH[3]);
2164  pPadH1D -> SetLogx(fLogX);
2165  pPadH1D -> SetLogy(fLogY);
2166  pPadR1D -> SetGrid(fGrid, fGrid);
2167  pPadR1D -> SetTicks(fTick, fTick);
2168  pPadR1D -> SetBorderMode(fMode);
2169  pPadR1D -> SetBorderSize(fBord);
2170  pPadR1D -> SetFrameBorderMode(fFrame);
2171  pPadR1D -> SetTopMargin(fMargin1DR[0]);
2172  pPadR1D -> SetRightMargin(fMargin1DR[1]);
2173  pPadR1D -> SetBottomMargin(fMargin1DR[2]);
2174  pPadR1D -> SetLeftMargin(fMargin1DR[3]);
2175  pPadR1D -> SetLogx(fLogX);
2176  pPadR1D -> SetLogy(fLogY);
2177  cPlot1D -> cd();
2178  pPadH1D -> Draw();
2179  pPadR1D -> Draw();
2180 
2181  // draw old vs. new histograms
2182  pPadH1D -> cd();
2183  m_vecOldHists1D[iHistRow1D][iHist1D] -> Draw();
2184  vecNewHists1D[iHistRow1D][iHist1D] -> Draw("same");
2185  leg -> Draw();
2186  txt -> Draw();
2187 
2188  // draw ratio
2189  pPadR1D -> cd();
2190  vecRatios1D[iHistRow1D][iHist1D] -> Draw();
2191 
2192  // draw text and save
2193  m_vecPlotDirs.at(iDir) -> cd();
2194  cPlot1D -> Write();
2195  cPlot1D -> Close();
2196  }
2197  }
2198 
2199  // make 2d plot without ratio
2200  for (size_t iHistRow2D = 0; iHistRow2D < nHistRows2D; iHistRow2D++) {
2201  for (size_t iHistType2D = 0; iHistType2D < nHistTypes2D; iHistType2D++) {
2202  for (size_t iHist2D = 0; iHist2D < nHist2D; iHist2D++) {
2203 
2204  // construct canvas
2205  TCanvas* cPlot2D = new TCanvas(vecCanvasNames2D[iHistRow2D][iHistType2D][iHist2D].data(), "", width2D, height);
2206  TPad* pOld = new TPad(sOldPadName.data(), "", xyOldPad[0], xyOldPad[1], xyOldPad[2], xyOldPad[3]);
2207  TPad* pNew = new TPad(sNewPadName.data(), "", xyNewPad[0], xyNewPad[1], xyNewPad[2], xyNewPad[3]);
2208  cPlot2D -> SetGrid(fGrid, fGrid);
2209  cPlot2D -> SetTicks(fTick, fTick);
2210  cPlot2D -> SetBorderMode(fMode);
2211  cPlot2D -> SetBorderSize(fBord);
2212  pOld -> SetGrid(fGrid, fGrid);
2213  pOld -> SetTicks(fTick, fTick);
2214  pOld -> SetBorderMode(fMode);
2215  pOld -> SetBorderSize(fBord);
2216  pOld -> SetFrameBorderMode(fFrame);
2217  pOld -> SetTopMargin(fMargin2D[0]);
2218  pOld -> SetRightMargin(fMargin2D[1]);
2219  pOld -> SetBottomMargin(fMargin2D[2]);
2220  pOld -> SetLeftMargin(fMargin2D[3]);
2221  pOld -> SetLogx(fLogX);
2222  pOld -> SetLogy(fLogY);
2223  pOld -> SetLogz(fLogZ);
2224  pNew -> SetGrid(fGrid, fGrid);
2225  pNew -> SetTicks(fTick, fTick);
2226  pNew -> SetBorderMode(fMode);
2227  pNew -> SetBorderSize(fBord);
2228  pNew -> SetFrameBorderMode(fFrame);
2229  pNew -> SetTopMargin(fMargin2D[0]);
2230  pNew -> SetRightMargin(fMargin2D[1]);
2231  pNew -> SetBottomMargin(fMargin2D[2]);
2232  pNew -> SetLeftMargin(fMargin2D[3]);
2233  pNew -> SetLogx(fLogX);
2234  pNew -> SetLogy(fLogY);
2235  pNew -> SetLogz(fLogZ);
2236  cPlot2D -> cd();
2237  pOld -> Draw();
2238  pNew -> Draw();
2239 
2240  // draw old vs. new histograms
2241  pOld -> cd();
2242  m_vecOldHists2D[iHistRow2D][iHistType2D][iHist2D] -> Draw("colz");
2243  pNew -> cd();
2244  vecNewHists2D[iHistRow2D][iHistType2D][iHist2D] -> Draw("colz");
2245 
2246  // draw text and save
2247  pNew -> cd();
2248  txt -> Draw();
2249  m_vecPlotDirs.at(iDir) -> cd();
2250  cPlot2D -> Write();
2251  cPlot2D -> Close();
2252  }
2253  }
2254  }
2255 
2256  // make 2d plot with ratio
2257  for (size_t iHistRow2D = 0; iHistRow2D < nHistRows2D; iHistRow2D++) {
2258  for (size_t iHistType2D = 0; iHistType2D < nHistTypes2D; iHistType2D++) {
2259  for (size_t iHist2D = 0; iHist2D < nHist2D; iHist2D++) {
2260 
2261  // make new name
2262  const string sNameWithRatio = vecCanvasNames2D[iHistRow2D][iHistType2D][iHist2D] + "_" + sLabel;
2263 
2264  // construct canvas
2265  TCanvas* cPlot2D = new TCanvas(sNameWithRatio.data(), "", width2DR, heightR2);
2266  TPad* pOld = new TPad(sOldPadName.data(), "", xyOldPadR[0], xyOldPadR[1], xyOldPadR[2], xyOldPadR[3]);
2267  TPad* pNew = new TPad(sNewPadName.data(), "", xyNewPadR[0], xyNewPadR[1], xyNewPadR[2], xyNewPadR[3]);
2268  TPad* pRat = new TPad(sRatPadName.data(), "", xyRatPadR2[0], xyRatPadR2[1], xyRatPadR2[2], xyRatPadR2[3]);
2269  cPlot2D -> SetGrid(fGrid, fGrid);
2270  cPlot2D -> SetTicks(fTick, fTick);
2271  cPlot2D -> SetBorderMode(fMode);
2272  cPlot2D -> SetBorderSize(fBord);
2273  pOld -> SetGrid(fGrid, fGrid);
2274  pOld -> SetTicks(fTick, fTick);
2275  pOld -> SetBorderMode(fMode);
2276  pOld -> SetBorderSize(fBord);
2277  pOld -> SetFrameBorderMode(fFrame);
2278  pOld -> SetTopMargin(fMargin2D[0]);
2279  pOld -> SetRightMargin(fMargin2D[1]);
2280  pOld -> SetBottomMargin(fMargin2D[2]);
2281  pOld -> SetLeftMargin(fMargin2D[3]);
2282  pOld -> SetLogx(fLogX);
2283  pOld -> SetLogy(fLogY);
2284  pOld -> SetLogz(fLogZ);
2285  pNew -> SetGrid(fGrid, fGrid);
2286  pNew -> SetTicks(fTick, fTick);
2287  pNew -> SetBorderMode(fMode);
2288  pNew -> SetBorderSize(fBord);
2289  pNew -> SetFrameBorderMode(fFrame);
2290  pNew -> SetTopMargin(fMargin2D[0]);
2291  pNew -> SetRightMargin(fMargin2D[1]);
2292  pNew -> SetBottomMargin(fMargin2D[2]);
2293  pNew -> SetLeftMargin(fMargin2D[3]);
2294  pNew -> SetLogx(fLogX);
2295  pNew -> SetLogy(fLogY);
2296  pNew -> SetLogz(fLogZ);
2297  pRat -> SetGrid(fGrid, fGrid);
2298  pRat -> SetTicks(fTick, fTick);
2299  pRat -> SetBorderMode(fMode);
2300  pRat -> SetBorderSize(fBord);
2301  pRat -> SetFrameBorderMode(fFrame);
2302  pRat -> SetTopMargin(fMargin2D[0]);
2303  pRat -> SetRightMargin(fMargin2D[1]);
2304  pRat -> SetBottomMargin(fMargin2D[2]);
2305  pRat -> SetLeftMargin(fMargin2D[3]);
2306  pRat -> SetLogx(fLogX);
2307  pRat -> SetLogy(fLogY);
2308  pRat -> SetLogz(fLogZ);
2309  cPlot2D -> cd();
2310  pOld -> Draw();
2311  pNew -> Draw();
2312  pRat -> Draw();
2313 
2314  // draw old vs. new vs. ratio histograms
2315  pOld -> cd();
2316  m_vecOldHists2D[iHistRow2D][iHistType2D][iHist2D] -> Draw("colz");
2317  txt -> Draw();
2318  pNew -> cd();
2319  vecNewHists2D[iHistRow2D][iHistType2D][iHist2D] -> Draw("colz");
2320  pRat -> cd();
2321  vecRatios2D[iHistRow2D][iHistType2D][iHist2D] -> Draw("colz");
2322 
2323  // draw text and save
2324  m_vecPlotDirs.at(iDir) -> cd();
2325  cPlot2D -> Write();
2326  cPlot2D -> Close();
2327  }
2328  }
2329  }
2330  cout << " Made 2D plot with ratio." << endl;
2331 
2332  // save ratio histograms
2333  m_vecRatioDirs.at(iDir) -> cd();
2334  for (auto histRow1D : vecRatios1D) {
2335  for (auto hRatio1D : histRow1D) {
2336  hRatio1D -> Write();
2337  }
2338  }
2339  for (auto histRow2D : vecRatios2D) {
2340  for (auto histTypes2D : histRow2D) {
2341  for (auto hRatio2D : histTypes2D) {
2342  hRatio2D -> Write();
2343  }
2344  }
2345  }
2346 
2347  // announce end and return
2348  cout << " Saved ratio histograms.\n"
2349  << " Finished making ratios and plots."
2350  << endl;
2351  return;
2352 
2353 } // end 'MakeRatiosAndPlots(vector<TH1D*>, vector<TH2D*>, int, string)'
2354 
2355 
2356 
2357 void STrackMatcherComparator::SaveHistograms() {
2358 
2359  // save 1d histograms
2360  for (size_t iHistRow = 0; iHistRow < m_vecTreeHists1D.size(); iHistRow++) {
2361  for (size_t iHistType = 0; iHistType < m_vecTreeHists1D[iHistRow].size(); iHistType++) {
2362  m_vecHistDirs[Src::NewTree] -> cd();
2363  m_vecTreeHists1D[iHistRow][iHistType] -> Write();
2364  m_vecHistDirs[Src::NewTuple] -> cd();
2365  m_vecTupleHists1D[iHistRow][iHistType] -> Write();
2366  m_vecHistDirs[Src::OldTuple] -> cd();
2367  m_vecOldHists1D[iHistRow][iHistType] -> Write();
2368  } // end type loop
2369  } // end row loop
2370  cout << " Saved 1d histograms." << endl;
2371 
2372  // save 2d histograms
2373  for (size_t iHistRow = 0; iHistRow < m_vecTreeHists2D.size(); iHistRow++) {
2374  for (size_t iHistType = 0; iHistType < m_vecTreeHists2D[iHistRow].size(); iHistType++) {
2375  for (size_t iVsHistMod = 0; iVsHistMod < m_vecTreeHists2D[iHistRow][iHistType].size(); iVsHistMod++) {
2376  m_vecHistDirs[Src::NewTree] -> cd();
2377  m_vecTreeHists2D[iHistRow][iHistType][iVsHistMod] -> Write();
2378  m_vecHistDirs[Src::NewTuple] -> cd();
2379  m_vecTupleHists2D[iHistRow][iHistType][iVsHistMod] -> Write();
2380  m_vecHistDirs[Src::OldTuple] -> cd();
2381  m_vecOldHists2D[iHistRow][iHistType][iVsHistMod] -> Write();
2382  } // end modifier loop
2383  } // end type loop
2384  } // end row loop
2385 
2386  // announce end and return
2387  cout << " Saved 2d histograms." << endl;
2388  return;
2389 
2390 } // end 'SaveHistograms()'
2391 
2392 
2393 
2394 void STrackMatcherComparator::CloseInput() {
2395 
2396  m_treeInFileTrue -> cd();
2397  m_treeInFileTrue -> Close();
2398  m_treeInFileReco -> cd();
2399  m_treeInFileReco -> Close();
2400  m_tupleInFileTrue -> cd();
2401  m_tupleInFileTrue -> Close();
2402  m_tupleInFileReco -> cd();
2403  m_tupleInFileReco -> Close();
2404  m_oldInFileTrue -> cd();
2405  m_oldInFileTrue -> Close();
2406  m_oldInFileReco -> cd();
2407  m_oldInFileReco -> Close();
2408 
2409  // announce end and return
2410  cout << " Closed input files." << endl;
2411  return;
2412 
2413 } // end 'CloseInput()'
2414 
2415 
2416 
2417 void STrackMatcherComparator::CloseOutput() {
2418 
2419  m_outFile -> cd();
2420  m_outFile -> Close();
2421 
2422  // announce end and return
2423  cout << " Closed output file." << endl;
2424  return;
2425 
2426 } // end 'CloseOutput()'
2427 
2428 
2429 
2430 // helper functions -----------------------------------------------------------
2431 
2432 void STrackMatcherComparator::FillHistogram1D(
2433  const STrackMatcherComparatorHistContent& content,
2434  const Type type,
2435  const vector<vector<TH1D*>> vecHist1D
2436 ) {
2437 
2438  vecHist1D.at(Var::NTot).at(type) -> Fill(content.nTot);
2439  vecHist1D.at(Var::NIntt).at(type) -> Fill(content.nIntt);
2440  vecHist1D.at(Var::NMvtx).at(type) -> Fill(content.nMvtx);
2441  vecHist1D.at(Var::NTpc).at(type) -> Fill(content.nTpc);
2442  vecHist1D.at(Var::RTot).at(type) -> Fill(content.rTot);
2443  vecHist1D.at(Var::RIntt).at(type) -> Fill(content.rIntt);
2444  vecHist1D.at(Var::RMvtx).at(type) -> Fill(content.rMvtx);
2445  vecHist1D.at(Var::RTpc).at(type) -> Fill(content.rTpc);
2446  vecHist1D.at(Var::Phi).at(type) -> Fill(content.phi);
2447  vecHist1D.at(Var::Eta).at(type) -> Fill(content.eta);
2448  vecHist1D.at(Var::Pt).at(type) -> Fill(content.pt);
2449  vecHist1D.at(Var::Frac).at(type) -> Fill(content.ptFrac);
2450  vecHist1D.at(Var::Qual).at(type) -> Fill(content.quality);
2451  vecHist1D.at(Var::PtErr).at(type) -> Fill(content.ptErr);
2452  vecHist1D.at(Var::EtaErr).at(type) -> Fill(content.etaErr);
2453  vecHist1D.at(Var::PhiErr).at(type) -> Fill(content.phiErr);
2454  vecHist1D.at(Var::PtRes).at(type) -> Fill(content.ptRes);
2455  vecHist1D.at(Var::EtaRes).at(type) -> Fill(content.etaRes);
2456  vecHist1D.at(Var::PhiRes).at(type) -> Fill(content.phiRes);
2457  return;
2458 
2459 } // end 'FillHistogram1D(STrackMatcherComparatorHistContent&, int, vector<vector<TH1D*>>)'
2460 
2461 
2462 
2463 
2464 void STrackMatcherComparator::FillHistogram2D(
2465  const STrackMatcherComparatorHistContent& content,
2466  const Type type,
2467  const Comp comparison,
2468  const double value,
2469  const vector<vector<vector<TH2D*>>> vecHist2D
2470 
2471 ) {
2472 
2473  vecHist2D.at(Var::NTot).at(type).at(comparison) -> Fill(value, content.nTot);
2474  vecHist2D.at(Var::NIntt).at(type).at(comparison) -> Fill(value, content.nIntt);
2475  vecHist2D.at(Var::NMvtx).at(type).at(comparison) -> Fill(value, content.nMvtx);
2476  vecHist2D.at(Var::NTpc).at(type).at(comparison) -> Fill(value, content.nTpc);
2477  vecHist2D.at(Var::RTot).at(type).at(comparison) -> Fill(value, content.rTot);
2478  vecHist2D.at(Var::RIntt).at(type).at(comparison) -> Fill(value, content.rIntt);
2479  vecHist2D.at(Var::RMvtx).at(type).at(comparison) -> Fill(value, content.rMvtx);
2480  vecHist2D.at(Var::RTpc).at(type).at(comparison) -> Fill(value, content.rTpc);
2481  vecHist2D.at(Var::Phi).at(type).at(comparison) -> Fill(value, content.phi);
2482  vecHist2D.at(Var::Eta).at(type).at(comparison) -> Fill(value, content.eta);
2483  vecHist2D.at(Var::Pt).at(type).at(comparison) -> Fill(value, content.pt);
2484  vecHist2D.at(Var::Frac).at(type).at(comparison) -> Fill(value, content.ptFrac);
2485  vecHist2D.at(Var::Qual).at(type).at(comparison) -> Fill(value, content.quality);
2486  vecHist2D.at(Var::PtErr).at(type).at(comparison) -> Fill(value, content.ptErr);
2487  vecHist2D.at(Var::EtaErr).at(type).at(comparison) -> Fill(value, content.etaErr);
2488  vecHist2D.at(Var::PhiErr).at(type).at(comparison) -> Fill(value, content.phiErr);
2489  vecHist2D.at(Var::PtRes).at(type).at(comparison) -> Fill(value, content.ptRes);
2490  vecHist2D.at(Var::EtaRes).at(type).at(comparison) -> Fill(value, content.etaRes);
2491  vecHist2D.at(Var::PhiRes).at(type).at(comparison) -> Fill(value, content.phiRes);
2492  return;
2493 
2494 } // end 'FillVsHistogram2D(STrackMatcherComparator&, int, int, double, vector<vector<vector<TH2D*>>>)'
2495 
2496 
2497 
2498 bool STrackMatcherComparator::IsNearSector(const float phi) {
2499 
2500  bool isNearSector = false;
2501  for (size_t iSector = 0; iSector < m_const.nSectors; iSector++) {
2502  const float cutVal = m_config.sigCutVal * m_config.phiSectors[iSector].second;
2503  const float minPhi = m_config.phiSectors[iSector].first - cutVal;
2504  const float maxPhi = m_config.phiSectors[iSector].first + cutVal;
2505  const bool isNear = ((phi >= minPhi) && (phi <= maxPhi));
2506  if (isNear) {
2507  isNearSector = true;
2508  break;
2509  }
2510  } // end sector loop
2511  return isNearSector;
2512 
2513 } // end 'IsNearSector(float)'
2514 
2515 // end ------------------------------------------------------------------------