Fork me on GitHub

source: git/validation/DelphesValidation.cpp@ 1388e73

ImprovedOutputFile
Last change on this file since 1388e73 was 341014c, checked in by Pavel Demin <pavel-demin@…>, 6 years ago

apply .clang-format to all .h, .cc and .cpp files

  • Property mode set to 100644
File size: 107.5 KB
RevLine 
[734b267]1/*
2 * Delphes: a framework for fast simulation of a generic collider experiment
3 * Copyright (C) 2012-2014 Universite catholique de Louvain (UCL), Belgium
4 *
5 * This program is free software: you can redistribute it and/or modify
6 * it under the terms of the GNU General Public License as published by
7 * the Free Software Foundation, either version 3 of the License, or
8 * (at your option) any later version.
9 *
10 * This program is distributed in the hope that it will be useful,
11 * but WITHOUT ANY WARRANTY; without even the implied warranty of
12 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13 * GNU General Public License for more details.
14 *
15 * You should have received a copy of the GNU General Public License
16 * along with this program. If not, see <http://www.gnu.org/licenses/>.
17 */
18
19#include <iostream>
[341014c]20#include <typeinfo>
[734b267]21#include <utility>
22#include <vector>
23
[341014c]24#include "TApplication.h"
[734b267]25#include "TROOT.h"
26#include "TSystem.h"
27
28#include "TString.h"
29
[341014c]30#include "TCanvas.h"
31#include "TClonesArray.h"
32#include "TGraph.h"
33#include "TGraphErrors.h"
[1408174]34#include "TH1.h"
[734b267]35#include "TH2.h"
36#include "THStack.h"
37#include "TLegend.h"
38#include "TLorentzVector.h"
[341014c]39#include "TMath.h"
[1408174]40#include "TMultiGraph.h"
[341014c]41#include "TPaveText.h"
42#include "TStyle.h"
[734b267]43
44#include "classes/DelphesClasses.h"
45
[341014c]46#include "ExRootAnalysis/ExRootResult.h"
47#include "ExRootAnalysis/ExRootTreeBranch.h"
[734b267]48#include "ExRootAnalysis/ExRootTreeReader.h"
49#include "ExRootAnalysis/ExRootTreeWriter.h"
50#include "ExRootAnalysis/ExRootUtilities.h"
51
52using namespace std;
53
54//------------------------------------------------------------------------------
55
[92237ca]56static const int Nbins = 50;
[1408174]57
58int objStyle = 1;
59int trackStyle = 7;
60int towerStyle = 3;
61
62Color_t objColor = kBlack;
63Color_t trackColor = kBlack;
64Color_t towerColor = kBlack;
65
66double effLegXmin = 0.22;
67double effLegXmax = 0.7;
68double effLegYmin = 0.22;
69double effLegYmax = 0.5;
70
71double resLegXmin = 0.62;
72double resLegXmax = 0.9;
73double resLegYmin = 0.52;
74double resLegYmax = 0.85;
75
76double topLeftLegXmin = 0.22;
77double topLeftLegXmax = 0.7;
78double topLeftLegYmin = 0.52;
79double topLeftLegYmax = 0.85;
80
[92237ca]81unsigned int k;
82
[1408174]83struct resolPlot
84{
[341014c]85 TH1 *resolHist;
[1408174]86 double ptmin;
87 double ptmax;
[92237ca]88 double etamin;
89 double etamax;
[1408174]90 double xmin;
91 double xmax;
92 TString obj;
93
94 resolPlot();
95 resolPlot(double ptdown, double ptup, TString object);
[92237ca]96 resolPlot(double etadown, double etaup, double ptdown, double ptup, TString object);
[1408174]97 void set(double ptdown, double ptup, TString object, double xmin = 0, double xmax = 2);
[92237ca]98 void set(double etadown, double etaup, double ptdown, double ptup, TString object, double xmin = 0, double xmax = 2);
[341014c]99 void print() { std::cout << ptmin << std::endl; }
[1408174]100};
101
102resolPlot::resolPlot()
103{
104}
105
106resolPlot::resolPlot(double ptdown, double ptup, TString object)
107{
[341014c]108 this->set(ptdown, ptup, object);
[1408174]109}
110
[92237ca]111resolPlot::resolPlot(double etadown, double etaup, double ptdown, double ptup, TString object)
112{
113 this->set(etadown, etaup, ptdown, ptup, object);
114}
115
[1408174]116void resolPlot::set(double ptdown, double ptup, TString object, double xmin, double xmax)
117{
118 ptmin = ptdown;
119 ptmax = ptup;
120 obj = object;
121
[341014c]122 resolHist = new TH1D(obj + "_delta_pt_" + Form("%4.2f", ptmin) + "_" + Form("%4.2f", ptmax), obj + "_delta_pt_" + Form("%4.2f", ptmin) + "_" + Form("%4.2f", ptmax), 1000, xmin, xmax);
[92237ca]123}
124
125void resolPlot::set(double etadown, double etaup, double ptdown, double ptup, TString object, double xmin, double xmax)
126{
127 etamin = etadown;
128 etamax = etaup;
129 ptmin = ptdown;
130 ptmax = ptup;
131 obj = object;
132
[341014c]133 resolHist = new TH1D(obj + "_delta_pt_" + Form("%4.2f", ptmin) + "_" + Form("%4.2f", ptmax) + "_" + Form("%4.2f", etamin) + "_" + Form("%4.2f", etamax), obj + "_delta_pt_" + Form("%4.2f", ptmin) + "_" + Form("%4.2f", ptmax) + "_" + Form("%4.2f", etamin) + "_" + Form("%4.2f", etamax), 1000, xmin, xmax);
[1408174]134}
135
136void HistogramsCollection(std::vector<resolPlot> *histos, double ptmin, double ptmax, TString obj, double xmin = 0, double xmax = 2)
137{
138 double width;
139 double ptdown;
140 double ptup;
141 resolPlot ptemp;
142
[341014c]143 for(int i = 0; i < Nbins; i++)
[1408174]144 {
145 width = (ptmax - ptmin) / Nbins;
[341014c]146 ptdown = TMath::Power(10, ptmin + i * width);
147 ptup = TMath::Power(10, ptmin + (i + 1) * width);
[1408174]148 ptemp.set(ptdown, ptup, obj, xmin, xmax);
149 histos->push_back(ptemp);
150 }
151}
152
[92237ca]153void HistogramsCollectionVsEta(std::vector<resolPlot> *histos, double etamin, double etamax, double ptmin, double ptmax, TString obj, double xmin = 0, double xmax = 2)
154{
155 resolPlot ptemp;
156 double width;
157 double etadown;
158 double etaup;
159
[341014c]160 for(int i = 0; i < Nbins; i++)
[92237ca]161 {
162 width = (etamax - etamin) / Nbins;
163 etadown = etamin + i * width;
[341014c]164 etaup = etamin + (i + 1) * width;
[92237ca]165
166 ptemp.set(etadown, etaup, ptmin, ptmax, obj, xmin, xmax);
167 histos->push_back(ptemp);
168 }
169}
170
[1408174]171//------------------------------------------------------------------------------
172
173class ExRootResult;
174class ExRootTreeReader;
175
176//------------------------------------------------------------------------------
177
[341014c]178void BinLogX(TH1 *h)
[1408174]179{
180 TAxis *axis = h->GetXaxis();
181 int bins = axis->GetNbins();
182
183 Axis_t from = axis->GetXmin();
184 Axis_t to = axis->GetXmax();
185 Axis_t width = (to - from) / bins;
186 Axis_t *new_bins = new Axis_t[bins + 1];
187
[341014c]188 for(int i = 0; i <= bins; i++)
[1408174]189 {
190 new_bins[i] = TMath::Power(10, from + i * width);
191 }
192 axis->Set(bins, new_bins);
193 delete new_bins;
194}
195
196//------------------------------------------------------------------------------
197
[341014c]198template <typename T>
199TH1D *GetEffPt(TClonesArray *branchReco, TClonesArray *branchParticle, TString name, int pdgID, double ptmin, double ptmax, double etamin, double etamax, ExRootTreeReader *treeReader)
[1408174]200{
201
[341014c]202 cout << "** Computing Efficiency of reconstructing " << branchReco->GetName() << " induced by " << branchParticle->GetName() << " with PID " << pdgID << endl;
[1408174]203
204 Long64_t allEntries = treeReader->GetEntries();
205
206 GenParticle *particle;
207 T *recoObj;
208
209 TLorentzVector recoMomentum, genMomentum, bestRecoMomentum;
210
211 Float_t deltaR;
212 Float_t pt, eta;
213 Long64_t entry;
214
215 Int_t i, j;
216
[341014c]217 TH1D *histGenPt = new TH1D(name + " gen spectra Pt", name + " gen spectra cen", Nbins, TMath::Log10(ptmin), TMath::Log10(ptmax));
218 TH1D *histRecoPt = new TH1D(name + " reco spectra Pt", name + " reco spectra cen", Nbins, TMath::Log10(ptmin), TMath::Log10(ptmax));
[1408174]219
[92237ca]220 histGenPt->SetDirectory(0);
221 histRecoPt->SetDirectory(0);
[1408174]222
[92237ca]223 BinLogX(histGenPt);
224 BinLogX(histRecoPt);
[1408174]225
226 // Loop over all events
227 for(entry = 0; entry < allEntries; ++entry)
228 {
229 // Load selected branches with data from specified event
230 treeReader->ReadEntry(entry);
231
232 // Loop over all generated particle in event
233 for(i = 0; i < branchParticle->GetEntriesFast(); ++i)
234 {
235
[341014c]236 particle = (GenParticle *)branchParticle->At(i);
[1408174]237 genMomentum = particle->P4();
238
239 deltaR = 999;
240
[341014c]241 pt = genMomentum.Pt();
[92237ca]242 eta = TMath::Abs(genMomentum.Eta());
243
[341014c]244 if(eta > etamax || eta < etamin) continue;
[cb159ed]245
[341014c]246 if(particle->PID == pdgID && genMomentum.Pt() > ptmin && genMomentum.Pt() < ptmax)
[534d945]247 //if (TMath::Abs(particle->PID) == pdgID && (particle->Status>20 && particle->Status <30) && genMomentum.Pt() > ptmin && genMomentum.Pt() < ptmax )
[92237ca]248 {
[1408174]249 // Loop over all reco object in event
250 for(j = 0; j < branchReco->GetEntriesFast(); ++j)
251 {
[341014c]252 recoObj = (T *)branchReco->At(j);
[1408174]253 recoMomentum = recoObj->P4();
[341014c]254 //if(Momentum.Px() == 0 && genMomentum.Py() == 0) continue;
[06a83b3]255
[1408174]256 // take the closest parton candidate
257 if(TMath::Abs(pdgID) == 5)
258 {
259 Jet *jet = (Jet *)recoObj;
[341014c]260 if(!(jet->BTag & (1 << 0))) continue;
[92237ca]261
262 //if(jet->BTag != ) continue;
263 }
264
[341014c]265 if(TMath::Abs(pdgID) == 4)
[92237ca]266 {
267 Jet *jet = (Jet *)recoObj;
[341014c]268 if(!(jet->BTag & (1 << 0))) continue;
[92237ca]269 }
[cb159ed]270
[341014c]271 if(TMath::Abs(pdgID) == 1)
[92237ca]272 {
273 Jet *jet = (Jet *)recoObj;
[341014c]274 if(!(jet->BTag & (1 << 0))) continue;
[1408174]275 }
[92237ca]276
[1408174]277 if(TMath::Abs(pdgID) == 15)
278 {
279 Jet *jet = (Jet *)recoObj;
280 if(jet->TauTag != 1) continue;
281 }
[92237ca]282
[1408174]283 if(genMomentum.DeltaR(recoMomentum) < deltaR)
284 {
285 deltaR = genMomentum.DeltaR(recoMomentum);
286 bestRecoMomentum = recoMomentum;
287 }
288 }
[92237ca]289 histGenPt->Fill(pt);
[341014c]290 if(deltaR < 0.3 && bestRecoMomentum.Pt() > 0.20 * pt)
291 {
292 histRecoPt->Fill(pt);
293 }
[92237ca]294 }
295 }
296 }
297
298 histRecoPt->Sumw2();
299 histGenPt->Sumw2();
300
301 histRecoPt->Divide(histGenPt);
302 histRecoPt->Scale(100.);
303
304 return histRecoPt;
305}
306
[341014c]307template <typename T>
308TH1D *GetEffEta(TClonesArray *branchReco, TClonesArray *branchParticle, TString name, int pdgID, double ptmin, double ptmax, double etamin, double etamax, ExRootTreeReader *treeReader)
[92237ca]309{
310
[341014c]311 cout << "** Computing Efficiency of reconstructing " << branchReco->GetName() << " induced by " << branchParticle->GetName() << " with PID " << pdgID << endl;
[92237ca]312
313 Long64_t allEntries = treeReader->GetEntries();
314
315 GenParticle *particle;
316 T *recoObj;
317
318 TLorentzVector recoMomentum, genMomentum, bestRecoMomentum;
319
320 Float_t deltaR;
321 Float_t pt, eta;
322 Long64_t entry;
323
324 Int_t i, j;
325
[341014c]326 TH1D *histGenEta = new TH1D(name + " gen spectra Eta", name + " gen spectra", Nbins, etamin, etamax);
327 TH1D *histRecoEta = new TH1D(name + " reco spectra Eta", name + " reco spectra", Nbins, etamin, etamax);
[92237ca]328
329 histGenEta->SetDirectory(0);
330 histRecoEta->SetDirectory(0);
331
332 // Loop over all events
333 for(entry = 0; entry < allEntries; ++entry)
334 {
335 // Load selected branches with data from specified event
336 treeReader->ReadEntry(entry);
337
338 // Loop over all generated particle in event
339 for(i = 0; i < branchParticle->GetEntriesFast(); ++i)
340 {
341
[341014c]342 particle = (GenParticle *)branchParticle->At(i);
[92237ca]343 genMomentum = particle->P4();
344
345 deltaR = 999;
346
[341014c]347 pt = genMomentum.Pt();
[92237ca]348 eta = genMomentum.Eta();
349
[341014c]350 if(pt > ptmax || pt < ptmin) continue;
[92237ca]351
[341014c]352 if(particle->PID == pdgID && genMomentum.Pt() > ptmin && genMomentum.Pt() < ptmax)
[534d945]353 //if (TMath::Abs(particle->PID) == pdgID && (particle->Status>20 && particle->Status <30) && genMomentum.Pt() > ptmin && genMomentum.Pt() < ptmax )
[92237ca]354 {
355 // Loop over all reco object in event
356 for(j = 0; j < branchReco->GetEntriesFast(); ++j)
[1408174]357 {
[341014c]358 recoObj = (T *)branchReco->At(j);
[92237ca]359 recoMomentum = recoObj->P4();
360 // this is simply to avoid warnings from initial state particle
361 // having infite rapidity ...
[341014c]362 //if(Momentum.Px() == 0 && genMomentum.Py() == 0) continue;
[1408174]363
[92237ca]364 // take the closest parton candidate
365 if(TMath::Abs(pdgID) == 5)
366 {
367 Jet *jet = (Jet *)recoObj;
[341014c]368 if(!(jet->BTag & (1 << 0))) continue;
[92237ca]369 }
370
371 if(TMath::Abs(pdgID) == 4)
372 {
373 Jet *jet = (Jet *)recoObj;
[341014c]374 if(!(jet->BTag & (1 << 0))) continue;
[92237ca]375 }
376
[341014c]377 if(TMath::Abs(pdgID) == 1)
[92237ca]378 {
379 Jet *jet = (Jet *)recoObj;
[341014c]380 if(!(jet->BTag & (1 << 0))) continue;
[92237ca]381 }
382
383 if(TMath::Abs(pdgID) == 15)
384 {
385 Jet *jet = (Jet *)recoObj;
386 if(jet->TauTag != 1) continue;
387 }
388 if(genMomentum.DeltaR(recoMomentum) < deltaR)
389 {
390 deltaR = genMomentum.DeltaR(recoMomentum);
391 bestRecoMomentum = recoMomentum;
392 }
[1408174]393 }
[92237ca]394
395 histGenEta->Fill(eta);
[341014c]396 if(deltaR < 0.3)
397 {
398 histRecoEta->Fill(eta);
399 }
[1408174]400 }
401 }
402 }
403
[92237ca]404 histRecoEta->Sumw2();
405 histGenEta->Sumw2();
[1408174]406
[92237ca]407 histRecoEta->Divide(histGenEta);
408 histRecoEta->Scale(100.);
[1408174]409
[92237ca]410 return histRecoEta;
411}
[1408174]412
[534d945]413//------------------------------------------------------------------------------
414
[341014c]415template <typename T>
416TH1D *GetJetEffPt(TClonesArray *branchJet, TString name, int pdgID, double ptmin, double ptmax, double etamin, double etamax, ExRootTreeReader *treeReader)
[534d945]417{
418
[341014c]419 cout << "** Computing Efficiency of reconstructing " << branchJet->GetName() << " with PID " << pdgID << endl;
[534d945]420
421 Long64_t allEntries = treeReader->GetEntries();
422
423 T *recoObj;
424
425 TLorentzVector recoMomentum, genMomentum, bestRecoMomentum;
426
427 Float_t pt, eta;
428 Long64_t entry;
429
430 Int_t j;
431
[341014c]432 TH1D *histGenPt = new TH1D(name + " gen spectra Pt", name + " gen spectra cen", Nbins, TMath::Log10(ptmin), TMath::Log10(ptmax));
433 TH1D *histRecoPt = new TH1D(name + " reco spectra Pt", name + " reco spectra cen", Nbins, TMath::Log10(ptmin), TMath::Log10(ptmax));
[534d945]434
435 histGenPt->SetDirectory(0);
436 histRecoPt->SetDirectory(0);
437
438 BinLogX(histGenPt);
439 BinLogX(histRecoPt);
440
441 // Loop over all events
442 for(entry = 0; entry < allEntries; ++entry)
443 {
444 // Load selected branches with data from specified event
445 treeReader->ReadEntry(entry);
446
447 // Loop over all reco object in event
448 for(j = 0; j < branchJet->GetEntriesFast(); ++j)
449 {
[341014c]450 recoObj = (T *)branchJet->At(j);
[534d945]451 recoMomentum = recoObj->P4();
[341014c]452 pt = recoMomentum.Pt();
[534d945]453 eta = TMath::Abs(recoMomentum.Eta());
454 Jet *jet = (Jet *)recoObj;
[341014c]455
456 if(eta > etamax || eta < etamin) continue;
[534d945]457 if(pt < ptmin || pt > ptmax) continue;
[341014c]458
[534d945]459 Int_t flavor = jet->Flavor;
460 if(flavor == 21) flavor = 0;
[341014c]461
[534d945]462 if(TMath::Abs(pdgID) == 1)
463 {
[341014c]464
[534d945]465 if(flavor < 4)
466 {
467 histGenPt->Fill(pt);
[341014c]468 if(jet->BTag & (1 << 0)) histRecoPt->Fill(pt);
[534d945]469 }
470 }
471 if(TMath::Abs(pdgID) == 4)
472 {
473 if(flavor == 4)
474 {
475 histGenPt->Fill(pt);
[341014c]476 if(jet->BTag & (1 << 0)) histRecoPt->Fill(pt);
[534d945]477 }
478 }
479 if(TMath::Abs(pdgID) == 5)
480 {
481 if(flavor == 5)
482 {
483 histGenPt->Fill(pt);
[341014c]484 if(jet->BTag & (1 << 0)) histRecoPt->Fill(pt);
[534d945]485 }
486 }
[341014c]487 }
[534d945]488 }
489
490 histRecoPt->Sumw2();
491 histGenPt->Sumw2();
492
493 histRecoPt->Divide(histGenPt);
494 histRecoPt->Scale(100.);
495
496 return histRecoPt;
497}
498
499// ------------------------------------------------------------------------------------------------------------------------------------------------------
500
[341014c]501template <typename T>
502TH1D *GetJetEffEta(TClonesArray *branchJet, TString name, int pdgID, double ptmin, double ptmax, double etamin, double etamax, ExRootTreeReader *treeReader)
[534d945]503{
504
[341014c]505 cout << "** Computing Efficiency of reconstructing " << branchJet->GetName() << " with PID " << pdgID << endl;
[534d945]506
507 Long64_t allEntries = treeReader->GetEntries();
508
509 T *recoObj;
510
511 TLorentzVector recoMomentum, genMomentum, bestRecoMomentum;
512
513 Float_t pt, eta;
514 Long64_t entry;
515
516 Int_t j;
517
[341014c]518 TH1D *histGenEta = new TH1D(name + " gen spectra Eta", name + " gen spectra", Nbins, etamin, etamax);
519 TH1D *histRecoEta = new TH1D(name + " reco spectra Eta", name + " reco spectra", Nbins, etamin, etamax);
[534d945]520
521 histGenEta->SetDirectory(0);
522 histRecoEta->SetDirectory(0);
523 // Loop over all events
524 for(entry = 0; entry < allEntries; ++entry)
525 {
526 // Load selected branches with data from specified event
527 treeReader->ReadEntry(entry);
528
529 // Loop over all reco object in event
530 for(j = 0; j < branchJet->GetEntriesFast(); ++j)
531 {
[341014c]532 recoObj = (T *)branchJet->At(j);
[534d945]533 recoMomentum = recoObj->P4();
[341014c]534 pt = recoMomentum.Pt();
[534d945]535 eta = recoMomentum.Eta();
536 Jet *jet = (Jet *)recoObj;
[341014c]537
538 if(eta > etamax || eta < etamin) continue;
[534d945]539 if(pt < ptmin || pt > ptmax) continue;
[341014c]540
[534d945]541 Int_t flavor = jet->Flavor;
542 if(flavor == 21) flavor = 0;
[341014c]543
[534d945]544 if(TMath::Abs(pdgID) == 1)
545 {
[341014c]546 if(flavor == 1 || flavor == 21)
[534d945]547 {
548 histGenEta->Fill(eta);
[341014c]549 if(jet->BTag & (1 << 0)) histRecoEta->Fill(eta);
[534d945]550 }
551 }
552 if(TMath::Abs(pdgID) == 4)
553 {
554 if(flavor == 4)
555 {
556 histGenEta->Fill(eta);
[341014c]557 if(jet->BTag & (1 << 0)) histRecoEta->Fill(eta);
[534d945]558 }
559 }
560 if(TMath::Abs(pdgID) == 5)
561 {
562 if(flavor == 5)
563 {
564 histGenEta->Fill(eta);
[341014c]565 if(jet->BTag & (1 << 0)) histRecoEta->Fill(eta);
[534d945]566 }
567 }
[341014c]568 }
[534d945]569 }
570
571 histRecoEta->Sumw2();
572 histGenEta->Sumw2();
573
574 histRecoEta->Divide(histGenEta);
575 histRecoEta->Scale(100.);
576
577 return histRecoEta;
578}
579
580// -----------------------------------------------------------------------------------------------------------------------------------------------------
[1408174]581
[341014c]582template <typename T>
583TH1D *GetTauEffPt(TClonesArray *branchReco, TClonesArray *branchParticle, TString name, int pdgID, double ptmin, double ptmax, double etamin, double etamax, ExRootTreeReader *treeReader)
[1408174]584{
585
[341014c]586 cout << "** Computing Efficiency of reconstructing " << branchReco->GetName() << " induced by " << branchParticle->GetName() << " with PID " << pdgID << endl;
[92237ca]587
588 Long64_t allEntries = treeReader->GetEntries();
[1408174]589
590 GenParticle *particle;
[92237ca]591 T *recoObj;
[1408174]592
[92237ca]593 TLorentzVector recoMomentum, genMomentum, bestRecoMomentum;
[1408174]594
595 Float_t deltaR;
596 Float_t pt, eta;
597 Long64_t entry;
598
[92237ca]599 Int_t i, j;
600
[341014c]601 TH1D *histGenPt = new TH1D(name + " gen spectra Pt", name + " gen spectra cen", Nbins, TMath::Log10(ptmin), TMath::Log10(ptmax));
602 TH1D *histRecoPt = new TH1D(name + " reco spectra Pt", name + " reco spectra cen", Nbins, TMath::Log10(ptmin), TMath::Log10(ptmax));
[92237ca]603
604 histGenPt->SetDirectory(0);
605 histRecoPt->SetDirectory(0);
606
607 BinLogX(histGenPt);
608 BinLogX(histRecoPt);
[1408174]609
610 // Loop over all events
611 for(entry = 0; entry < allEntries; ++entry)
612 {
613 // Load selected branches with data from specified event
614 treeReader->ReadEntry(entry);
615
[92237ca]616 // Loop over all generated particle in event
617 for(i = 0; i < branchParticle->GetEntriesFast(); ++i)
[1408174]618 {
[92237ca]619
[341014c]620 particle = (GenParticle *)branchParticle->At(i);
[92237ca]621 genMomentum = particle->P4();
[1408174]622
623 deltaR = 999;
624
[341014c]625 pt = genMomentum.Pt();
[92237ca]626 eta = TMath::Abs(genMomentum.Eta());
[1408174]627
[341014c]628 if(eta > etamax || eta < etamin) continue;
[92237ca]629
[341014c]630 if(particle->PID == pdgID && genMomentum.Pt() > ptmin && genMomentum.Pt() < ptmax)
[92237ca]631 {
632 // Loop over all reco object in event
633 for(j = 0; j < branchReco->GetEntriesFast(); ++j)
634 {
[341014c]635 recoObj = (T *)branchReco->At(j);
[92237ca]636 recoMomentum = recoObj->P4();
[1408174]637 // this is simply to avoid warnings from initial state particle
638 // having infite rapidity ...
[92237ca]639 //if(Momentum.Px() == 0 && genMomentum.Py() == 0) continue;
640
641 if(TMath::Abs(pdgID) == 1)
642 {
643 Jet *jet = (Jet *)recoObj;
[341014c]644 if(jet->TauTag != 1) continue;
[92237ca]645 }
646
647 if(TMath::Abs(pdgID) == 15)
648 {
649 Jet *jet = (Jet *)recoObj;
650 if(jet->TauTag != 1) continue;
651 }
[1408174]652
653 if(genMomentum.DeltaR(recoMomentum) < deltaR)
654 {
[92237ca]655 deltaR = genMomentum.DeltaR(recoMomentum);
656 bestRecoMomentum = recoMomentum;
[1408174]657 }
658 }
[92237ca]659
660 histGenPt->Fill(pt);
[341014c]661 if(deltaR < 0.3)
662 {
663 histRecoPt->Fill(pt);
664 }
[1408174]665 }
[92237ca]666 }
667 }
[1408174]668
[92237ca]669 histRecoPt->Sumw2();
670 histGenPt->Sumw2();
[1408174]671
[92237ca]672 histRecoPt->Divide(histGenPt);
673 histRecoPt->Scale(100.);
[341014c]674 if(TMath::Abs(pdgID) == 15) histRecoPt->Scale(1 / 0.648);
[92237ca]675
676 return histRecoPt;
677}
678
[341014c]679template <typename T>
680TH1D *GetTauEffEta(TClonesArray *branchReco, TClonesArray *branchParticle, TString name, int pdgID, double ptmin, double ptmax, double etamin, double etamax, ExRootTreeReader *treeReader)
[92237ca]681{
682
[341014c]683 cout << "** Computing Efficiency of reconstructing " << branchReco->GetName() << " induced by " << branchParticle->GetName() << " with PID " << pdgID << endl;
[92237ca]684
685 Long64_t allEntries = treeReader->GetEntries();
686
687 GenParticle *particle;
688 T *recoObj;
689
690 TLorentzVector recoMomentum, genMomentum, bestRecoMomentum;
691
692 Float_t deltaR;
693 Float_t pt, eta;
694 Long64_t entry;
695
696 Int_t i, j;
697
[341014c]698 TH1D *histGenEta = new TH1D(name + " gen spectra Eta", name + " gen spectra", Nbins, etamin, etamax);
699 TH1D *histRecoEta = new TH1D(name + " reco spectra Eta", name + " reco spectra", Nbins, etamin, etamax);
[92237ca]700
701 histGenEta->SetDirectory(0);
702 histRecoEta->SetDirectory(0);
703
704 // Loop over all events
705 for(entry = 0; entry < allEntries; ++entry)
706 {
707 // Load selected branches with data from specified event
708 treeReader->ReadEntry(entry);
709
710 // Loop over all generated particle in event
711 for(i = 0; i < branchParticle->GetEntriesFast(); ++i)
712 {
713
[341014c]714 particle = (GenParticle *)branchParticle->At(i);
[92237ca]715 genMomentum = particle->P4();
716
717 deltaR = 999;
718
[341014c]719 pt = genMomentum.Pt();
[92237ca]720 eta = genMomentum.Eta();
721
[341014c]722 if(pt > ptmax || pt < ptmin) continue;
[92237ca]723
[341014c]724 if(particle->PID == pdgID && genMomentum.Pt() > ptmin && genMomentum.Pt() < ptmax)
[92237ca]725 {
726 // Loop over all reco object in event
727 for(j = 0; j < branchReco->GetEntriesFast(); ++j)
[1408174]728 {
[341014c]729 recoObj = (T *)branchReco->At(j);
[92237ca]730 recoMomentum = recoObj->P4();
731 // this is simply to avoid warnings from initial state particle
732 // having infite rapidity ...
[341014c]733 //if(Momentum.Px() == 0 && genMomentum.Py() == 0) continue;
[92237ca]734
[341014c]735 if(TMath::Abs(pdgID) == 1)
[92237ca]736 {
737 Jet *jet = (Jet *)recoObj;
[341014c]738 if(jet->TauTag != 1) continue;
[92237ca]739 }
740
741 if(TMath::Abs(pdgID) == 15)
[1408174]742 {
[92237ca]743 Jet *jet = (Jet *)recoObj;
744 if(jet->TauTag != 1) continue;
745 }
746
747 if(genMomentum.DeltaR(recoMomentum) < deltaR)
748 {
749 deltaR = genMomentum.DeltaR(recoMomentum);
750 bestRecoMomentum = recoMomentum;
[1408174]751 }
752 }
[92237ca]753
754 histGenEta->Fill(eta);
[341014c]755 if(deltaR < 0.3)
756 {
757 histRecoEta->Fill(eta);
758 }
[1408174]759 }
760 }
761 }
[92237ca]762
763 histRecoEta->Sumw2();
764 histGenEta->Sumw2();
765
766 histRecoEta->Divide(histGenEta);
767 histRecoEta->Scale(100.);
[341014c]768 if(TMath::Abs(pdgID) == 15) histRecoEta->Scale(1 / 0.648);
[92237ca]769
770 return histRecoEta;
[1408174]771}
[32677fb]772
[341014c]773template <typename T>
[92237ca]774void GetPtres(std::vector<resolPlot> *histos, TClonesArray *branchReco, TClonesArray *branchParticle, int pdgID, Double_t etaMin, Double_t etaMax, ExRootTreeReader *treeReader)
[32677fb]775{
776 Long64_t allEntries = treeReader->GetEntries();
777
778 cout << "** Computing pt resolution of " << branchReco->GetName() << " induced by " << branchParticle->GetName() << " with PID " << pdgID << endl;
779
780 GenParticle *particle;
[341014c]781 T *recoObj;
[32677fb]782
783 TLorentzVector recoMomentum, genMomentum, bestGenMomentum;
784
785 Float_t deltaR;
786 Float_t pt, eta;
787 Long64_t entry;
788
789 Int_t i, j, bin;
790
791 // Loop over all events
792 for(entry = 0; entry < allEntries; ++entry)
793 {
794 // Load selected branches with data from specified event
795 treeReader->ReadEntry(entry);
796
797 // Loop over all reconstructed jets in event
798 for(i = 0; i < branchReco->GetEntriesFast(); ++i)
799 {
[341014c]800 recoObj = (T *)branchReco->At(i);
[32677fb]801 recoMomentum = recoObj->P4();
802
803 deltaR = 999;
804
[341014c]805 // Loop over all hard partons in event
806 for(j = 0; j < branchParticle->GetEntriesFast(); ++j)
807 {
808 particle = (GenParticle *)branchParticle->At(j);
809 if(particle->PID == pdgID && particle->Status == 1)
[32677fb]810 {
811 genMomentum = particle->P4();
812
813 // this is simply to avoid warnings from initial state particle
814 // having infite rapidity ...
815 if(genMomentum.Px() == 0 && genMomentum.Py() == 0) continue;
816
817 // take the closest parton candidate
818 if(genMomentum.DeltaR(recoMomentum) < deltaR)
819 {
[341014c]820 deltaR = genMomentum.DeltaR(recoMomentum);
821 bestGenMomentum = genMomentum;
[32677fb]822 }
823 }
824 }
825
826 if(deltaR < 0.3)
827 {
[341014c]828 pt = bestGenMomentum.Pt();
[32677fb]829 eta = TMath::Abs(bestGenMomentum.Eta());
830
[341014c]831 for(bin = 0; bin < Nbins; bin++)
[32677fb]832 {
[92237ca]833 if(pt > histos->at(bin).ptmin && pt < histos->at(bin).ptmax && eta > etaMin && eta < etaMax)
[32677fb]834 {
[341014c]835 histos->at(bin).resolHist->Fill(recoMomentum.Pt() / bestGenMomentum.Pt());
[32677fb]836 }
837 }
838 }
839 }
840 }
841}
842
[341014c]843template <typename T>
[92237ca]844void GetEres(std::vector<resolPlot> *histos, TClonesArray *branchReco, TClonesArray *branchParticle, int pdgID, Double_t etaMin, Double_t etaMax, ExRootTreeReader *treeReader)
[1408174]845{
846 Long64_t allEntries = treeReader->GetEntries();
847
[92237ca]848 cout << "** Computing e resolution of " << branchReco->GetName() << " induced by " << branchParticle->GetName() << " with PID " << pdgID << endl;
849
850 GenParticle *particle;
[341014c]851 T *recoObj;
[92237ca]852
853 TLorentzVector recoMomentum, genMomentum, bestGenMomentum;
854
855 Float_t deltaR;
856 Float_t e, eta;
857 Long64_t entry;
858
859 Int_t i, j, bin;
860
861 // Loop over all events
862 for(entry = 0; entry < allEntries; ++entry)
863 {
864 // Load selected branches with data from specified event
865 treeReader->ReadEntry(entry);
866
[341014c]867 // Loop over all reconstructed jets in event
[92237ca]868 for(i = 0; i < branchReco->GetEntriesFast(); ++i)
869 {
[341014c]870 recoObj = (T *)branchReco->At(i);
[92237ca]871 recoMomentum = recoObj->P4();
872
873 deltaR = 999;
874
[341014c]875 // Loop over all hard partons in event
876 for(j = 0; j < branchParticle->GetEntriesFast(); ++j)
877 {
878 particle = (GenParticle *)branchParticle->At(j);
879 if(particle->PID == pdgID && particle->Status == 1)
[92237ca]880 {
881 genMomentum = particle->P4();
882
883 // this is simply to avoid warnings from initial state particle
884 // having infite rapidity ...
885 if(genMomentum.Px() == 0 && genMomentum.Py() == 0) continue;
886
887 // take the closest parton candidate
888 if(genMomentum.DeltaR(recoMomentum) < deltaR)
889 {
[341014c]890 deltaR = genMomentum.DeltaR(recoMomentum);
891 bestGenMomentum = genMomentum;
[92237ca]892 }
893 }
894 }
895
896 if(deltaR < 0.3)
897 {
[341014c]898 e = bestGenMomentum.E();
[92237ca]899 eta = TMath::Abs(bestGenMomentum.Eta());
900
[341014c]901 for(bin = 0; bin < Nbins; bin++)
[92237ca]902 {
903 if(e > histos->at(bin).ptmin && e < histos->at(bin).ptmax && eta > etaMin && eta < etaMax)
904 {
[341014c]905 histos->at(bin).resolHist->Fill(recoMomentum.E() / bestGenMomentum.E());
[92237ca]906 }
907 }
908 }
909 }
910 }
911}
912
[341014c]913template <typename T>
[92237ca]914void GetPtresVsEta(std::vector<resolPlot> *histos, TClonesArray *branchReco, TClonesArray *branchParticle, int pdgID, Double_t ptMin, Double_t ptMax, ExRootTreeReader *treeReader)
915{
916 Long64_t allEntries = treeReader->GetEntries();
917
918 cout << "** Computing pt resolution of " << branchReco->GetName() << " induced by " << branchParticle->GetName() << " with PID " << pdgID << endl;
919
920 GenParticle *particle;
[341014c]921 T *recoObj;
[92237ca]922
923 TLorentzVector recoMomentum, genMomentum, bestGenMomentum;
924
925 Float_t deltaR;
926 Float_t pt, eta;
927 Long64_t entry;
928
929 Int_t i, j, bin;
930
931 // Loop over all events
932 for(entry = 0; entry < allEntries; ++entry)
933 {
934 // Load selected branches with data from specified event
935 treeReader->ReadEntry(entry);
936
937 // Loop over all reconstructed jets in event
938 for(i = 0; i < branchReco->GetEntriesFast(); ++i)
939 {
[341014c]940 recoObj = (T *)branchReco->At(i);
[92237ca]941 recoMomentum = recoObj->P4();
942
943 deltaR = 999;
944
[341014c]945 // Loop over all hard partons in event
946 for(j = 0; j < branchParticle->GetEntriesFast(); ++j)
947 {
948 particle = (GenParticle *)branchParticle->At(j);
949 if(particle->PID == pdgID && particle->Status == 1)
[92237ca]950 {
951 genMomentum = particle->P4();
952
953 // this is simply to avoid warnings from initial state particle
954 // having infite rapidity ...
955 if(genMomentum.Px() == 0 && genMomentum.Py() == 0) continue;
956
957 // take the closest parton candidate
958 if(genMomentum.DeltaR(recoMomentum) < deltaR)
959 {
[341014c]960 deltaR = genMomentum.DeltaR(recoMomentum);
961 bestGenMomentum = genMomentum;
[92237ca]962 }
963 }
964 }
965
966 if(deltaR < 0.3)
967 {
[341014c]968 pt = bestGenMomentum.Pt();
[92237ca]969 eta = bestGenMomentum.Eta();
970
[341014c]971 for(bin = 0; bin < Nbins; bin++)
[92237ca]972 {
973 if(eta > histos->at(bin).etamin && eta < histos->at(bin).etamax && pt > ptMin && pt < ptMax)
974 {
[341014c]975 histos->at(bin).resolHist->Fill(recoMomentum.Pt() / bestGenMomentum.Pt());
[92237ca]976 }
977 }
978 }
979 }
980 }
981}
982
[341014c]983template <typename T>
[92237ca]984void GetEresVsEta(std::vector<resolPlot> *histos, TClonesArray *branchReco, TClonesArray *branchParticle, int pdgID, Double_t eMin, Double_t eMax, ExRootTreeReader *treeReader)
985{
986 Long64_t allEntries = treeReader->GetEntries();
987
988 cout << "** Computing E resolution of " << branchReco->GetName() << " induced by " << branchParticle->GetName() << " with PID " << pdgID << endl;
989
990 GenParticle *particle;
[341014c]991 T *recoObj;
[92237ca]992
993 TLorentzVector recoMomentum, genMomentum, bestGenMomentum;
994
995 Float_t deltaR;
996 Float_t e, eta;
997 Long64_t entry;
998
999 Int_t i, j, bin;
1000
1001 // Loop over all events
1002 for(entry = 0; entry < allEntries; ++entry)
1003 {
1004 // Load selected branches with data from specified event
1005 treeReader->ReadEntry(entry);
1006
1007 // Loop over all reconstructed jets in event
1008 for(i = 0; i < branchReco->GetEntriesFast(); ++i)
1009 {
[341014c]1010 recoObj = (T *)branchReco->At(i);
[92237ca]1011 recoMomentum = recoObj->P4();
1012
1013 deltaR = 999;
1014
[341014c]1015 // Loop over all hard partons in event
1016 for(j = 0; j < branchParticle->GetEntriesFast(); ++j)
1017 {
1018 particle = (GenParticle *)branchParticle->At(j);
1019 if(particle->PID == pdgID && particle->Status == 1)
[92237ca]1020 {
1021 genMomentum = particle->P4();
1022
1023 // this is simply to avoid warnings from initial state particle
1024 // having infite rapidity ...
1025 if(genMomentum.Px() == 0 && genMomentum.Py() == 0) continue;
1026
1027 // take the closest parton candidate
1028 if(genMomentum.DeltaR(recoMomentum) < deltaR)
1029 {
[341014c]1030 deltaR = genMomentum.DeltaR(recoMomentum);
1031 bestGenMomentum = genMomentum;
[92237ca]1032 }
1033 }
1034 }
1035
1036 if(deltaR < 0.3)
1037 {
[341014c]1038 e = bestGenMomentum.E();
[92237ca]1039 eta = bestGenMomentum.Eta();
1040
[341014c]1041 for(bin = 0; bin < Nbins; bin++)
[92237ca]1042 {
1043 if(eta > histos->at(bin).etamin && eta < histos->at(bin).etamax && e > eMin && e < eMax)
1044 {
[341014c]1045 histos->at(bin).resolHist->Fill(recoMomentum.E() / bestGenMomentum.E());
[92237ca]1046 }
1047 }
1048 }
1049 }
1050 }
1051}
1052
1053void GetJetsEres(std::vector<resolPlot> *histos, TClonesArray *branchJet, TClonesArray *branchGenJet, ExRootTreeReader *treeReader, Double_t etaMin, Double_t etaMax)
1054{
1055
1056 Long64_t allEntries = treeReader->GetEntries();
1057
1058 cout << "** Computing resolution of " << branchJet->GetName() << " induced by " << branchGenJet->GetName() << endl;
[1408174]1059
1060 Jet *jet, *genjet;
1061
[92237ca]1062 TLorentzVector jetMomentum, genJetMomentum, bestGenJetMomentum;
1063
1064 Float_t deltaR;
1065 Float_t pt, eta;
1066 Long64_t entry;
1067
1068 Int_t i, j, bin;
1069
1070 // Loop over all events
1071 for(entry = 0; entry < allEntries; ++entry)
1072 {
1073 // Load selected branches with data from specified event
1074 treeReader->ReadEntry(entry);
1075
[341014c]1076 if(entry % 10000 == 0) cout << "Event number: " << entry << endl;
[92237ca]1077
1078 // Loop over all reconstructed jets in event
[341014c]1079 for(i = 0; i < TMath::Min(2, branchJet->GetEntriesFast()); ++i) //branchJet->GetEntriesFast(); ++i)
[92237ca]1080 {
1081
[341014c]1082 jet = (Jet *)branchJet->At(i);
[92237ca]1083 jetMomentum = jet->P4();
1084
1085 deltaR = 999;
1086
[341014c]1087 // Loop over all hard partons in event
1088 for(j = 0; j < TMath::Min(2, branchGenJet->GetEntriesFast()); ++j)
1089 {
1090 genjet = (Jet *)branchGenJet->At(j);
[92237ca]1091
1092 genJetMomentum = genjet->P4();
1093
1094 // this is simply to avoid warnings from initial state particle
1095 // having infite rapidity ...
1096 if(genJetMomentum.Px() == 0 && genJetMomentum.Py() == 0) continue;
1097
1098 // take the closest parton candidate
1099 if(genJetMomentum.DeltaR(jetMomentum) < deltaR)
1100 {
[341014c]1101 deltaR = genJetMomentum.DeltaR(jetMomentum);
1102 bestGenJetMomentum = genJetMomentum;
[92237ca]1103 }
1104 }
1105
1106 if(deltaR < 0.3)
1107 {
[341014c]1108 pt = genJetMomentum.E();
[92237ca]1109 eta = genJetMomentum.Eta();
1110
[341014c]1111 for(bin = 0; bin < Nbins; bin++)
[92237ca]1112 {
1113 if(pt > histos->at(bin).ptmin && pt < histos->at(bin).ptmax && eta < etaMax && eta > etaMin)
1114 {
[341014c]1115 histos->at(bin).resolHist->Fill(jetMomentum.E() / bestGenJetMomentum.E());
[92237ca]1116 }
1117 }
1118 }
1119 }
1120 }
1121}
1122
1123void GetJetsEresVsEta(std::vector<resolPlot> *histos, TClonesArray *branchJet, TClonesArray *branchGenJet, ExRootTreeReader *treeReader, Double_t eMin, Double_t eMax)
1124{
1125
1126 Long64_t allEntries = treeReader->GetEntries();
1127
1128 cout << "** Computing resolution of " << branchJet->GetName() << " induced by " << branchGenJet->GetName() << endl;
1129
1130 Jet *jet, *genjet;
1131
1132 TLorentzVector jetMomentum, genJetMomentum, bestGenJetMomentum;
1133
1134 Float_t deltaR;
1135 Float_t pt, eta;
1136 Long64_t entry;
1137
1138 Int_t i, j, bin;
1139
1140 // Loop over all events
1141 for(entry = 0; entry < allEntries; ++entry)
1142 {
1143 // Load selected branches with data from specified event
1144 treeReader->ReadEntry(entry);
1145
[341014c]1146 if(entry % 10000 == 0) cout << "Event number: " << entry << endl;
[92237ca]1147
1148 // Loop over all reconstructed jets in event
[341014c]1149 for(i = 0; i < TMath::Min(2, branchJet->GetEntriesFast()); ++i) //branchJet->GetEntriesFast(); ++i)
[92237ca]1150 {
1151
[341014c]1152 jet = (Jet *)branchJet->At(i);
[92237ca]1153 jetMomentum = jet->P4();
1154
1155 deltaR = 999;
1156
[341014c]1157 // Loop over all hard partons in event
1158 for(j = 0; j < TMath::Min(2, branchGenJet->GetEntriesFast()); ++j)
1159 {
1160 genjet = (Jet *)branchGenJet->At(j);
[92237ca]1161
1162 genJetMomentum = genjet->P4();
1163
1164 // this is simply to avoid warnings from initial state particle
1165 // having infite rapidity ...
1166 if(genJetMomentum.Px() == 0 && genJetMomentum.Py() == 0) continue;
1167
1168 // take the closest parton candidate
1169 if(genJetMomentum.DeltaR(jetMomentum) < deltaR)
1170 {
[341014c]1171 deltaR = genJetMomentum.DeltaR(jetMomentum);
1172 bestGenJetMomentum = genJetMomentum;
[92237ca]1173 }
1174 }
1175
1176 if(deltaR < 0.3)
1177 {
1178
[341014c]1179 pt = genJetMomentum.E();
[92237ca]1180 eta = genJetMomentum.Eta();
1181
[341014c]1182 for(bin = 0; bin < Nbins; bin++)
[92237ca]1183 {
1184 if(eta > histos->at(bin).etamin && eta < histos->at(bin).etamax && pt < eMax && pt > eMin)
1185 {
[341014c]1186 histos->at(bin).resolHist->Fill(jetMomentum.E() / bestGenJetMomentum.E());
[92237ca]1187 }
1188 }
1189 }
1190 }
1191 }
1192}
1193
[341014c]1194std::pair<Double_t, Double_t> GausFit(TH1 *hist)
[92237ca]1195{
[341014c]1196 TF1 *f1 = new TF1("f1", "gaus", hist->GetMean() - 2 * hist->GetRMS(), hist->GetMean() + 2 * hist->GetRMS());
1197 hist->Fit("f1", "RQ");
[92237ca]1198
[341014c]1199 TF1 *f2 = new TF1("f2", "gaus", f1->GetParameter(1) - 2 * f1->GetParameter(2), f1->GetParameter(1) + 2 * f1->GetParameter(2));
1200 hist->Fit("f2", "RQ");
[92237ca]1201
1202 Double_t sig = f2->GetParameter(2);
1203 Double_t sigErr = f2->GetParError(2);
1204
1205 delete f1;
1206 delete f2;
[341014c]1207 return make_pair(sig, sigErr);
[92237ca]1208}
1209
1210TGraphErrors EresGraph(std::vector<resolPlot> *histos, bool rms = false)
1211{
1212 Int_t bin;
1213 Int_t count = 0;
[341014c]1214 TGraphErrors gr = TGraphErrors(Nbins / 2);
[92237ca]1215 double val, error;
[341014c]1216 for(bin = 0; bin < Nbins; bin++)
[92237ca]1217 {
[341014c]1218 std::pair<Double_t, Double_t> sigvalues = GausFit(histos->at(bin).resolHist);
1219 if(rms == true)
1220 {
1221 gr.SetPoint(count, (histos->at(bin).ptmin + histos->at(bin).ptmax) / 2.0, 100 * histos->at(bin).resolHist->GetRMS());
1222 //gr.SetPointError(count,0, 100*sigvalues.second); // to correct
1223 error = 100 * histos->at(bin).resolHist->GetRMSError();
1224 val = 100 * histos->at(bin).resolHist->GetRMS();
1225 if(error > 0.2 * val) error = 0.2 * val;
1226 gr.SetPointError(count, 0, error); // to correct
1227 }
1228 else
1229 {
1230
1231 gr.SetPoint(count, (histos->at(bin).ptmin + histos->at(bin).ptmax) / 2.0, 100 * sigvalues.first);
1232 error = 100 * sigvalues.second;
1233 val = 100 * sigvalues.first;
1234 if(error > 0.2 * val) error = 0.2 * val;
1235 gr.SetPointError(count, 0, error); // to correct
1236 //gr.SetPointError(count,0, 100*sigvalues.second);
1237 }
1238 count++;
[92237ca]1239 }
1240
1241 return gr;
1242}
1243
1244TGraphErrors MetResGraph(std::vector<resolPlot> *histos, bool rms = false)
1245{
1246 Int_t bin;
1247 Int_t count = 0;
[341014c]1248 TGraphErrors gr = TGraphErrors(Nbins / 2);
[92237ca]1249 double val, error;
[341014c]1250 for(bin = 0; bin < Nbins; bin++)
[92237ca]1251 {
[341014c]1252 std::pair<Double_t, Double_t> sigvalues = GausFit(histos->at(bin).resolHist);
1253 if(rms == true)
1254 {
1255 gr.SetPoint(count, (histos->at(bin).ptmin + histos->at(bin).ptmax) / 2.0, histos->at(bin).resolHist->GetRMS());
1256 error = histos->at(bin).resolHist->GetRMSError();
1257 val = histos->at(bin).resolHist->GetRMS();
1258 if(error > 0.2 * val) error = 0.2 * val;
1259 gr.SetPointError(count, 0, error); // to correct
1260 }
1261 else
1262 {
1263
1264 gr.SetPoint(count, (histos->at(bin).ptmin + histos->at(bin).ptmax) / 2.0, sigvalues.first);
1265 val = sigvalues.first;
1266 error = sigvalues.second;
1267 if(error > 0.2 * val) error = 0.2 * val;
1268 gr.SetPointError(count, 0, error);
1269 //gr.SetPointError(count,0, 100*sigvalues.second);
1270 }
1271 count++;
[92237ca]1272 }
1273
1274 return gr;
1275}
1276
1277TGraphErrors EresGraphVsEta(std::vector<resolPlot> *histos, bool rms = false)
1278{
1279 Int_t bin;
1280 Int_t count = 0;
[341014c]1281 TGraphErrors gr = TGraphErrors(Nbins / 2);
[92237ca]1282 double val, error;
[341014c]1283 for(bin = 0; bin < Nbins; bin++)
[92237ca]1284 {
1285
[341014c]1286 std::pair<Double_t, Double_t> sigvalues = GausFit(histos->at(bin).resolHist);
1287 if(rms == true)
1288 {
1289 gr.SetPoint(count, (histos->at(bin).etamin + histos->at(bin).etamax) / 2.0, histos->at(bin).resolHist->GetRMS());
1290 error = 100 * histos->at(bin).resolHist->GetRMSError();
1291 val = 100 * histos->at(bin).resolHist->GetRMS();
1292 if(error > 0.2 * val) error = 0.2 * val;
1293 gr.SetPointError(count, 0, error); // to correct
1294 }
1295 else
1296 {
1297 gr.SetPoint(count, (histos->at(bin).etamin + histos->at(bin).etamax) / 2.0, 100 * sigvalues.first);
1298 val = 100 * sigvalues.first;
1299 error = 100 * sigvalues.second;
1300 if(error > 0.2 * val) error = 0.2 * val;
1301 gr.SetPointError(count, 0, error);
1302 //gr.SetPointError(count,0, 100*sigvalues.second);
1303 }
1304 count++;
[92237ca]1305 }
1306
1307 return gr;
1308}
1309
1310void GetMetres(std::vector<resolPlot> *histos, TClonesArray *branchScalarHT, TClonesArray *branchMet, TClonesArray *branchJet, ExRootTreeReader *treeReader)
1311{
1312
1313 Long64_t allEntries = treeReader->GetEntries();
1314
1315 cout << "** Computing resolution of " << branchMet->GetName() << " vs " << branchScalarHT->GetName() << endl;
1316
1317 MissingET *met;
1318 ScalarHT *scalarHT;
1319
1320 Long64_t entry;
1321
1322 Int_t bin;
1323 Double_t ht;
1324
1325 Jet *jet;
1326 TLorentzVector p1, p2;
1327
1328 // Loop over all events
1329 for(entry = 0; entry < allEntries; ++entry)
1330 {
1331 // Load selected branches with data from specified event
1332 treeReader->ReadEntry(entry);
1333
[341014c]1334 if(entry % 10000 == 0) cout << "Event number: " << entry << endl;
[92237ca]1335
[341014c]1336 if(branchJet->GetEntriesFast() > 1)
[92237ca]1337 {
1338
[341014c]1339 jet = (Jet *)branchJet->At(0);
[92237ca]1340 p1 = jet->P4();
[341014c]1341 jet = (Jet *)branchJet->At(1);
[92237ca]1342 p2 = jet->P4();
1343
[341014c]1344 met = (MissingET *)branchMet->At(0);
1345 scalarHT = (ScalarHT *)branchScalarHT->At(0);
[92237ca]1346 ht = scalarHT->HT;
1347
[341014c]1348 if(p1.Pt() < 0.75 * ht / 2) continue;
1349 if(p2.Pt() < 0.75 * ht / 2) continue;
[92237ca]1350
[341014c]1351 for(bin = 0; bin < Nbins; bin++)
[92237ca]1352 {
[341014c]1353 if(ht > histos->at(bin).ptmin && ht < histos->at(bin).ptmax)
[92237ca]1354 {
1355 histos->at(bin).resolHist->Fill(met->P4().Px());
1356 }
1357 }
1358 }
1359 }
1360}
1361
1362//------------------------------------------------------------------------------
1363
1364void addResoGraph(TMultiGraph *mg, TGraphErrors *gr, TLegend *leg, int style, Color_t color, TString text)
1365{
1366
1367 gr->SetLineWidth(2);
1368 gr->SetLineColor(color);
1369 gr->SetMarkerStyle(style);
1370 gr->SetMarkerColor(color);
1371 gr->SetMarkerSize(1.5);
1372
1373 std::cout << "Adding " << gr->GetName() << std::endl;
1374 mg->Add(gr);
[341014c]1375 leg->AddEntry(gr, text, "p");
[92237ca]1376}
1377
1378void DrawAxis(TMultiGraph *mg, TLegend *leg, double xmin, double xmax, double ymin, double ymax, TString tx, TString ty, bool logx = 0, bool logy = 0)
1379{
1380 mg->SetMinimum(ymin);
1381 mg->SetMaximum(ymax);
[341014c]1382 mg->GetXaxis()->SetLimits(xmin, xmax);
[92237ca]1383
1384 mg->GetXaxis()->SetTitle(tx);
1385 mg->GetYaxis()->SetTitle(ty);
1386
1387 mg->GetYaxis()->SetTitleSize(0.07);
1388 mg->GetXaxis()->SetTitleSize(0.07);
1389 mg->GetYaxis()->SetLabelSize(0.06);
1390 mg->GetXaxis()->SetLabelSize(0.06);
1391 mg->GetYaxis()->SetLabelOffset(0.03);
1392 mg->GetYaxis()->SetTitleOffset(1.4);
1393 mg->GetXaxis()->SetTitleOffset(1.4);
1394
1395 mg->GetXaxis()->SetTitleFont(132);
1396 mg->GetYaxis()->SetTitleFont(132);
1397 mg->GetXaxis()->SetLabelFont(132);
1398 mg->GetYaxis()->SetLabelFont(132);
1399
1400 mg->GetYaxis()->SetNdivisions(505);
1401
1402 leg->SetTextFont(132);
1403 leg->SetBorderSize(0);
1404 leg->SetShadowColor(0);
1405 leg->SetFillColor(0);
1406 leg->SetFillStyle(0);
1407
1408 gStyle->SetOptTitle(0);
1409
1410 if(logx) gPad->SetLogx();
1411 if(logy) gPad->SetLogy();
1412
1413 //gPad->SetGridx();
1414 //gPad->SetGridy();
1415 gPad->SetBottomMargin(0.2);
1416 gPad->SetLeftMargin(0.2);
1417 gPad->Modified();
1418 gPad->Update();
1419}
1420
[181c061]1421void DelphesValidation(
1422 const char *inputFilePion,
1423 const char *inputFileElectron,
1424 const char *inputFileMuon,
1425 const char *inputFilePhoton,
1426 const char *inputFileNeutralHadron,
1427 const char *inputFileJet,
1428 const char *inputFileBJet,
1429 const char *inputFileCJet,
1430 const char *inputFileTauJet,
1431 const char *outputFile,
1432 const char *version)
[92237ca]1433{
1434
1435 TChain *chainPion = new TChain("Delphes");
1436 chainPion->Add(inputFilePion);
1437 ExRootTreeReader *treeReaderPion = new ExRootTreeReader(chainPion);
1438
1439 TChain *chainElectron = new TChain("Delphes");
1440 chainElectron->Add(inputFileElectron);
1441 ExRootTreeReader *treeReaderElectron = new ExRootTreeReader(chainElectron);
1442
1443 TChain *chainMuon = new TChain("Delphes");
1444 chainMuon->Add(inputFileMuon);
1445 ExRootTreeReader *treeReaderMuon = new ExRootTreeReader(chainMuon);
1446
1447 TChain *chainPhoton = new TChain("Delphes");
1448 chainPhoton->Add(inputFilePhoton);
1449 ExRootTreeReader *treeReaderPhoton = new ExRootTreeReader(chainPhoton);
1450
1451 TChain *chainNeutralHadron = new TChain("Delphes");
1452 chainNeutralHadron->Add(inputFileNeutralHadron);
1453 ExRootTreeReader *treeReaderNeutralHadron = new ExRootTreeReader(chainNeutralHadron);
1454
1455 TChain *chainJet = new TChain("Delphes");
1456 chainJet->Add(inputFileJet);
1457 ExRootTreeReader *treeReaderJet = new ExRootTreeReader(chainJet);
1458
1459 TChain *chainBJet = new TChain("Delphes");
1460 chainBJet->Add(inputFileBJet);
1461 ExRootTreeReader *treeReaderBJet = new ExRootTreeReader(chainBJet);
1462
1463 TChain *chainCJet = new TChain("Delphes");
1464 chainCJet->Add(inputFileCJet);
1465 ExRootTreeReader *treeReaderCJet = new ExRootTreeReader(chainCJet);
1466
1467 TChain *chainTauJet = new TChain("Delphes");
1468 chainTauJet->Add(inputFileTauJet);
1469 ExRootTreeReader *treeReaderTauJet = new ExRootTreeReader(chainTauJet);
1470
1471 TClonesArray *branchParticleElectron = treeReaderElectron->UseBranch("Particle");
1472 TClonesArray *branchTrackElectron = treeReaderElectron->UseBranch("Track");
1473 TClonesArray *branchElectron = treeReaderElectron->UseBranch("Electron");
[b4ec6ac]1474 TClonesArray *branchElectronPF = treeReaderElectron->UseBranch("ElectronPF");
[92237ca]1475
1476 TClonesArray *branchParticleMuon = treeReaderMuon->UseBranch("Particle");
1477 TClonesArray *branchTrackMuon = treeReaderMuon->UseBranch("Track");
1478 TClonesArray *branchMuon = treeReaderMuon->UseBranch("Muon");
1479
1480 TClonesArray *branchParticlePion = treeReaderPion->UseBranch("Particle");
1481 TClonesArray *branchTrackPion = treeReaderPion->UseBranch("Track");
1482 TClonesArray *branchPion = treeReaderPion->UseBranch("Pion");
1483
1484 TClonesArray *branchParticlePhoton = treeReaderPhoton->UseBranch("Particle");
1485 TClonesArray *branchTowerPhoton = treeReaderPhoton->UseBranch("Tower");
1486 TClonesArray *branchPhoton = treeReaderPhoton->UseBranch("Photon");
1487
1488 TClonesArray *branchParticleNeutralHadron = treeReaderNeutralHadron->UseBranch("Particle");
1489 TClonesArray *branchTowerNeutralHadron = treeReaderNeutralHadron->UseBranch("Tower");
1490
1491 TClonesArray *branchGenJet = treeReaderJet->UseBranch("GenJet");
1492 TClonesArray *branchParticleJet = treeReaderJet->UseBranch("Particle");
[13cb1b6]1493 TClonesArray *branchPFJet = treeReaderJet->UseBranch("PFJet");
[92237ca]1494 TClonesArray *branchCaloJet = treeReaderJet->UseBranch("CaloJet");
[06a83b3]1495 TClonesArray *branchJet = treeReaderJet->UseBranch("Jet");
[cb159ed]1496
[92237ca]1497 TClonesArray *branchParticleBJet = treeReaderBJet->UseBranch("Particle");
1498 TClonesArray *branchPFBJet = treeReaderBJet->UseBranch("Jet");
1499
1500 TClonesArray *branchParticleCJet = treeReaderCJet->UseBranch("Particle");
1501 TClonesArray *branchPFCJet = treeReaderCJet->UseBranch("Jet");
1502
1503 TClonesArray *branchParticleTauJet = treeReaderTauJet->UseBranch("Particle");
1504 TClonesArray *branchPFTauJet = treeReaderTauJet->UseBranch("Jet");
1505
1506 TClonesArray *branchGenScalarHT = treeReaderJet->UseBranch("GenScalarHT");
[13cb1b6]1507 TClonesArray *branchMet = treeReaderJet->UseBranch("PFMissingET");
[92237ca]1508 TClonesArray *branchCaloMet = treeReaderJet->UseBranch("CaloMissingET");
1509
1510 std::vector<Color_t> colors;
1511
1512 colors.push_back(kBlack);
1513 colors.push_back(kBlue);
1514 colors.push_back(kRed);
[341014c]1515 colors.push_back(kGreen + 1);
1516 colors.push_back(kMagenta + 1);
[92237ca]1517 colors.push_back(kOrange);
1518
1519 std::vector<Int_t> markerStyles;
1520
1521 markerStyles.push_back(20);
1522 markerStyles.push_back(21);
1523 markerStyles.push_back(22);
1524 markerStyles.push_back(23);
1525 markerStyles.push_back(33);
1526 markerStyles.push_back(34);
1527
1528 TString pdfOutput(outputFile);
1529 pdfOutput.ReplaceAll(".root", ".pdf");
1530
1531 TString figPath = inputFilePion;
1532 figPath.ReplaceAll(".root", "");
1533 figPath.ReplaceAll("root", "www/fig");
1534 Int_t lastSlash = figPath.Last('/');
1535 Int_t sizePath = figPath.Length();
[341014c]1536 figPath.Remove(lastSlash + 1, sizePath);
[92237ca]1537
[85ad2b9]1538 TString header = pdfOutput;
1539 header.ReplaceAll(".pdf", "");
1540 header.ReplaceAll("validation_", "");
1541 lastSlash = header.Last('/');
1542 sizePath = header.Length();
[341014c]1543 header.Remove(0, lastSlash + 1);
[85ad2b9]1544
1545 TString vrs(version);
1546
[341014c]1547 TPaveText *pave = new TPaveText(0.0, 0.89, 0.94, 0.94, "NDC");
[bba3653]1548 pave->SetTextAlign(30);
[85ad2b9]1549 pave->SetTextFont(132);
1550 pave->SetBorderSize(0);
1551 pave->SetShadowColor(0);
1552 pave->SetFillColor(0);
1553 pave->SetFillStyle(0);
[341014c]1554 pave->AddText("Delphes " + vrs + " - " + header);
[85ad2b9]1555
[92237ca]1556 TString s_etaMin, s_etaMax, s_eta, s_pt, s_e;
1557
1558 Double_t ptMin = 1.;
1559 Double_t ptMax = 50000.;
1560
1561 Double_t etaMin = -6.;
1562 Double_t etaMax = 6.;
1563
1564 std::vector<Double_t> ptVals;
1565
1566 ptVals.push_back(1.);
1567 ptVals.push_back(10.);
1568 ptVals.push_back(100.);
1569 ptVals.push_back(1000.);
1570 ptVals.push_back(10000.);
1571
1572 std::vector<Double_t> etaVals;
1573
1574 etaVals.push_back(0.);
1575 etaVals.push_back(1.5);
1576 etaVals.push_back(2.5);
1577 etaVals.push_back(4.0);
1578 etaVals.push_back(6.0);
1579
[341014c]1580 const int n_etabins = etaVals.size() - 1;
[92237ca]1581 const int n_ptbins = ptVals.size();
1582
1583 //////////////////////////
1584 // Tracking performance //
1585 //////////////////////////
1586
1587 // --------- Pion Tracks --------- //
1588
[341014c]1589 TMultiGraph *mg_trkpi_res_pt = new TMultiGraph("", "");
1590 TMultiGraph *mg_trkpi_eff_pt = new TMultiGraph("", "");
1591 TMultiGraph *mg_trkpi_res_eta = new TMultiGraph("", "");
1592 TMultiGraph *mg_trkpi_eff_eta = new TMultiGraph("", "");
[92237ca]1593
[341014c]1594 TLegend *leg_trkpi_res_pt = new TLegend(0.55, 0.22, 0.90, 0.48);
1595 TLegend *leg_trkpi_eff_pt = (TLegend *)leg_trkpi_res_pt->Clone();
1596 TLegend *leg_trkpi_res_eta = (TLegend *)leg_trkpi_res_pt->Clone();
1597 TLegend *leg_trkpi_eff_eta = (TLegend *)leg_trkpi_res_eta->Clone();
[92237ca]1598
[341014c]1599 TGraphErrors *gr_trkpi_res_pt = new TGraphErrors[n_etabins];
1600 TGraphErrors *gr_trkpi_eff_pt = new TGraphErrors[n_etabins];
[dd5e213]1601 TGraphErrors *gr_trkpi_res_eta = new TGraphErrors[n_ptbins];
1602 TGraphErrors *gr_trkpi_eff_eta = new TGraphErrors[n_ptbins];
[341014c]1603 TH1D *h_trkpi_eff_pt, *h_trkpi_eff_eta;
[92237ca]1604
[341014c]1605 std::vector<resolPlot> *plots_trkpi_res_pt = new std::vector<resolPlot>[n_etabins];
1606 std::vector<resolPlot> *plots_trkpi_res_eta = new std::vector<resolPlot>[n_ptbins];
[92237ca]1607
1608 // loop over eta bins
[341014c]1609 for(k = 0; k < etaVals.size() - 1; k++)
[92237ca]1610 {
[341014c]1611 HistogramsCollection(&plots_trkpi_res_pt[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "trkpi");
1612 GetPtres<Track>(&plots_trkpi_res_pt[k], branchTrackPion, branchParticlePion, 211, etaVals.at(k), etaVals.at(k + 1), treeReaderPion);
1613 gr_trkpi_res_pt[k] = EresGraph(&plots_trkpi_res_pt[k]);
[92237ca]1614
[341014c]1615 h_trkpi_eff_pt = GetEffPt<Track>(branchTrackPion, branchParticlePion, "Pion", 211, ptMin, ptMax, etaVals.at(k), etaVals.at(k + 1), treeReaderPion);
1616 gr_trkpi_eff_pt[k] = TGraphErrors(h_trkpi_eff_pt);
[92237ca]1617
[341014c]1618 s_etaMin = Form("%.1f", etaVals.at(k));
1619 s_etaMax = Form("%.1f", etaVals.at(k + 1));
[92237ca]1620
[341014c]1621 s_eta = "#pi^{ #pm} , " + s_etaMin + " < | #eta | < " + s_etaMax;
[92237ca]1622
[341014c]1623 gr_trkpi_res_pt[k].SetName("trkRes_" + s_etaMin + "_" + s_etaMax);
1624 gr_trkpi_eff_pt[k].SetName("trkEff_" + s_etaMin + "_" + s_etaMax);
[92237ca]1625
[341014c]1626 addResoGraph(mg_trkpi_res_pt, &gr_trkpi_res_pt[k], leg_trkpi_res_pt, markerStyles.at(k), colors.at(k), s_eta);
1627 addResoGraph(mg_trkpi_eff_pt, &gr_trkpi_eff_pt[k], leg_trkpi_eff_pt, markerStyles.at(k), colors.at(k), s_eta);
[92237ca]1628 }
1629
1630 // loop over pt
[341014c]1631 for(k = 0; k < ptVals.size(); k++)
[92237ca]1632 {
[341014c]1633 HistogramsCollectionVsEta(&plots_trkpi_res_eta[k], etaMin, etaMax, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), "trkpi", 0.0, 2.0);
1634 GetPtresVsEta<Track>(&plots_trkpi_res_eta[k], branchTrackPion, branchParticlePion, 211, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), treeReaderPion);
1635 gr_trkpi_res_eta[k] = EresGraphVsEta(&plots_trkpi_res_eta[k]);
[92237ca]1636
[341014c]1637 h_trkpi_eff_eta = GetEffEta<Track>(branchTrackPion, branchParticlePion, "Pion", 211, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), etaMin, etaMax, treeReaderPion);
1638 gr_trkpi_eff_eta[k] = TGraphErrors(h_trkpi_eff_eta);
[92237ca]1639
[341014c]1640 s_pt = Form("#pi^{ #pm} , p_{T} = %.0f GeV", ptVals.at(k));
1641 if(ptVals.at(k) >= 1000.) s_pt = Form("#pi^{ #pm} , p_{T} = %.0f TeV", ptVals.at(k) / 1000.);
[92237ca]1642
[341014c]1643 addResoGraph(mg_trkpi_res_eta, &gr_trkpi_res_eta[k], leg_trkpi_res_eta, markerStyles.at(k), colors.at(k), s_pt);
1644 addResoGraph(mg_trkpi_eff_eta, &gr_trkpi_eff_eta[k], leg_trkpi_eff_eta, markerStyles.at(k), colors.at(k), s_pt);
[92237ca]1645 }
1646
[341014c]1647 TCanvas *c_trkpi_res_pt = new TCanvas("", "", 800, 600);
[92237ca]1648
1649 mg_trkpi_res_pt->Draw("APE");
1650 DrawAxis(mg_trkpi_res_pt, leg_trkpi_res_pt, ptMin, ptMax, 0.01, 100, "p_{T} [GeV]", "(track resolution in p_{T})/p_{T} (%)", true, true);
1651 leg_trkpi_res_pt->Draw();
[85ad2b9]1652 pave->Draw();
[92237ca]1653
[341014c]1654 c_trkpi_res_pt->Print(pdfOutput + "(", "pdf");
1655 c_trkpi_res_pt->Print(figPath + "img_trkpi_res_pt.pdf", "pdf");
1656 c_trkpi_res_pt->Print(figPath + "img_trkpi_res_pt.png", "png");
[92237ca]1657
[341014c]1658 TCanvas *c_trkpi_res_eta = new TCanvas("", "", 800, 600);
[92237ca]1659
1660 mg_trkpi_res_eta->Draw("APE");
1661 DrawAxis(mg_trkpi_res_eta, leg_trkpi_res_eta, etaMin, etaMax, 0.01, 100, " #eta ", "(track resolution in p_{T})/p_{T} (%)", false, true);
1662 leg_trkpi_res_eta->Draw();
[85ad2b9]1663 pave->Draw();
[92237ca]1664
[341014c]1665 c_trkpi_res_eta->Print(pdfOutput, "pdf");
1666 c_trkpi_res_eta->Print(figPath + "img_trkpi_res_eta.pdf", "pdf");
1667 c_trkpi_res_eta->Print(figPath + "img_trkpi_res_eta.png", "png");
[92237ca]1668
[341014c]1669 TCanvas *c_trkpi_eff_pt = new TCanvas("", "", 800, 600);
[92237ca]1670
1671 mg_trkpi_eff_pt->Draw("APE");
1672 DrawAxis(mg_trkpi_eff_pt, leg_trkpi_eff_pt, ptMin, ptMax, 0.0, 100, "p_{T} [GeV]", "tracking efficiency (%)", true, false);
1673 leg_trkpi_eff_pt->Draw();
[85ad2b9]1674 pave->Draw();
[92237ca]1675
[341014c]1676 c_trkpi_eff_pt->Print(pdfOutput, "pdf");
1677 c_trkpi_eff_pt->Print(figPath + "img_trkpi_eff_pt.pdf", "pdf");
1678 c_trkpi_eff_pt->Print(figPath + "img_trkpi_eff_pt.png", "png");
[92237ca]1679
[341014c]1680 TCanvas *c_trkpi_eff_eta = new TCanvas("", "", 800, 600);
[92237ca]1681
1682 mg_trkpi_eff_eta->Draw("APE");
1683 DrawAxis(mg_trkpi_eff_eta, leg_trkpi_eff_eta, etaMin, etaMax, 0.0, 100, " #eta ", "tracking efficiency (%)", false, false);
1684 leg_trkpi_eff_eta->Draw();
[85ad2b9]1685 pave->Draw();
1686
[341014c]1687 c_trkpi_eff_eta->Print(pdfOutput, "pdf");
1688 c_trkpi_eff_eta->Print(figPath + "img_trkpi_eff_eta.pdf", "pdf");
1689 c_trkpi_eff_eta->Print(figPath + "img_trkpi_eff_eta.png", "png");
[92237ca]1690
1691 // --------- Electron Tracks --------- //
1692
[341014c]1693 TMultiGraph *mg_trkele_res_pt = new TMultiGraph("", "");
1694 TMultiGraph *mg_trkele_eff_pt = new TMultiGraph("", "");
1695 TMultiGraph *mg_trkele_res_eta = new TMultiGraph("", "");
1696 TMultiGraph *mg_trkele_eff_eta = new TMultiGraph("", "");
[92237ca]1697
[341014c]1698 TLegend *leg_trkele_res_pt = new TLegend(0.55, 0.22, 0.90, 0.48);
1699 TLegend *leg_trkele_eff_pt = (TLegend *)leg_trkele_res_pt->Clone();
1700 TLegend *leg_trkele_res_eta = (TLegend *)leg_trkele_res_pt->Clone();
1701 TLegend *leg_trkele_eff_eta = (TLegend *)leg_trkele_res_eta->Clone();
[92237ca]1702
[341014c]1703 TGraphErrors *gr_trkele_res_pt = new TGraphErrors[n_etabins];
1704 TGraphErrors *gr_trkele_eff_pt = new TGraphErrors[n_etabins];
[dd5e213]1705 TGraphErrors *gr_trkele_res_eta = new TGraphErrors[n_ptbins];
1706 TGraphErrors *gr_trkele_eff_eta = new TGraphErrors[n_ptbins];
1707
[341014c]1708 TH1D *h_trkele_eff_pt, *h_trkele_eff_eta;
[92237ca]1709
[341014c]1710 std::vector<resolPlot> *plots_trkele_res_pt = new std::vector<resolPlot>[n_etabins];
1711 std::vector<resolPlot> *plots_trkele_res_eta = new std::vector<resolPlot>[n_ptbins];
[92237ca]1712
1713 // loop over eta bins
[341014c]1714 for(k = 0; k < etaVals.size() - 1; k++)
[92237ca]1715 {
[341014c]1716 HistogramsCollection(&plots_trkele_res_pt[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "trkele");
1717 GetPtres<Track>(&plots_trkele_res_pt[k], branchTrackElectron, branchParticleElectron, 11, etaVals.at(k), etaVals.at(k + 1), treeReaderElectron);
1718 gr_trkele_res_pt[k] = EresGraph(&plots_trkele_res_pt[k]);
[92237ca]1719
[341014c]1720 h_trkele_eff_pt = GetEffPt<Track>(branchTrackElectron, branchParticleElectron, "Electron", 11, ptMin, ptMax, etaVals.at(k), etaVals.at(k + 1), treeReaderElectron);
1721 gr_trkele_eff_pt[k] = TGraphErrors(h_trkele_eff_pt);
[92237ca]1722
[341014c]1723 s_etaMin = Form("%.1f", etaVals.at(k));
1724 s_etaMax = Form("%.1f", etaVals.at(k + 1));
[92237ca]1725
[341014c]1726 s_eta = "e^{ #pm} , " + s_etaMin + " < | #eta | < " + s_etaMax;
[92237ca]1727
[341014c]1728 gr_trkele_res_pt[k].SetName("trkRes_" + s_etaMin + "_" + s_etaMax);
1729 gr_trkele_eff_pt[k].SetName("trkEff_" + s_etaMin + "_" + s_etaMax);
[92237ca]1730
[341014c]1731 addResoGraph(mg_trkele_res_pt, &gr_trkele_res_pt[k], leg_trkele_res_pt, markerStyles.at(k), colors.at(k), s_eta);
1732 addResoGraph(mg_trkele_eff_pt, &gr_trkele_eff_pt[k], leg_trkele_eff_pt, markerStyles.at(k), colors.at(k), s_eta);
[92237ca]1733 }
1734
1735 // loop over pt
[341014c]1736 for(k = 0; k < ptVals.size(); k++)
[92237ca]1737 {
[341014c]1738 HistogramsCollectionVsEta(&plots_trkele_res_eta[k], etaMin, etaMax, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), "trkele", 0.0, 2.0);
1739 GetPtresVsEta<Track>(&plots_trkele_res_eta[k], branchTrackElectron, branchParticleElectron, 11, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), treeReaderElectron);
1740 gr_trkele_res_eta[k] = EresGraphVsEta(&plots_trkele_res_eta[k]);
[92237ca]1741
[341014c]1742 h_trkele_eff_eta = GetEffEta<Track>(branchTrackElectron, branchParticleElectron, "Electron", 11, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), etaMin, etaMax, treeReaderElectron);
1743 gr_trkele_eff_eta[k] = TGraphErrors(h_trkele_eff_eta);
[92237ca]1744
[341014c]1745 s_pt = Form("e^{ #pm} , p_{T} = %.0f GeV", ptVals.at(k));
1746 if(ptVals.at(k) >= 1000.) s_pt = Form("e^{ #pm} , p_{T} = %.0f TeV", ptVals.at(k) / 1000.);
[92237ca]1747
[341014c]1748 addResoGraph(mg_trkele_res_eta, &gr_trkele_res_eta[k], leg_trkele_res_eta, markerStyles.at(k), colors.at(k), s_pt);
1749 addResoGraph(mg_trkele_eff_eta, &gr_trkele_eff_eta[k], leg_trkele_eff_eta, markerStyles.at(k), colors.at(k), s_pt);
[92237ca]1750 }
1751
[341014c]1752 TCanvas *c_trkele_res_pt = new TCanvas("", "", 800, 600);
[92237ca]1753
1754 mg_trkele_res_pt->Draw("APE");
1755 DrawAxis(mg_trkele_res_pt, leg_trkele_res_pt, ptMin, ptMax, 0.01, 100, "p_{T} [GeV]", "(track resolution in p_{T})/p_{T} (%)", true, true);
1756 leg_trkele_res_pt->Draw();
[85ad2b9]1757 pave->Draw();
[92237ca]1758
[341014c]1759 c_trkele_res_pt->Print(pdfOutput, "pdf");
1760 c_trkele_res_pt->Print(figPath + "img_trkele_res_pt.pdf", "pdf");
1761 c_trkele_res_pt->Print(figPath + "img_trkele_res_pt.png", "png");
[92237ca]1762
[341014c]1763 TCanvas *c_trkele_res_eta = new TCanvas("", "", 800, 600);
[92237ca]1764
1765 mg_trkele_res_eta->Draw("APE");
1766 DrawAxis(mg_trkele_res_eta, leg_trkele_res_eta, etaMin, etaMax, 0.01, 100, " #eta ", "(track resolution in p_{T})/p_{T} (%)", false, true);
1767 leg_trkele_res_eta->Draw();
[85ad2b9]1768 pave->Draw();
[92237ca]1769
[341014c]1770 c_trkele_res_eta->Print(pdfOutput, "pdf");
1771 c_trkele_res_eta->Print(figPath + "img_trkele_res_eta.pdf", "pdf");
1772 c_trkele_res_eta->Print(figPath + "img_trkele_res_eta.png", "png");
[92237ca]1773
[341014c]1774 TCanvas *c_trkele_eff_pt = new TCanvas("", "", 800, 600);
[92237ca]1775
1776 mg_trkele_eff_pt->Draw("APE");
1777 DrawAxis(mg_trkele_eff_pt, leg_trkele_eff_pt, ptMin, ptMax, 0.0, 100, "p_{T} [GeV]", "tracking efficiency (%)", true, false);
1778 leg_trkele_eff_pt->Draw();
[85ad2b9]1779 pave->Draw();
[92237ca]1780
[341014c]1781 c_trkele_eff_pt->Print(pdfOutput, "pdf");
1782 c_trkele_eff_pt->Print(figPath + "img_trkele_eff_pt.pdf", "pdf");
1783 c_trkele_eff_pt->Print(figPath + "img_trkele_eff_pt.png", "png");
[92237ca]1784
[341014c]1785 TCanvas *c_trkele_eff_eta = new TCanvas("", "", 800, 600);
[92237ca]1786
1787 mg_trkele_eff_eta->Draw("APE");
1788 DrawAxis(mg_trkele_eff_eta, leg_trkele_eff_eta, etaMin, etaMax, 0.0, 100, " #eta ", "tracking efficiency (%)", false, false);
1789 leg_trkele_eff_eta->Draw();
[85ad2b9]1790 pave->Draw();
[92237ca]1791
[341014c]1792 c_trkele_eff_eta->Print(pdfOutput, "pdf");
1793 c_trkele_eff_eta->Print(figPath + "img_trkele_eff_eta.pdf", "pdf");
1794 c_trkele_eff_eta->Print(figPath + "img_trkele_eff_eta.png", "png");
[92237ca]1795
1796 // --------- Muon Tracks --------- //
1797
[341014c]1798 TMultiGraph *mg_trkmu_res_pt = new TMultiGraph("", "");
1799 TMultiGraph *mg_trkmu_eff_pt = new TMultiGraph("", "");
1800 TMultiGraph *mg_trkmu_res_eta = new TMultiGraph("", "");
1801 TMultiGraph *mg_trkmu_eff_eta = new TMultiGraph("", "");
[92237ca]1802
[341014c]1803 TLegend *leg_trkmu_res_pt = new TLegend(0.55, 0.22, 0.90, 0.48);
1804 TLegend *leg_trkmu_eff_pt = (TLegend *)leg_trkmu_res_pt->Clone();
[92237ca]1805
[341014c]1806 TLegend *leg_trkmu_res_eta = (TLegend *)leg_trkmu_res_pt->Clone();
1807 TLegend *leg_trkmu_eff_eta = (TLegend *)leg_trkmu_res_eta->Clone();
[92237ca]1808
[341014c]1809 TGraphErrors *gr_trkmu_res_pt = new TGraphErrors[n_etabins];
1810 TGraphErrors *gr_trkmu_eff_pt = new TGraphErrors[n_etabins];
[dd5e213]1811 TGraphErrors *gr_trkmu_res_eta = new TGraphErrors[n_ptbins];
1812 TGraphErrors *gr_trkmu_eff_eta = new TGraphErrors[n_ptbins];
1813
[341014c]1814 TH1D *h_trkmu_eff_pt, *h_trkmu_eff_eta;
[92237ca]1815
[341014c]1816 std::vector<resolPlot> *plots_trkmu_res_pt = new std::vector<resolPlot>[n_etabins];
1817 std::vector<resolPlot> *plots_trkmu_res_eta = new std::vector<resolPlot>[n_ptbins];
[92237ca]1818
1819 // loop over eta bins
[341014c]1820 for(k = 0; k < etaVals.size() - 1; k++)
[92237ca]1821 {
[341014c]1822 HistogramsCollection(&plots_trkmu_res_pt[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "trkmu");
1823 GetPtres<Track>(&plots_trkmu_res_pt[k], branchTrackMuon, branchParticleMuon, 13, etaVals.at(k), etaVals.at(k + 1), treeReaderMuon);
1824 gr_trkmu_res_pt[k] = EresGraph(&plots_trkmu_res_pt[k]);
[92237ca]1825
[341014c]1826 h_trkmu_eff_pt = GetEffPt<Track>(branchTrackMuon, branchParticleMuon, "Muon", 13, ptMin, ptMax, etaVals.at(k), etaVals.at(k + 1), treeReaderMuon);
1827 gr_trkmu_eff_pt[k] = TGraphErrors(h_trkmu_eff_pt);
[92237ca]1828
[341014c]1829 s_etaMin = Form("%.1f", etaVals.at(k));
1830 s_etaMax = Form("%.1f", etaVals.at(k + 1));
[92237ca]1831
[341014c]1832 s_eta = "#mu^{ #pm} , " + s_etaMin + " < | #eta | < " + s_etaMax;
[92237ca]1833
[341014c]1834 gr_trkmu_res_pt[k].SetName("trkRes_" + s_etaMin + "_" + s_etaMax);
1835 gr_trkmu_eff_pt[k].SetName("trkEff_" + s_etaMin + "_" + s_etaMax);
[92237ca]1836
[341014c]1837 addResoGraph(mg_trkmu_res_pt, &gr_trkmu_res_pt[k], leg_trkmu_res_pt, markerStyles.at(k), colors.at(k), s_eta);
1838 addResoGraph(mg_trkmu_eff_pt, &gr_trkmu_eff_pt[k], leg_trkmu_eff_pt, markerStyles.at(k), colors.at(k), s_eta);
[92237ca]1839 }
1840
1841 // loop over pt
[341014c]1842 for(k = 0; k < ptVals.size(); k++)
[92237ca]1843 {
[341014c]1844 HistogramsCollectionVsEta(&plots_trkmu_res_eta[k], etaMin, etaMax, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), "trkmu", 0.0, 2.0);
1845 GetPtresVsEta<Track>(&plots_trkmu_res_eta[k], branchTrackMuon, branchParticleMuon, 13, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), treeReaderMuon);
1846 gr_trkmu_res_eta[k] = EresGraphVsEta(&plots_trkmu_res_eta[k]);
[92237ca]1847
[341014c]1848 h_trkmu_eff_eta = GetEffEta<Track>(branchTrackMuon, branchParticleMuon, "Muon", 13, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), etaMin, etaMax, treeReaderMuon);
1849 gr_trkmu_eff_eta[k] = TGraphErrors(h_trkmu_eff_eta);
[92237ca]1850
[341014c]1851 s_pt = Form("#mu^{ #pm} , p_{T} = %.0f GeV", ptVals.at(k));
1852 if(ptVals.at(k) >= 1000.) s_pt = Form("#mu^{ #pm} , p_{T} = %.0f TeV", ptVals.at(k) / 1000.);
[92237ca]1853
[341014c]1854 addResoGraph(mg_trkmu_res_eta, &gr_trkmu_res_eta[k], leg_trkmu_res_eta, markerStyles.at(k), colors.at(k), s_pt);
1855 addResoGraph(mg_trkmu_eff_eta, &gr_trkmu_eff_eta[k], leg_trkmu_eff_eta, markerStyles.at(k), colors.at(k), s_pt);
[92237ca]1856 }
1857
[341014c]1858 TCanvas *c_trkmu_res_pt = new TCanvas("", "", 800, 600);
[92237ca]1859
1860 mg_trkmu_res_pt->Draw("APE");
1861 DrawAxis(mg_trkmu_res_pt, leg_trkmu_res_pt, ptMin, ptMax, 0.01, 100, "p_{T} [GeV]", "(track resolution in p_{T})/p_{T} (%)", true, true);
1862 leg_trkmu_res_pt->Draw();
[85ad2b9]1863 pave->Draw();
[92237ca]1864
[341014c]1865 c_trkmu_res_pt->Print(pdfOutput, "pdf");
1866 c_trkmu_res_pt->Print(figPath + "img_trkmu_res_pt.pdf", "pdf");
1867 c_trkmu_res_pt->Print(figPath + "img_trkmu_res_pt.png", "png");
[92237ca]1868
[341014c]1869 TCanvas *c_trkmu_res_eta = new TCanvas("", "", 800, 600);
[92237ca]1870
1871 mg_trkmu_res_eta->Draw("APE");
1872 DrawAxis(mg_trkmu_res_eta, leg_trkmu_res_eta, etaMin, etaMax, 0.01, 100, " #eta ", "(track resolution in p_{T})/p_{T} (%)", false, true);
1873 leg_trkmu_res_eta->Draw();
[85ad2b9]1874 pave->Draw();
[92237ca]1875
[341014c]1876 c_trkmu_res_eta->Print(pdfOutput, "pdf");
1877 c_trkmu_res_eta->Print(figPath + "img_trkmu_res_eta.pdf", "pdf");
1878 c_trkmu_res_eta->Print(figPath + "img_trkmu_res_eta.png", "png");
[92237ca]1879
[341014c]1880 TCanvas *c_trkmu_eff_pt = new TCanvas("", "", 800, 600);
[92237ca]1881
1882 mg_trkmu_eff_pt->Draw("APE");
1883 DrawAxis(mg_trkmu_eff_pt, leg_trkmu_eff_pt, ptMin, ptMax, 0.0, 100, "p_{T} [GeV]", "tracking efficiency (%)", true, false);
1884 leg_trkmu_eff_pt->Draw();
[85ad2b9]1885 pave->Draw();
[92237ca]1886
[341014c]1887 c_trkmu_eff_pt->Print(pdfOutput, "pdf");
1888 c_trkmu_eff_pt->Print(figPath + "img_trkmu_eff_pt.pdf", "pdf");
1889 c_trkmu_eff_pt->Print(figPath + "img_trkmu_eff_pt.png", "png");
[92237ca]1890
[341014c]1891 TCanvas *c_trkmu_eff_eta = new TCanvas("", "", 800, 600);
[92237ca]1892
1893 mg_trkmu_eff_eta->Draw("APE");
1894 DrawAxis(mg_trkmu_eff_eta, leg_trkmu_eff_eta, etaMin, etaMax, 0.0, 100, " #eta ", "tracking efficiency (%)", false, false);
1895 leg_trkmu_eff_eta->Draw();
[85ad2b9]1896 pave->Draw();
[92237ca]1897
[341014c]1898 c_trkmu_eff_eta->Print(pdfOutput, "pdf");
1899 c_trkmu_eff_eta->Print(figPath + "img_trkmu_eff_eta.pdf", "pdf");
1900 c_trkmu_eff_eta->Print(figPath + "img_trkmu_eff_eta.png", "png");
[92237ca]1901
1902 //////////////////////
1903 // Ecal performance //
1904 //////////////////////
1905
[341014c]1906 TMultiGraph *mg_ecal_res_e = new TMultiGraph("", "");
1907 TMultiGraph *mg_ecal_res_eta = new TMultiGraph("", "");
[92237ca]1908
[341014c]1909 TLegend *leg_ecal_res_e = new TLegend(0.55, 0.64, 0.90, 0.90);
1910 TLegend *leg_ecal_res_eta = new TLegend(0.60, 0.59, 0.95, 0.90);
[92237ca]1911
[341014c]1912 TGraphErrors *gr_ecal_res_e = new TGraphErrors[n_etabins];
1913 TGraphErrors *gr_ecal_res_eta = new TGraphErrors[n_ptbins];
[92237ca]1914
[341014c]1915 std::vector<resolPlot> *plots_ecal_res_e = new std::vector<resolPlot>[n_etabins];
[dd5e213]1916 std::vector<resolPlot> *plots_ecal_res_eta = new std::vector<resolPlot>[n_ptbins];
[92237ca]1917
1918 // loop over eta bins
[341014c]1919 for(k = 0; k < etaVals.size() - 1; k++)
[92237ca]1920 {
[341014c]1921 HistogramsCollection(&plots_ecal_res_e[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "ecal");
1922 GetEres<Tower>(&plots_ecal_res_e[k], branchTowerPhoton, branchParticlePhoton, 22, etaVals.at(k), etaVals.at(k + 1), treeReaderPhoton);
1923 gr_ecal_res_e[k] = EresGraph(&plots_ecal_res_e[k]);
[92237ca]1924
[341014c]1925 s_etaMin = Form("%.1f", etaVals.at(k));
1926 s_etaMax = Form("%.1f", etaVals.at(k + 1));
[92237ca]1927
[341014c]1928 s_eta = "#gamma , " + s_etaMin + " < | #eta | < " + s_etaMax;
[92237ca]1929
[341014c]1930 gr_ecal_res_e[k].SetName("trkRes_" + s_etaMin + "_" + s_etaMax);
[92237ca]1931
[341014c]1932 addResoGraph(mg_ecal_res_e, &gr_ecal_res_e[k], leg_ecal_res_e, markerStyles.at(k), colors.at(k), s_eta);
[92237ca]1933 }
1934
1935 // loop over pt
[341014c]1936 for(k = 0; k < ptVals.size(); k++)
[92237ca]1937 {
[341014c]1938 HistogramsCollectionVsEta(&plots_ecal_res_eta[k], etaMin, etaMax, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), "ecal", 0.0, 2.0);
1939 GetEresVsEta<Tower>(&plots_ecal_res_eta[k], branchTowerPhoton, branchParticlePhoton, 22, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), treeReaderPhoton);
1940 gr_ecal_res_eta[k] = EresGraphVsEta(&plots_ecal_res_eta[k]);
[92237ca]1941
[341014c]1942 s_e = Form("#gamma , E = %.0f GeV", ptVals.at(k));
1943 if(ptVals.at(k) >= 1000.) s_e = Form("#gamma , E = %.0f TeV", ptVals.at(k) / 1000.);
[92237ca]1944
[341014c]1945 addResoGraph(mg_ecal_res_eta, &gr_ecal_res_eta[k], leg_ecal_res_eta, markerStyles.at(k), colors.at(k), s_e);
[92237ca]1946 }
1947
[341014c]1948 TCanvas *c_ecal_res_e = new TCanvas("", "", 800, 600);
[92237ca]1949
1950 mg_ecal_res_e->Draw("APE");
[341014c]1951 // DrawAxis(mg_ecal_res_e, leg_ecal_res_e, ptMin, ptMax, 0.5, 100, "E [GeV]", "(ECAL resolution in E)/E (%)", true, true);
[92237ca]1952 DrawAxis(mg_ecal_res_e, leg_ecal_res_e, ptMin, ptMax, 0.0, 20, "E [GeV]", "(ECAL resolution in E)/E (%)", true, false);
1953 leg_ecal_res_e->Draw();
[85ad2b9]1954 pave->Draw();
[92237ca]1955
[341014c]1956 c_ecal_res_e->Print(pdfOutput, "pdf");
1957 c_ecal_res_e->Print(figPath + "img_ecal_res_e.pdf", "pdf");
1958 c_ecal_res_e->Print(figPath + "img_ecal_res_e.png", "png");
[92237ca]1959
[341014c]1960 TCanvas *c_ecal_res_eta = new TCanvas("", "", 800, 600);
[92237ca]1961
1962 mg_ecal_res_eta->Draw("APE");
1963 //DrawAxis(mg_ecal_res_eta, leg_ecal_res_eta, etaMin, etaMax, 0.5, 100, " #eta ", "(ECAL resolution in E)/E (%)", false, true);
1964 DrawAxis(mg_ecal_res_eta, leg_ecal_res_eta, etaMin, etaMax, 0.0, 20, " #eta ", "(ECAL resolution in E)/E (%)", false, false);
1965 leg_ecal_res_eta->Draw();
[85ad2b9]1966 pave->Draw();
[92237ca]1967
[341014c]1968 c_ecal_res_eta->Print(pdfOutput, "pdf");
1969 c_ecal_res_eta->Print(figPath + "img_ecal_res_eta.pdf", "pdf");
1970 c_ecal_res_eta->Print(figPath + "img_ecal_res_eta.png", "png");
[92237ca]1971
1972 //////////////////////
1973 // Hcal performance //
1974 //////////////////////
1975
[341014c]1976 TMultiGraph *mg_hcal_res_e = new TMultiGraph("", "");
1977 TMultiGraph *mg_hcal_res_eta = new TMultiGraph("", "");
[92237ca]1978
[341014c]1979 TLegend *leg_hcal_res_e = new TLegend(0.55, 0.64, 0.90, 0.90);
1980 TLegend *leg_hcal_res_eta = new TLegend(0.60, 0.59, 0.95, 0.90);
[92237ca]1981
[341014c]1982 TGraphErrors *gr_hcal_res_e = new TGraphErrors[n_etabins];
[dd5e213]1983 TGraphErrors *gr_hcal_res_eta = new TGraphErrors[n_ptbins];
[92237ca]1984
[341014c]1985 std::vector<resolPlot> *plots_hcal_res_e = new std::vector<resolPlot>[n_etabins];
[dd5e213]1986 std::vector<resolPlot> *plots_hcal_res_eta = new std::vector<resolPlot>[n_ptbins];
[92237ca]1987
1988 // loop over eta bins
[341014c]1989 for(k = 0; k < etaVals.size() - 1; k++)
[92237ca]1990 {
[341014c]1991 HistogramsCollection(&plots_hcal_res_e[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "hcal");
1992 GetEres<Tower>(&plots_hcal_res_e[k], branchTowerNeutralHadron, branchParticleNeutralHadron, 2112, etaVals.at(k), etaVals.at(k + 1), treeReaderNeutralHadron);
[92237ca]1993
[341014c]1994 gr_hcal_res_e[k] = EresGraph(&plots_hcal_res_e[k]);
[92237ca]1995
[341014c]1996 s_etaMin = Form("%.1f", etaVals.at(k));
1997 s_etaMax = Form("%.1f", etaVals.at(k + 1));
[92237ca]1998
[341014c]1999 s_eta = "n , " + s_etaMin + " < | #eta | < " + s_etaMax;
[92237ca]2000
[341014c]2001 gr_hcal_res_e[k].SetName("trkRes_" + s_etaMin + "_" + s_etaMax);
[92237ca]2002
[341014c]2003 addResoGraph(mg_hcal_res_e, &gr_hcal_res_e[k], leg_hcal_res_e, markerStyles.at(k), colors.at(k), s_eta);
[92237ca]2004 }
2005
2006 // loop over pt
[341014c]2007 for(k = 0; k < ptVals.size(); k++)
[92237ca]2008 {
[341014c]2009 HistogramsCollectionVsEta(&plots_hcal_res_eta[k], etaMin, etaMax, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), "hcal", 0.0, 2.0);
2010 GetEresVsEta<Tower>(&plots_hcal_res_eta[k], branchTowerNeutralHadron, branchParticleNeutralHadron, 2112, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), treeReaderNeutralHadron);
2011 gr_hcal_res_eta[k] = EresGraphVsEta(&plots_hcal_res_eta[k]);
[92237ca]2012
[341014c]2013 s_e = Form("n , E = %.0f GeV", ptVals.at(k));
2014 if(ptVals.at(k) >= 1000.) s_e = Form("n , E = %.0f TeV", ptVals.at(k) / 1000.);
[92237ca]2015
[341014c]2016 addResoGraph(mg_hcal_res_eta, &gr_hcal_res_eta[k], leg_hcal_res_eta, markerStyles.at(k), colors.at(k), s_e);
[92237ca]2017 }
2018
[341014c]2019 TCanvas *c_hcal_res_e = new TCanvas("", "", 800, 600);
[92237ca]2020
2021 mg_hcal_res_e->Draw("APE");
2022 //DrawAxis(mg_hcal_res_e, leg_hcal_res_e, ptMin, ptMax, 1, 100, "E [GeV]", "(HCAL resolution in E)/E (%)", true, true);
2023 DrawAxis(mg_hcal_res_e, leg_hcal_res_e, ptMin, ptMax, 0.0, 50, "E [GeV]", "(HCAL resolution in E)/E (%)", true, false);
2024 leg_hcal_res_e->Draw();
[85ad2b9]2025 pave->Draw();
[92237ca]2026
[341014c]2027 c_hcal_res_e->Print(pdfOutput, "pdf");
2028 c_hcal_res_e->Print(figPath + "img_hcal_res_e.pdf", "pdf");
2029 c_hcal_res_e->Print(figPath + "img_hcal_res_e.png", "png");
[92237ca]2030
[341014c]2031 TCanvas *c_hcal_res_eta = new TCanvas("", "", 800, 600);
[92237ca]2032
2033 mg_hcal_res_eta->Draw("APE");
2034 //DrawAxis(mg_hcal_res_eta, leg_hcal_res_eta, etaMin, etaMax, 1, 100, " #eta ", "(HCAL resolution in E)/E (%)", false, true);
2035 DrawAxis(mg_hcal_res_eta, leg_hcal_res_eta, etaMin, etaMax, 0.0, 50, " #eta ", "(HCAL resolution in E)/E (%)", false, false);
2036 leg_hcal_res_eta->Draw();
[85ad2b9]2037 pave->Draw();
[92237ca]2038
[341014c]2039 c_hcal_res_eta->Print(pdfOutput, "pdf");
2040 c_hcal_res_eta->Print(figPath + "img_hcal_res_eta.pdf", "pdf");
2041 c_hcal_res_eta->Print(figPath + "img_hcal_res_eta.png", "png");
[92237ca]2042
2043 ////////////////////
2044 // PF - electrons //
2045 ////////////////////
2046
2047 TMultiGraph *mg_pfele_res_e[n_etabins];
2048 TMultiGraph *mg_pfele_res_eta[n_ptbins];
2049
2050 TLegend *leg_pfele_res_e[n_etabins];
2051 TLegend *leg_pfele_res_eta[n_ptbins];
2052
[dd5e213]2053 TGraphErrors *gr_pfele_res_e = new TGraphErrors[n_etabins];
2054 TGraphErrors *gr_pfele_res_eta = new TGraphErrors[n_ptbins];
2055 TGraphErrors *gr_trkele_res_e = new TGraphErrors[n_etabins];
2056 TGraphErrors *gr_trkele_res_eeta = new TGraphErrors[n_ptbins];
[cb159ed]2057
[341014c]2058 std::vector<resolPlot> *plots_pfele_res_e = new std::vector<resolPlot>[n_etabins];
2059 std::vector<resolPlot> *plots_pfele_res_eta = new std::vector<resolPlot>[n_ptbins];
2060 std::vector<resolPlot> *plots_trkele_res_e = new std::vector<resolPlot>[n_etabins];
[dd5e213]2061 std::vector<resolPlot> *plots_trkele_res_eeta = new std::vector<resolPlot>[n_ptbins];
[92237ca]2062
2063 TCanvas *c_pfele_res_e[n_etabins];
2064 TCanvas *c_pfele_res_eta[n_ptbins];
2065
2066 // loop over eta bins
[341014c]2067 for(k = 0; k < etaVals.size() - 1; k++)
[92237ca]2068 {
[341014c]2069 mg_pfele_res_e[k] = new TMultiGraph("", "");
2070 leg_pfele_res_e[k] = new TLegend(0.40, 0.60, 0.75, 0.90);
[92237ca]2071
[341014c]2072 HistogramsCollection(&plots_pfele_res_e[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "pfele");
2073 GetEres<Electron>(&plots_pfele_res_e[k], branchElectronPF, branchParticleElectron, 11, etaVals.at(k), etaVals.at(k + 1), treeReaderElectron);
2074 gr_pfele_res_e[k] = EresGraph(&plots_pfele_res_e[k]);
[92237ca]2075
[341014c]2076 HistogramsCollection(&plots_trkele_res_e[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "trkele");
2077 GetEres<Track>(&plots_trkele_res_e[k], branchTrackElectron, branchParticleElectron, 11, etaVals.at(k), etaVals.at(k + 1), treeReaderElectron);
2078 gr_trkele_res_e[k] = EresGraph(&plots_trkele_res_e[k]);
[cb159ed]2079
[341014c]2080 s_etaMin = Form("%.1f", etaVals.at(k));
2081 s_etaMax = Form("%.1f", etaVals.at(k + 1));
2082 s_eta = "e^{ #pm}, " + s_etaMin + " < | #eta | < " + s_etaMax;
[92237ca]2083
[341014c]2084 leg_pfele_res_e[k]->SetTextFont(132);
2085 leg_pfele_res_e[k]->SetHeader(s_eta);
[92237ca]2086
[341014c]2087 addResoGraph(mg_pfele_res_e[k], &gr_ecal_res_e[k], leg_pfele_res_e[k], markerStyles.at(0), colors.at(0), "ECAL");
2088 addResoGraph(mg_pfele_res_e[k], &gr_trkele_res_e[k], leg_pfele_res_e[k], markerStyles.at(1), colors.at(1), "Track");
2089 addResoGraph(mg_pfele_res_e[k], &gr_pfele_res_e[k], leg_pfele_res_e[k], markerStyles.at(2), colors.at(2), "Particle-flow");
[92237ca]2090
[341014c]2091 c_pfele_res_e[k] = new TCanvas("", "", 800, 600);
[92237ca]2092
[341014c]2093 mg_pfele_res_e[k]->Draw("APE");
2094 //DrawAxis(mg_pfele_res_e[k], leg_pfele_res_e[k], ptMin, ptMax, 0.1, 100, "E [GeV]", "(resolution in E)/E (%)", true, true);
2095 DrawAxis(mg_pfele_res_e[k], leg_pfele_res_e[k], ptMin, ptMax, 0.0, 20, "E [GeV]", "(resolution in E)/E (%)", true, false);
2096 leg_pfele_res_e[k]->Draw();
2097 pave->Draw();
[92237ca]2098
[341014c]2099 TString s_etarange = "eta_" + s_etaMin + "_" + s_etaMax + "_";
[92237ca]2100
[341014c]2101 c_pfele_res_e[k]->Print(pdfOutput, "pdf");
2102 c_pfele_res_e[k]->Print(figPath + "img_pfele_res_" + s_etarange + "e.pdf", "pdf");
2103 c_pfele_res_e[k]->Print(figPath + "img_pfele_res_" + s_etarange + "e.png", "png");
[92237ca]2104 }
2105
[341014c]2106 // loop over eta bins
2107 for(k = 0; k < ptVals.size(); k++)
[92237ca]2108 {
2109
[341014c]2110 mg_pfele_res_eta[k] = new TMultiGraph("", "");
2111 leg_pfele_res_eta[k] = new TLegend(0.40, 0.60, 0.75, 0.90);
[cb159ed]2112
[341014c]2113 HistogramsCollectionVsEta(&plots_pfele_res_eta[k], etaMin, etaMax, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), "pfele", 0.0, 2.0);
2114 GetEresVsEta<Electron>(&plots_pfele_res_eta[k], branchElectronPF, branchParticleElectron, 11, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), treeReaderElectron);
2115 gr_pfele_res_eta[k] = EresGraphVsEta(&plots_pfele_res_eta[k]);
[92237ca]2116
[341014c]2117 HistogramsCollectionVsEta(&plots_trkele_res_eeta[k], etaMin, etaMax, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), "trkele", 0.0, 2.0);
2118 GetEresVsEta<Track>(&plots_trkele_res_eeta[k], branchTrackElectron, branchParticleElectron, 11, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), treeReaderElectron);
2119 gr_trkele_res_eeta[k] = EresGraphVsEta(&plots_trkele_res_eeta[k]);
[92237ca]2120
[341014c]2121 s_e = Form("e^{ #pm}, E = %.0f GeV", ptVals.at(k));
2122 if(ptVals.at(k) >= 1000.) s_e = Form("e^{ #pm}, E = %.0f TeV", ptVals.at(k) / 1000.);
[92237ca]2123
[341014c]2124 leg_pfele_res_eta[k]->SetTextFont(132);
2125 leg_pfele_res_eta[k]->SetHeader(s_e);
[92237ca]2126
[341014c]2127 addResoGraph(mg_pfele_res_eta[k], &gr_ecal_res_eta[k], leg_pfele_res_eta[k], markerStyles.at(0), colors.at(0), "ECAL");
2128 addResoGraph(mg_pfele_res_eta[k], &gr_trkele_res_eeta[k], leg_pfele_res_eta[k], markerStyles.at(1), colors.at(1), "Track");
2129 addResoGraph(mg_pfele_res_eta[k], &gr_pfele_res_eta[k], leg_pfele_res_eta[k], markerStyles.at(2), colors.at(2), "Particle-flow");
[92237ca]2130
[341014c]2131 c_pfele_res_eta[k] = new TCanvas("", "", 800, 600);
[92237ca]2132
[341014c]2133 mg_pfele_res_eta[k]->Draw("APE");
2134 //DrawAxis(mg_pfele_res_eta[k], leg_pfele_res_eta[k], etaMin, etaMax, 0.1, 1000, "#eta", "(resolution in E)/E (%)", false, true);
2135 DrawAxis(mg_pfele_res_eta[k], leg_pfele_res_eta[k], etaMin, etaMax, 0.0, 50, "#eta", "(resolution in E)/E (%)", false, false);
2136 leg_pfele_res_eta[k]->Draw();
2137 pave->Draw();
[92237ca]2138
[341014c]2139 TString s_ptrange = Form("pt_%.0f_", ptVals.at(k));
[92237ca]2140
[341014c]2141 c_pfele_res_eta[k]->Print(pdfOutput, "pdf");
2142 c_pfele_res_eta[k]->Print(figPath + "img_pfele_res_" + s_ptrange + "eta.pdf", "pdf");
2143 c_pfele_res_eta[k]->Print(figPath + "img_pfele_res_" + s_ptrange + "eta.png", "png");
[92237ca]2144 }
2145
2146 /////////////////
2147 // PF - Pions //
2148 /////////////////
2149
2150 TMultiGraph *mg_pfpi_res_e[n_etabins];
2151 TMultiGraph *mg_pfpi_res_eta[n_ptbins];
2152
2153 TLegend *leg_pfpi_res_e[n_etabins];
2154 TLegend *leg_pfpi_res_eta[n_ptbins];
2155
[341014c]2156 TGraphErrors *gr_pfpi_res_e = new TGraphErrors[n_etabins];
2157 TGraphErrors *gr_pfpi_res_eta = new TGraphErrors[n_ptbins];
[92237ca]2158
[341014c]2159 TGraphErrors *gr_trkpi_res_e = new TGraphErrors[n_etabins];
[dd5e213]2160 TGraphErrors *gr_trkpi_res_eeta = new TGraphErrors[n_ptbins];
[cb159ed]2161
[341014c]2162 std::vector<resolPlot> *plots_pfpi_res_e = new std::vector<resolPlot>[n_etabins];
2163 std::vector<resolPlot> *plots_pfpi_res_eta = new std::vector<resolPlot>[n_ptbins];
2164 std::vector<resolPlot> *plots_trkpi_res_e = new std::vector<resolPlot>[n_etabins];
[dd5e213]2165 std::vector<resolPlot> *plots_trkpi_res_eeta = new std::vector<resolPlot>[n_ptbins];
[92237ca]2166
2167 TCanvas *c_pfpi_res_e[n_etabins];
2168 TCanvas *c_pfpi_res_eta[n_ptbins];
2169
2170 // loop over eta bins
[341014c]2171 for(k = 0; k < etaVals.size() - 1; k++)
[92237ca]2172 {
[341014c]2173 mg_pfpi_res_e[k] = new TMultiGraph("", "");
2174 leg_pfpi_res_e[k] = new TLegend(0.40, 0.60, 0.75, 0.90);
[cb159ed]2175
[341014c]2176 HistogramsCollection(&plots_pfpi_res_e[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "pfpi");
2177 GetEres<Track>(&plots_pfpi_res_e[k], branchPion, branchParticlePion, 211, etaVals.at(k), etaVals.at(k + 1), treeReaderPion);
2178 gr_pfpi_res_e[k] = EresGraph(&plots_pfpi_res_e[k]);
[cb159ed]2179
[341014c]2180 HistogramsCollection(&plots_trkpi_res_e[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "trkpi");
2181 GetEres<Track>(&plots_trkpi_res_e[k], branchTrackPion, branchParticlePion, 211, etaVals.at(k), etaVals.at(k + 1), treeReaderPion);
2182 gr_trkpi_res_e[k] = EresGraph(&plots_trkpi_res_e[k]);
[92237ca]2183
[341014c]2184 s_etaMin = Form("%.1f", etaVals.at(k));
2185 s_etaMax = Form("%.1f", etaVals.at(k + 1));
2186 s_eta = "#pi^{ #pm}, " + s_etaMin + " < | #eta | < " + s_etaMax;
[92237ca]2187
[341014c]2188 leg_pfpi_res_e[k]->SetTextFont(132);
2189 leg_pfpi_res_e[k]->SetHeader(s_eta);
[92237ca]2190
[341014c]2191 addResoGraph(mg_pfpi_res_e[k], &gr_hcal_res_e[k], leg_pfpi_res_e[k], markerStyles.at(0), colors.at(0), "HCAL");
2192 addResoGraph(mg_pfpi_res_e[k], &gr_trkpi_res_e[k], leg_pfpi_res_e[k], markerStyles.at(1), colors.at(1), "Track");
2193 addResoGraph(mg_pfpi_res_e[k], &gr_pfpi_res_e[k], leg_pfpi_res_e[k], markerStyles.at(2), colors.at(2), "Particle-flow");
[92237ca]2194
[341014c]2195 c_pfpi_res_e[k] = new TCanvas("", "", 800, 600);
[92237ca]2196
[341014c]2197 mg_pfpi_res_e[k]->Draw("APE");
2198 //DrawAxis(mg_pfpi_res_e[k], leg_pfpi_res_e[k], ptMin, ptMax, 0.1, 100, "E [GeV]", "(resolution in E)/E (%)", true, true);
2199 DrawAxis(mg_pfpi_res_e[k], leg_pfpi_res_e[k], ptMin, ptMax, 0.1, 50, "E [GeV]", "(resolution in E)/E (%)", true, false);
2200 leg_pfpi_res_e[k]->Draw();
2201 pave->Draw();
[92237ca]2202
[341014c]2203 TString s_etarange = "eta_" + s_etaMin + "_" + s_etaMax + "_";
[92237ca]2204
[341014c]2205 c_pfpi_res_e[k]->Print(pdfOutput, "pdf");
2206 c_pfpi_res_e[k]->Print(figPath + "img_pfpi_res_" + s_etarange + "e.pdf", "pdf");
2207 c_pfpi_res_e[k]->Print(figPath + "img_pfpi_res_" + s_etarange + "e.png", "png");
[92237ca]2208 }
2209
[341014c]2210 // loop over eta bins
2211 for(k = 0; k < ptVals.size(); k++)
[92237ca]2212 {
2213
[341014c]2214 mg_pfpi_res_eta[k] = new TMultiGraph("", "");
2215 leg_pfpi_res_eta[k] = new TLegend(0.40, 0.60, 0.75, 0.90);
[cb159ed]2216
[341014c]2217 HistogramsCollectionVsEta(&plots_pfpi_res_eta[k], etaMin, etaMax, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), "pfpi", 0.0, 2.0);
2218 GetEresVsEta<Track>(&plots_pfpi_res_eta[k], branchPion, branchParticlePion, 211, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), treeReaderPion);
2219 gr_pfpi_res_eta[k] = EresGraphVsEta(&plots_pfpi_res_eta[k]);
[cb159ed]2220
[341014c]2221 HistogramsCollectionVsEta(&plots_trkpi_res_eeta[k], etaMin, etaMax, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), "trkpi", 0.0, 2.0);
2222 GetEresVsEta<Track>(&plots_trkpi_res_eeta[k], branchPion, branchParticlePion, 211, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), treeReaderPion);
2223 gr_trkpi_res_eeta[k] = EresGraphVsEta(&plots_trkpi_res_eeta[k]);
[92237ca]2224
[341014c]2225 s_e = Form("#pi^{ #pm}, E = %.0f GeV", ptVals.at(k));
2226 if(ptVals.at(k) >= 1000.) s_e = Form("#pi^{ #pm}, E = %.0f TeV", ptVals.at(k) / 1000.);
[92237ca]2227
[341014c]2228 leg_pfpi_res_eta[k]->SetTextFont(132);
2229 leg_pfpi_res_eta[k]->SetHeader(s_e);
[92237ca]2230
[341014c]2231 addResoGraph(mg_pfpi_res_eta[k], &gr_hcal_res_eta[k], leg_pfpi_res_eta[k], markerStyles.at(0), colors.at(0), "HCAL");
2232 addResoGraph(mg_pfpi_res_eta[k], &gr_trkpi_res_eeta[k], leg_pfpi_res_eta[k], markerStyles.at(1), colors.at(1), "Track");
2233 addResoGraph(mg_pfpi_res_eta[k], &gr_pfpi_res_eta[k], leg_pfpi_res_eta[k], markerStyles.at(2), colors.at(2), "Particle-flow");
[92237ca]2234
[341014c]2235 c_pfpi_res_eta[k] = new TCanvas("", "", 800, 600);
[92237ca]2236
[341014c]2237 mg_pfpi_res_eta[k]->Draw("APE");
2238 //DrawAxis(mg_pfpi_res_eta[k], leg_pfpi_res_eta[k], etaMin, etaMax, 0.1, 1000, "#eta", "(resolution in E)/E (%)", false, true);
2239 DrawAxis(mg_pfpi_res_eta[k], leg_pfpi_res_eta[k], etaMin, etaMax, 0.0, 50, "#eta", "(resolution in E)/E (%)", false, false);
2240 leg_pfpi_res_eta[k]->Draw();
2241 pave->Draw();
[92237ca]2242
[341014c]2243 TString s_ptrange = Form("pt_%.0f_", ptVals.at(k));
[92237ca]2244
[341014c]2245 c_pfpi_res_eta[k]->Print(pdfOutput, "pdf");
2246 c_pfpi_res_eta[k]->Print(figPath + "img_pfpi_res_" + s_ptrange + "eta.pdf", "pdf");
2247 c_pfpi_res_eta[k]->Print(figPath + "img_pfpi_res_" + s_ptrange + "eta.png", "png");
[92237ca]2248 }
2249
2250 /////////////////
2251 // PF - jets //
2252 /////////////////
2253
2254 TMultiGraph *mg_pfjet_res_e[n_etabins];
2255 TMultiGraph *mg_pfjet_res_eta[n_ptbins];
2256
2257 TLegend *leg_pfjet_res_e[n_etabins];
2258 TLegend *leg_pfjet_res_eta[n_ptbins];
2259
[341014c]2260 TGraphErrors *gr_pfjet_res_e = new TGraphErrors[n_etabins];
[dd5e213]2261 TGraphErrors *gr_pfjet_res_eta = new TGraphErrors[n_ptbins];
[92237ca]2262
[341014c]2263 TGraphErrors *gr_cajet_res_e = new TGraphErrors[n_etabins];
[dd5e213]2264 TGraphErrors *gr_cajet_res_eta = new TGraphErrors[n_ptbins];
[92237ca]2265
[341014c]2266 std::vector<resolPlot> *plots_pfjet_res_e = new std::vector<resolPlot>[n_etabins];
[dd5e213]2267 std::vector<resolPlot> *plots_pfjet_res_eta = new std::vector<resolPlot>[n_ptbins];
[341014c]2268 std::vector<resolPlot> *plots_cajet_res_e = new std::vector<resolPlot>[n_etabins];
[dd5e213]2269 std::vector<resolPlot> *plots_cajet_res_eta = new std::vector<resolPlot>[n_ptbins];
[92237ca]2270
2271 TCanvas *c_pfjet_res_e[n_etabins];
2272 TCanvas *c_pfjet_res_eta[n_ptbins];
2273
2274 // loop over eta bins
[341014c]2275 for(k = 0; k < etaVals.size() - 1; k++)
[92237ca]2276 {
2277
[341014c]2278 mg_pfjet_res_e[k] = new TMultiGraph("", "");
2279 leg_pfjet_res_e[k] = new TLegend(0.40, 0.70, 0.90, 0.90);
[92237ca]2280
[341014c]2281 HistogramsCollection(&plots_pfjet_res_e[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "pfjet");
2282 HistogramsCollection(&plots_cajet_res_e[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "cajet");
[92237ca]2283
[341014c]2284 GetJetsEres(&plots_pfjet_res_e[k], branchPFJet, branchGenJet, treeReaderJet, etaVals.at(k), etaVals.at(k + 1));
2285 GetJetsEres(&plots_cajet_res_e[k], branchCaloJet, branchGenJet, treeReaderJet, etaVals.at(k), etaVals.at(k + 1));
[92237ca]2286
[341014c]2287 gr_pfjet_res_e[k] = EresGraph(&plots_pfjet_res_e[k]);
2288 gr_cajet_res_e[k] = EresGraph(&plots_cajet_res_e[k]);
[92237ca]2289
[341014c]2290 s_etaMin = Form("%.1f", etaVals.at(k));
2291 s_etaMax = Form("%.1f", etaVals.at(k + 1));
2292 s_eta = "anti-k_{T}, R = 0.4, " + s_etaMin + " < | #eta | < " + s_etaMax;
[92237ca]2293
[341014c]2294 leg_pfjet_res_e[k]->SetTextFont(132);
2295 leg_pfjet_res_e[k]->SetHeader(s_eta);
[92237ca]2296
[341014c]2297 addResoGraph(mg_pfjet_res_e[k], &gr_cajet_res_e[k], leg_pfjet_res_e[k], markerStyles.at(0), colors.at(0), "Calorimeter Jets");
2298 addResoGraph(mg_pfjet_res_e[k], &gr_pfjet_res_e[k], leg_pfjet_res_e[k], markerStyles.at(1), colors.at(1), "Particle-flow Jets");
[92237ca]2299
[341014c]2300 c_pfjet_res_e[k] = new TCanvas("", "", 800, 600);
[92237ca]2301
[341014c]2302 mg_pfjet_res_e[k]->Draw("APE");
2303 //DrawAxis(mg_pfjet_res_e[k], leg_pfjet_res_e[k], 10, ptMax, 0.5, 100, "E [GeV]", "(resolution in E)/E (%)", true, true);
2304 DrawAxis(mg_pfjet_res_e[k], leg_pfjet_res_e[k], 10, ptMax, 0.0, 30, "E [GeV]", "(resolution in E)/E (%)", true, false);
2305 leg_pfjet_res_e[k]->Draw();
2306 pave->Draw();
[92237ca]2307
[341014c]2308 TString s_etarange = "eta_" + s_etaMin + "_" + s_etaMax + "_";
[92237ca]2309
[341014c]2310 c_pfjet_res_e[k]->Print(pdfOutput, "pdf");
2311 c_pfjet_res_e[k]->Print(figPath + "img_pfjet_res_" + s_etarange + "e.pdf", "pdf");
2312 c_pfjet_res_e[k]->Print(figPath + "img_pfjet_res_" + s_etarange + "e.png", "png");
[92237ca]2313 }
2314
[341014c]2315 // loop over eta bins
2316 for(k = 0; k < ptVals.size(); k++)
[92237ca]2317 {
2318
[341014c]2319 mg_pfjet_res_eta[k] = new TMultiGraph("", "");
2320 leg_pfjet_res_eta[k] = new TLegend(0.30, 0.70, 0.85, 0.90);
[92237ca]2321
[341014c]2322 HistogramsCollectionVsEta(&plots_pfjet_res_eta[k], etaMin, etaMax, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), "pfjet", 0.0, 2.0);
2323 HistogramsCollectionVsEta(&plots_cajet_res_eta[k], etaMin, etaMax, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), "cajet", 0.0, 2.0);
[92237ca]2324
[341014c]2325 GetJetsEresVsEta(&plots_pfjet_res_eta[k], branchPFJet, branchGenJet, treeReaderJet, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k));
2326 GetJetsEresVsEta(&plots_cajet_res_eta[k], branchCaloJet, branchGenJet, treeReaderJet, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k));
[92237ca]2327
[341014c]2328 gr_pfjet_res_eta[k] = EresGraphVsEta(&plots_pfjet_res_eta[k]);
2329 gr_cajet_res_eta[k] = EresGraphVsEta(&plots_cajet_res_eta[k]);
[92237ca]2330
[341014c]2331 s_e = Form("anti-k_{T}, R = 0.4, jets, E = %.0f GeV", ptVals.at(k));
2332 if(ptVals.at(k) >= 1000.) s_e = Form("anti-k_{T}, R = 0.4, E = %.0f TeV", ptVals.at(k) / 1000.);
[92237ca]2333
[341014c]2334 leg_pfjet_res_eta[k]->SetTextFont(132);
2335 leg_pfjet_res_eta[k]->SetHeader(s_e);
[92237ca]2336
[341014c]2337 addResoGraph(mg_pfjet_res_eta[k], &gr_cajet_res_eta[k], leg_pfjet_res_eta[k], markerStyles.at(0), colors.at(0), "Calorimeter Jets");
2338 addResoGraph(mg_pfjet_res_eta[k], &gr_pfjet_res_eta[k], leg_pfjet_res_eta[k], markerStyles.at(1), colors.at(1), "Particle-flow Jets");
[92237ca]2339
[341014c]2340 c_pfjet_res_eta[k] = new TCanvas("", "", 800, 600);
[92237ca]2341
[341014c]2342 mg_pfjet_res_eta[k]->Draw("APE");
2343 //DrawAxis(mg_pfjet_res_eta[k], leg_pfjet_res_eta[k], etaMin, etaMax, 0.1, 1000, "#eta", "(resolution in E)/E (%)", false, true);
2344 DrawAxis(mg_pfjet_res_eta[k], leg_pfjet_res_eta[k], etaMin, etaMax, 0.0, 50, "#eta", "(resolution in E)/E (%)", false, false);
2345 leg_pfjet_res_eta[k]->Draw();
2346 pave->Draw();
[92237ca]2347
[341014c]2348 TString s_ptrange = Form("pt_%.0f_", ptVals.at(k));
[92237ca]2349
[341014c]2350 c_pfjet_res_eta[k]->Print(pdfOutput, "pdf");
2351 c_pfjet_res_eta[k]->Print(figPath + "img_pfjet_res_" + s_ptrange + "eta.pdf", "pdf");
2352 c_pfjet_res_eta[k]->Print(figPath + "img_pfjet_res_" + s_ptrange + "eta.png", "png");
[92237ca]2353 }
2354
[341014c]2355 /////////////////////
2356 // PF Missing ET ///
2357 /////////////////////
[92237ca]2358
[341014c]2359 TMultiGraph *mg_met_res_ht = new TMultiGraph("", "");
2360 TLegend *leg_met_res_ht = new TLegend(0.60, 0.22, 0.90, 0.42);
[92237ca]2361
[341014c]2362 std::vector<resolPlot> plots_pfmet, plots_camet;
[92237ca]2363
[341014c]2364 HistogramsCollection(&plots_pfmet, TMath::Log10(ptMin), TMath::Log10(ptMax), "pfMET", -500, 500);
2365 HistogramsCollection(&plots_camet, TMath::Log10(ptMin), TMath::Log10(ptMax), "caMET", -500, 500);
[92237ca]2366
[341014c]2367 GetMetres(&plots_pfmet, branchGenScalarHT, branchMet, branchPFJet, treeReaderJet);
2368 GetMetres(&plots_camet, branchGenScalarHT, branchCaloMet, branchCaloJet, treeReaderJet);
[92237ca]2369
[341014c]2370 TGraphErrors gr_pfmet_res_ht = MetResGraph(&plots_pfmet, true);
2371 TGraphErrors gr_camet_res_ht = MetResGraph(&plots_camet, true);
[92237ca]2372
[341014c]2373 addResoGraph(mg_met_res_ht, &gr_camet_res_ht, leg_met_res_ht, markerStyles.at(0), colors.at(0), "Calorimeter E_{T}^{miss}");
2374 addResoGraph(mg_met_res_ht, &gr_pfmet_res_ht, leg_met_res_ht, markerStyles.at(1), colors.at(1), "Particle-flow E_{T}^{miss}");
[92237ca]2375
[341014c]2376 TCanvas *c_met_res_ht = new TCanvas("", "", 800, 600);
[92237ca]2377
[341014c]2378 mg_met_res_ht->Draw("APE");
2379 DrawAxis(mg_met_res_ht, leg_met_res_ht, 1000, 100000, 0.1, 1000, " #sum p_{T} [GeV]", "resolution in E_{x,y}^{miss} [GeV]", true, true);
[92237ca]2380
[341014c]2381 leg_met_res_ht->Draw();
2382 pave->Draw();
2383 c_met_res_ht->Print(pdfOutput, "pdf");
2384 c_met_res_ht->Print(figPath + "img_met_res_ht.pdf", "pdf");
2385 c_met_res_ht->Print(figPath + "img_met_res_ht.png", "png");
[92237ca]2386
[341014c]2387 /////////////////////////////////////////
2388 // Electron Reconstruction Efficiency ///
2389 /////////////////////////////////////////
[92237ca]2390
[341014c]2391 TMultiGraph *mg_recele_eff_pt = new TMultiGraph("", "");
2392 TMultiGraph *mg_recele_eff_eta = new TMultiGraph("", "");
[92237ca]2393
[341014c]2394 TLegend *leg_recele_eff_pt = new TLegend(0.55, 0.22, 0.90, 0.48);
2395 TLegend *leg_recele_eff_eta = new TLegend(0.55, 0.22, 0.90, 0.48);
[92237ca]2396
[341014c]2397 TGraphErrors *gr_recele_eff_pt = new TGraphErrors[n_etabins];
2398 TGraphErrors *gr_recele_eff_eta = new TGraphErrors[n_ptbins];
2399 TH1D *h_recele_eff_pt, *h_recele_eff_eta;
[92237ca]2400
[341014c]2401 // loop over eta bins
2402 for(k = 0; k < etaVals.size() - 1; k++)
2403 {
[92237ca]2404
[341014c]2405 h_recele_eff_pt = GetEffPt<Electron>(branchElectron, branchParticleElectron, "Electron", 11, ptMin, ptMax, etaVals.at(k), etaVals.at(k + 1), treeReaderElectron);
2406 gr_recele_eff_pt[k] = TGraphErrors(h_recele_eff_pt);
[92237ca]2407
[341014c]2408 s_etaMin = Form("%.1f", etaVals.at(k));
2409 s_etaMax = Form("%.1f", etaVals.at(k + 1));
[92237ca]2410
[341014c]2411 s_eta = "e^{ #pm} , " + s_etaMin + " < | #eta | < " + s_etaMax;
[92237ca]2412
[341014c]2413 gr_recele_eff_pt[k].SetName("recEff_" + s_etaMin + "_" + s_etaMax);
[92237ca]2414
[341014c]2415 addResoGraph(mg_recele_eff_pt, &gr_recele_eff_pt[k], leg_recele_eff_pt, markerStyles.at(k), colors.at(k), s_eta);
2416 }
[92237ca]2417
[341014c]2418 // loop over pt
2419 for(k = 0; k < ptVals.size(); k++)
2420 {
2421 h_recele_eff_eta = GetEffEta<Electron>(branchElectron, branchParticleElectron, "Electron", 11, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), etaMin, etaMax, treeReaderElectron);
2422 gr_recele_eff_eta[k] = TGraphErrors(h_recele_eff_eta);
[92237ca]2423
[341014c]2424 s_pt = Form("e^{ #pm} , p_{T} = %.0f GeV", ptVals.at(k));
2425 if(ptVals.at(k) >= 1000.) s_pt = Form("e^{ #pm} , p_{T} = %.0f TeV", ptVals.at(k) / 1000.);
[92237ca]2426
[341014c]2427 addResoGraph(mg_recele_eff_eta, &gr_recele_eff_eta[k], leg_recele_eff_eta, markerStyles.at(k), colors.at(k), s_pt);
2428 }
[92237ca]2429
[341014c]2430 TCanvas *c_recele_eff_pt = new TCanvas("", "", 800, 600);
[1408174]2431
[341014c]2432 mg_recele_eff_pt->Draw("APE");
2433 DrawAxis(mg_recele_eff_pt, leg_recele_eff_pt, ptMin, ptMax, 0.0, 100, "p_{T} [GeV]", "reconstruction efficiency (%)", true, false);
2434 leg_recele_eff_pt->Draw();
2435 pave->Draw();
[1408174]2436
[341014c]2437 c_recele_eff_pt->Print(pdfOutput, "pdf");
2438 c_recele_eff_pt->Print(figPath + "img_recele_eff_pt.pdf", "pdf");
2439 c_recele_eff_pt->Print(figPath + "img_recele_eff_pt.png", "png");
[1408174]2440
[341014c]2441 TCanvas *c_recele_eff_eta = new TCanvas("", "", 800, 600);
[1408174]2442
[341014c]2443 mg_recele_eff_eta->Draw("APE");
2444 DrawAxis(mg_recele_eff_eta, leg_recele_eff_eta, etaMin, etaMax, 0.0, 100, " #eta ", "reconstruction efficiency (%)", false, false);
2445 leg_recele_eff_eta->Draw();
2446 pave->Draw();
[1408174]2447
[341014c]2448 c_recele_eff_eta->Print(pdfOutput, "pdf");
2449 c_recele_eff_eta->Print(figPath + "img_recele_eff_eta.pdf", "pdf");
2450 c_recele_eff_eta->Print(figPath + "img_recele_eff_eta.png", "png");
[92237ca]2451
[341014c]2452 /////////////////////////////////////////
2453 // Muon Reconstruction Efficiency ///
2454 /////////////////////////////////////////
[92237ca]2455
[341014c]2456 TMultiGraph *mg_recmu_eff_pt = new TMultiGraph("", "");
2457 TMultiGraph *mg_recmu_eff_eta = new TMultiGraph("", "");
[92237ca]2458
[341014c]2459 TLegend *leg_recmu_eff_pt = new TLegend(0.55, 0.22, 0.90, 0.48);
2460 TLegend *leg_recmu_eff_eta = new TLegend(0.55, 0.22, 0.90, 0.48);
[1408174]2461
[341014c]2462 TGraphErrors *gr_recmu_eff_pt = new TGraphErrors[n_etabins];
2463 TGraphErrors *gr_recmu_eff_eta = new TGraphErrors[n_ptbins];
2464 TH1D *h_recmu_eff_pt, *h_recmu_eff_eta;
[1408174]2465
[341014c]2466 // loop over eta bins
2467 for(k = 0; k < etaVals.size() - 1; k++)
2468 {
[1408174]2469
[341014c]2470 h_recmu_eff_pt = GetEffPt<Muon>(branchMuon, branchParticleMuon, "muon", 13, ptMin, ptMax, etaVals.at(k), etaVals.at(k + 1), treeReaderMuon);
2471 gr_recmu_eff_pt[k] = TGraphErrors(h_recmu_eff_pt);
[1408174]2472
[341014c]2473 s_etaMin = Form("%.1f", etaVals.at(k));
2474 s_etaMax = Form("%.1f", etaVals.at(k + 1));
[1408174]2475
[341014c]2476 s_eta = "#mu^{ #pm} , " + s_etaMin + " < | #eta | < " + s_etaMax;
[1408174]2477
[341014c]2478 gr_recmu_eff_pt[k].SetName("recEff_" + s_etaMin + "_" + s_etaMax);
[1408174]2479
[341014c]2480 addResoGraph(mg_recmu_eff_pt, &gr_recmu_eff_pt[k], leg_recmu_eff_pt, markerStyles.at(k), colors.at(k), s_eta);
2481 }
[1408174]2482
[341014c]2483 // loop over pt
2484 for(k = 0; k < ptVals.size(); k++)
2485 {
2486 h_recmu_eff_eta = GetEffEta<Muon>(branchMuon, branchParticleMuon, "muon", 13, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), etaMin, etaMax, treeReaderMuon);
2487 gr_recmu_eff_eta[k] = TGraphErrors(h_recmu_eff_eta);
[1408174]2488
[341014c]2489 s_pt = Form("#mu^{ #pm} , p_{T} = %.0f GeV", ptVals.at(k));
2490 if(ptVals.at(k) >= 1000.) s_pt = Form("#mu^{ #pm} , p_{T} = %.0f TeV", ptVals.at(k) / 1000.);
[1408174]2491
[341014c]2492 addResoGraph(mg_recmu_eff_eta, &gr_recmu_eff_eta[k], leg_recmu_eff_eta, markerStyles.at(k), colors.at(k), s_pt);
2493 }
[1408174]2494
[341014c]2495 TCanvas *c_recmu_eff_pt = new TCanvas("", "", 800, 600);
[1408174]2496
[341014c]2497 mg_recmu_eff_pt->Draw("APE");
2498 DrawAxis(mg_recmu_eff_pt, leg_recmu_eff_pt, ptMin, ptMax, 0.0, 100, "p_{T} [GeV]", "reconstruction efficiency (%)", true, false);
2499 leg_recmu_eff_pt->Draw();
2500 pave->Draw();
[1408174]2501
[341014c]2502 c_recmu_eff_pt->Print(pdfOutput, "pdf");
2503 c_recmu_eff_pt->Print(figPath + "img_recmu_eff_pt.pdf", "pdf");
2504 c_recmu_eff_pt->Print(figPath + "img_recmu_eff_pt.png", "png");
[1408174]2505
[341014c]2506 TCanvas *c_recmu_eff_eta = new TCanvas("", "", 800, 600);
[1408174]2507
[341014c]2508 mg_recmu_eff_eta->Draw("APE");
2509 DrawAxis(mg_recmu_eff_eta, leg_recmu_eff_eta, etaMin, etaMax, 0.0, 100, " #eta ", "reconstruction efficiency (%)", false, false);
2510 leg_recmu_eff_eta->Draw();
2511 pave->Draw();
[1408174]2512
[341014c]2513 c_recmu_eff_eta->Print(pdfOutput, "pdf");
2514 c_recmu_eff_eta->Print(figPath + "img_recmu_eff_eta.pdf", "pdf");
2515 c_recmu_eff_eta->Print(figPath + "img_recmu_eff_eta.png", "png");
[1408174]2516
[341014c]2517 /////////////////////////////////////////
2518 // Photon Reconstruction Efficiency ///
2519 /////////////////////////////////////////
[1408174]2520
[341014c]2521 TMultiGraph *mg_recpho_eff_pt = new TMultiGraph("", "");
2522 TMultiGraph *mg_recpho_eff_eta = new TMultiGraph("", "");
[92237ca]2523
[341014c]2524 TLegend *leg_recpho_eff_pt = new TLegend(0.55, 0.22, 0.90, 0.48);
2525 TLegend *leg_recpho_eff_eta = new TLegend(0.55, 0.22, 0.90, 0.48);
[92237ca]2526
[341014c]2527 TGraphErrors *gr_recpho_eff_pt = new TGraphErrors[n_etabins];
2528 TGraphErrors *gr_recpho_eff_eta = new TGraphErrors[n_ptbins];
2529 TH1D *h_recpho_eff_pt, *h_recpho_eff_eta;
[1408174]2530
[341014c]2531 // loop over eta bins
2532 for(k = 0; k < etaVals.size() - 1; k++)
2533 {
[1408174]2534
[341014c]2535 h_recpho_eff_pt = GetEffPt<Photon>(branchPhoton, branchParticlePhoton, "Photon", 22, ptMin, ptMax, etaVals.at(k), etaVals.at(k + 1), treeReaderPhoton);
2536 gr_recpho_eff_pt[k] = TGraphErrors(h_recpho_eff_pt);
[1408174]2537
[341014c]2538 s_etaMin = Form("%.1f", etaVals.at(k));
2539 s_etaMax = Form("%.1f", etaVals.at(k + 1));
[1408174]2540
[341014c]2541 s_eta = "#gamma , " + s_etaMin + " < | #eta | < " + s_etaMax;
[92237ca]2542
[341014c]2543 gr_recpho_eff_pt[k].SetName("recEff_" + s_etaMin + "_" + s_etaMax);
[1408174]2544
[341014c]2545 addResoGraph(mg_recpho_eff_pt, &gr_recpho_eff_pt[k], leg_recpho_eff_pt, markerStyles.at(k), colors.at(k), s_eta);
2546 }
[1160e4f]2547
[341014c]2548 // loop over pt
2549 for(k = 0; k < ptVals.size(); k++)
2550 {
2551 h_recpho_eff_eta = GetEffEta<Photon>(branchPhoton, branchParticlePhoton, "Photon", 22, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), etaMin, etaMax, treeReaderPhoton);
2552 gr_recpho_eff_eta[k] = TGraphErrors(h_recpho_eff_eta);
[1160e4f]2553
[341014c]2554 s_pt = Form("#gamma , p_{T} = %.0f GeV", ptVals.at(k));
2555 if(ptVals.at(k) >= 1000.) s_pt = Form("#gamma , p_{T} = %.0f TeV", ptVals.at(k) / 1000.);
[1160e4f]2556
[341014c]2557 addResoGraph(mg_recpho_eff_eta, &gr_recpho_eff_eta[k], leg_recpho_eff_eta, markerStyles.at(k), colors.at(k), s_pt);
2558 }
[1408174]2559
[341014c]2560 TCanvas *c_recpho_eff_pt = new TCanvas("", "", 800, 600);
[1408174]2561
[341014c]2562 mg_recpho_eff_pt->Draw("APE");
2563 DrawAxis(mg_recpho_eff_pt, leg_recpho_eff_pt, ptMin, ptMax, 0.0, 100, "p_{T} [GeV]", "reconstruction efficiency (%)", true, false);
2564 leg_recpho_eff_pt->Draw();
2565 pave->Draw();
[92237ca]2566
[341014c]2567 c_recpho_eff_pt->Print(pdfOutput, "pdf");
2568 c_recpho_eff_pt->Print(figPath + "img_recpho_eff_pt.pdf", "pdf");
2569 c_recpho_eff_pt->Print(figPath + "img_recpho_eff_pt.png", "png");
[92237ca]2570
[341014c]2571 TCanvas *c_recpho_eff_eta = new TCanvas("", "", 800, 600);
[92237ca]2572
[341014c]2573 mg_recpho_eff_eta->Draw("APE");
2574 DrawAxis(mg_recpho_eff_eta, leg_recpho_eff_eta, etaMin, etaMax, 0.0, 100, " #eta ", "reconstruction efficiency (%)", false, false);
2575 leg_recpho_eff_eta->Draw();
2576 pave->Draw();
[92237ca]2577
[341014c]2578 c_recpho_eff_eta->Print(pdfOutput, "pdf");
2579 c_recpho_eff_eta->Print(figPath + "img_recpho_eff_eta.pdf", "pdf");
2580 c_recpho_eff_eta->Print(figPath + "img_recpho_eff_eta.png", "png");
[92237ca]2581
[341014c]2582 /////////////////////////////////////////
2583 // B-jets Efficiency/ mistag rates ///
2584 /////////////////////////////////////////
[92237ca]2585
[341014c]2586 TMultiGraph *mg_recbjet_eff_pt = new TMultiGraph("", "");
2587 TMultiGraph *mg_recbjet_eff_eta = new TMultiGraph("", "");
[92237ca]2588
[341014c]2589 TLegend *leg_recbjet_eff_pt = new TLegend(0.50, 0.22, 0.90, 0.48);
2590 TLegend *leg_recbjet_eff_eta = new TLegend(0.50, 0.22, 0.90, 0.48);
[92237ca]2591
[341014c]2592 TGraphErrors *gr_recbjet_eff_pt = new TGraphErrors[n_etabins];
2593 TGraphErrors *gr_recbjet_eff_eta = new TGraphErrors[n_ptbins];
2594 TH1D *h_recbjet_eff_pt, *h_recbjet_eff_eta;
[92237ca]2595
[341014c]2596 // loop over eta bins
2597 for(k = 0; k < etaVals.size() - 1; k++)
2598 {
[92237ca]2599
[341014c]2600 h_recbjet_eff_pt = GetJetEffPt<Jet>(branchPFBJet, "BJet", 5, ptMin, ptMax, etaVals.at(k), etaVals.at(k + 1), treeReaderBJet);
2601 //h_recbjet_eff_pt = GetEffPt<Jet>(branchPFBJet, branchParticleBJet, "BJet", 5, ptMin, ptMax, etaVals.at(k), etaVals.at(k+1), treeReaderBJet);
2602 gr_recbjet_eff_pt[k] = TGraphErrors(h_recbjet_eff_pt);
[92237ca]2603
[341014c]2604 s_etaMin = Form("%.1f", etaVals.at(k));
2605 s_etaMax = Form("%.1f", etaVals.at(k + 1));
[1408174]2606
[341014c]2607 s_eta = "b-jet , " + s_etaMin + " < | #eta | < " + s_etaMax;
[92237ca]2608
[341014c]2609 gr_recbjet_eff_pt[k].SetName("recEff_" + s_etaMin + "_" + s_etaMax);
[92237ca]2610
[341014c]2611 addResoGraph(mg_recbjet_eff_pt, &gr_recbjet_eff_pt[k], leg_recbjet_eff_pt, markerStyles.at(k), colors.at(k), s_eta);
2612 }
[1408174]2613
[341014c]2614 // loop over pt
2615 for(k = 0; k < ptVals.size(); k++)
2616 {
2617 h_recbjet_eff_eta = GetJetEffEta<Jet>(branchPFBJet, "BJet", 5, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), etaMin, etaMax, treeReaderBJet);
2618 //h_recbjet_eff_eta = GetEffEta<Jet>(branchPFBJet, branchParticleBJet, "BJet", 5, 0.5*ptVals.at(k), 2.0*ptVals.at(k) ,etaMin, etaMax , treeReaderBJet);
2619 gr_recbjet_eff_eta[k] = TGraphErrors(h_recbjet_eff_eta);
[734b267]2620
[341014c]2621 s_pt = Form("b-jet , p_{T} = %.0f GeV", ptVals.at(k));
2622 if(ptVals.at(k) >= 1000.) s_pt = Form("b-jet , p_{T} = %.0f TeV", ptVals.at(k) / 1000.);
[1408174]2623
[341014c]2624 addResoGraph(mg_recbjet_eff_eta, &gr_recbjet_eff_eta[k], leg_recbjet_eff_eta, markerStyles.at(k), colors.at(k), s_pt);
2625 }
[1408174]2626
[341014c]2627 TCanvas *c_recbjet_eff_pt = new TCanvas("", "", 800, 600);
[1408174]2628
[341014c]2629 mg_recbjet_eff_pt->Draw("APE");
2630 DrawAxis(mg_recbjet_eff_pt, leg_recbjet_eff_pt, ptMin, ptMax, 0.0, 100, "p_{T} [GeV]", "b - tag efficiency (%)", true, false);
2631 leg_recbjet_eff_pt->Draw();
2632 pave->Draw();
[1408174]2633
[341014c]2634 c_recbjet_eff_pt->Print(pdfOutput, "pdf");
2635 c_recbjet_eff_pt->Print(figPath + "img_recbjet_eff_pt.pdf", "pdf");
2636 c_recbjet_eff_pt->Print(figPath + "img_recbjet_eff_pt.png", "png");
[1408174]2637
[341014c]2638 TCanvas *c_recbjet_eff_eta = new TCanvas("", "", 800, 600);
[1408174]2639
[341014c]2640 mg_recbjet_eff_eta->Draw("APE");
2641 DrawAxis(mg_recbjet_eff_eta, leg_recbjet_eff_eta, etaMin, etaMax, 0.0, 100, " #eta ", "b - tag efficiency (%)", false, false);
2642 leg_recbjet_eff_eta->Draw();
2643 pave->Draw();
[1408174]2644
[341014c]2645 c_recbjet_eff_eta->Print(pdfOutput, "pdf");
2646 c_recbjet_eff_eta->Print(figPath + "img_recbjet_eff_eta.pdf", "pdf");
2647 c_recbjet_eff_eta->Print(figPath + "img_recbjet_eff_eta.png", "png");
[1408174]2648
[341014c]2649 // ------ c - mistag ------
[1408174]2650
[341014c]2651 TMultiGraph *mg_recbjet_cmis_pt = new TMultiGraph("", "");
2652 TMultiGraph *mg_recbjet_cmis_eta = new TMultiGraph("", "");
[1408174]2653
[341014c]2654 TLegend *leg_recbjet_cmis_pt = new TLegend(0.50, 0.64, 0.90, 0.90);
2655 TLegend *leg_recbjet_cmis_eta = new TLegend(0.50, 0.64, 0.90, 0.90);
[1408174]2656
[341014c]2657 TGraphErrors *gr_recbjet_cmis_pt = new TGraphErrors[n_etabins];
2658 TGraphErrors *gr_recbjet_cmis_eta = new TGraphErrors[n_ptbins];
2659 TH1D *h_recbjet_cmis_pt, *h_recbjet_cmis_eta;
[1408174]2660
[341014c]2661 // loop over eta bins
2662 for(k = 0; k < etaVals.size() - 1; k++)
2663 {
[1408174]2664
[341014c]2665 h_recbjet_cmis_pt = GetJetEffPt<Jet>(branchPFCJet, "CJet", 4, ptMin, ptMax, etaVals.at(k), etaVals.at(k + 1), treeReaderCJet);
2666 //h_recbjet_cmis_pt = GetEffPt<Jet>(branchPFCJet, branchParticleCJet, "CJet", 4, ptMin, ptMax, etaVals.at(k), etaVals.at(k+1), treeReaderCJet);
2667 gr_recbjet_cmis_pt[k] = TGraphErrors(h_recbjet_cmis_pt);
[1408174]2668
[341014c]2669 s_etaMin = Form("%.1f", etaVals.at(k));
2670 s_etaMax = Form("%.1f", etaVals.at(k + 1));
[1408174]2671
[341014c]2672 s_eta = "c-jet , " + s_etaMin + " < | #eta | < " + s_etaMax;
[1408174]2673
[341014c]2674 gr_recbjet_cmis_pt[k].SetName("recEff_" + s_etaMin + "_" + s_etaMax);
[1408174]2675
[341014c]2676 addResoGraph(mg_recbjet_cmis_pt, &gr_recbjet_cmis_pt[k], leg_recbjet_cmis_pt, markerStyles.at(k), colors.at(k), s_eta);
2677 }
[1408174]2678
[341014c]2679 // loop over pt
2680 for(k = 0; k < ptVals.size(); k++)
2681 {
2682 h_recbjet_cmis_eta = GetJetEffEta<Jet>(branchPFCJet, "CJet", 4, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), etaMin, etaMax, treeReaderCJet);
2683 //h_recbjet_cmis_eta = GetEffEta<Jet>(branchPFCJet, branchParticleCJet, "CJet", 4, 0.5*ptVals.at(k), 2.0*ptVals.at(k) ,etaMin, etaMax , treeReaderCJet);
2684 gr_recbjet_cmis_eta[k] = TGraphErrors(h_recbjet_cmis_eta);
[1408174]2685
[341014c]2686 s_pt = Form("c-jet , p_{T} = %.0f GeV", ptVals.at(k));
2687 if(ptVals.at(k) >= 1000.) s_pt = Form("c-jet , p_{T} = %.0f TeV", ptVals.at(k) / 1000.);
[1408174]2688
[341014c]2689 addResoGraph(mg_recbjet_cmis_eta, &gr_recbjet_cmis_eta[k], leg_recbjet_cmis_eta, markerStyles.at(k), colors.at(k), s_pt);
2690 }
[1408174]2691
[341014c]2692 TCanvas *c_recbjet_cmis_pt = new TCanvas("", "", 800, 600);
[1408174]2693
[341014c]2694 mg_recbjet_cmis_pt->Draw("APE");
2695 DrawAxis(mg_recbjet_cmis_pt, leg_recbjet_cmis_pt, ptMin, ptMax, 0.0, 20, "p_{T} [GeV]", "c - mistag rate (%)", true, false);
2696 leg_recbjet_cmis_pt->Draw();
2697 pave->Draw();
[1408174]2698
[341014c]2699 c_recbjet_cmis_pt->Print(pdfOutput, "pdf");
2700 c_recbjet_cmis_pt->Print(figPath + "img_recbjet_cmis_pt.pdf", "pdf");
2701 c_recbjet_cmis_pt->Print(figPath + "img_recbjet_cmis_pt.png", "png");
[1408174]2702
[341014c]2703 TCanvas *c_recbjet_cmis_eta = new TCanvas("", "", 800, 600);
[1408174]2704
[341014c]2705 mg_recbjet_cmis_eta->Draw("APE");
2706 DrawAxis(mg_recbjet_cmis_eta, leg_recbjet_cmis_eta, etaMin, etaMax, 0.0, 20, " #eta ", "c - mistag rate (%)", false, false);
2707 leg_recbjet_cmis_eta->Draw();
2708 pave->Draw();
[1408174]2709
[341014c]2710 c_recbjet_cmis_eta->Print(pdfOutput, "pdf");
2711 c_recbjet_cmis_eta->Print(figPath + "img_recbjet_cmis_eta.pdf", "pdf");
2712 c_recbjet_cmis_eta->Print(figPath + "img_recbjet_cmis_eta.png", "png");
[1408174]2713
[341014c]2714 // ------ light - mistag ------
[1408174]2715
[341014c]2716 TMultiGraph *mg_recbjet_lmis_pt = new TMultiGraph("", "");
2717 TMultiGraph *mg_recbjet_lmis_eta = new TMultiGraph("", "");
[cb159ed]2718
[341014c]2719 TLegend *leg_recbjet_lmis_pt = new TLegend(0.50, 0.64, 0.90, 0.90);
2720 TLegend *leg_recbjet_lmis_eta = new TLegend(0.50, 0.64, 0.90, 0.90);
[1408174]2721
[341014c]2722 TGraphErrors *gr_recbjet_lmis_pt = new TGraphErrors[n_etabins];
2723 TGraphErrors *gr_recbjet_lmis_eta = new TGraphErrors[n_ptbins];
2724 TH1D *h_recbjet_lmis_pt, *h_recbjet_lmis_eta;
[1408174]2725
[341014c]2726 // loop over eta bins
2727 for(k = 0; k < etaVals.size() - 1; k++)
2728 {
[1408174]2729
[341014c]2730 h_recbjet_lmis_pt = GetJetEffPt<Jet>(branchJet, "Jet", 1, ptMin, ptMax, etaVals.at(k), etaVals.at(k + 1), treeReaderJet);
2731 //h_recbjet_lmis_pt = GetEffPt<Jet>(branchJet, branchParticleJet, "Jet", 1, ptMin, ptMax, etaVals.at(k), etaVals.at(k+1), treeReaderJet);
2732 gr_recbjet_lmis_pt[k] = TGraphErrors(h_recbjet_lmis_pt);
[1408174]2733
[341014c]2734 s_etaMin = Form("%.1f", etaVals.at(k));
2735 s_etaMax = Form("%.1f", etaVals.at(k + 1));
[1408174]2736
[341014c]2737 s_eta = "uds-jet , " + s_etaMin + " < | #eta | < " + s_etaMax;
[92237ca]2738
[341014c]2739 gr_recbjet_lmis_pt[k].SetName("recEff_" + s_etaMin + "_" + s_etaMax);
[92237ca]2740
[341014c]2741 addResoGraph(mg_recbjet_lmis_pt, &gr_recbjet_lmis_pt[k], leg_recbjet_lmis_pt, markerStyles.at(k), colors.at(k), s_eta);
2742 }
[92237ca]2743
[341014c]2744 // loop over pt
2745 for(k = 0; k < ptVals.size(); k++)
2746 {
2747 h_recbjet_lmis_eta = GetJetEffEta<Jet>(branchJet, "Jet", 1, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), etaMin, etaMax, treeReaderJet);
2748 //h_recbjet_lmis_eta = GetEffEta<Jet>(branchJet, branchParticleJet, "Jet", 1, 0.5*ptVals.at(k), 2.0*ptVals.at(k) ,etaMin, etaMax , treeReaderJet);
2749 gr_recbjet_lmis_eta[k] = TGraphErrors(h_recbjet_lmis_eta);
[92237ca]2750
[341014c]2751 s_pt = Form("uds-jet , p_{T} = %.0f GeV", ptVals.at(k));
2752 if(ptVals.at(k) >= 1000.) s_pt = Form("uds-jet , p_{T} = %.0f TeV", ptVals.at(k) / 1000.);
[92237ca]2753
[341014c]2754 addResoGraph(mg_recbjet_lmis_eta, &gr_recbjet_lmis_eta[k], leg_recbjet_lmis_eta, markerStyles.at(k), colors.at(k), s_pt);
2755 }
[92237ca]2756
[341014c]2757 TCanvas *c_recbjet_lmis_pt = new TCanvas("", "", 800, 600);
[92237ca]2758
[341014c]2759 mg_recbjet_lmis_pt->Draw("APE");
[92237ca]2760
[341014c]2761 DrawAxis(mg_recbjet_lmis_pt, leg_recbjet_lmis_pt, ptMin, ptMax, 0.0, 1.0, "p_{T} [GeV]", "light - mistag rate (%)", true, false);
[92237ca]2762
[341014c]2763 leg_recbjet_lmis_pt->Draw();
2764 pave->Draw();
[92237ca]2765
[341014c]2766 c_recbjet_lmis_pt->Print(pdfOutput, "pdf");
2767 c_recbjet_lmis_pt->Print(figPath + "img_recbjet_lmis_pt.pdf", "pdf");
2768 c_recbjet_lmis_pt->Print(figPath + "img_recbjet_lmis_pt.png", "png");
[92237ca]2769
[341014c]2770 TCanvas *c_recbjet_lmis_eta = new TCanvas("", "", 800, 600);
[fa068d3]2771
[341014c]2772 mg_recbjet_lmis_eta->Draw("APE");
2773 DrawAxis(mg_recbjet_lmis_eta, leg_recbjet_lmis_eta, etaMin, etaMax, 0.0, 1.0, " #eta ", "light - mistag rate (%)", false, false);
2774 leg_recbjet_lmis_eta->Draw();
2775 pave->Draw();
[1408174]2776
[341014c]2777 c_recbjet_lmis_eta->Print(pdfOutput, "pdf");
2778 c_recbjet_lmis_eta->Print(figPath + "img_recbjet_lmis_eta.pdf", "pdf");
2779 c_recbjet_lmis_eta->Print(figPath + "img_recbjet_lmis_eta.png", "png");
[1408174]2780
[341014c]2781 ///////////////////////////////////////////
2782 // tau-jets Efficiency/ mistag rates ///
2783 ///////////////////////////////////////////
[92237ca]2784
[341014c]2785 TMultiGraph *mg_rectaujet_eff_pt = new TMultiGraph("", "");
2786 TMultiGraph *mg_rectaujet_eff_eta = new TMultiGraph("", "");
[92237ca]2787
[341014c]2788 TLegend *leg_rectaujet_eff_pt = new TLegend(0.50, 0.22, 0.90, 0.48);
2789 TLegend *leg_rectaujet_eff_eta = new TLegend(0.50, 0.22, 0.90, 0.48);
[92237ca]2790
[341014c]2791 TGraphErrors *gr_rectaujet_eff_pt = new TGraphErrors[n_etabins];
2792 TGraphErrors *gr_rectaujet_eff_eta = new TGraphErrors[n_ptbins];
2793 TH1D *h_rectaujet_eff_pt, *h_rectaujet_eff_eta;
[92237ca]2794
[341014c]2795 // loop over eta bins
2796 for(k = 0; k < etaVals.size() - 1; k++)
2797 {
[92237ca]2798
[341014c]2799 h_rectaujet_eff_pt = GetTauEffPt<Jet>(branchPFTauJet, branchParticleTauJet, "TauJet", 15, ptMin, ptMax, etaVals.at(k), etaVals.at(k + 1), treeReaderTauJet);
2800 gr_rectaujet_eff_pt[k] = TGraphErrors(h_rectaujet_eff_pt);
[92237ca]2801
[341014c]2802 s_etaMin = Form("%.1f", etaVals.at(k));
2803 s_etaMax = Form("%.1f", etaVals.at(k + 1));
[92237ca]2804
[341014c]2805 s_eta = "#tau-jet , " + s_etaMin + " < | #eta | < " + s_etaMax;
[92237ca]2806
[341014c]2807 gr_rectaujet_eff_pt[k].SetName("recEff_" + s_etaMin + "_" + s_etaMax);
[92237ca]2808
[341014c]2809 addResoGraph(mg_rectaujet_eff_pt, &gr_rectaujet_eff_pt[k], leg_rectaujet_eff_pt, markerStyles.at(k), colors.at(k), s_eta);
2810 }
[92237ca]2811
[341014c]2812 // loop over pt
2813 for(k = 0; k < ptVals.size(); k++)
2814 {
2815 h_rectaujet_eff_eta = GetTauEffEta<Jet>(branchPFTauJet, branchParticleTauJet, "TauJet", 15, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), etaMin, etaMax, treeReaderTauJet);
2816 gr_rectaujet_eff_eta[k] = TGraphErrors(h_rectaujet_eff_eta);
[92237ca]2817
[341014c]2818 s_pt = Form("#tau-jet , p_{T} = %.0f GeV", ptVals.at(k));
2819 if(ptVals.at(k) >= 1000.) s_pt = Form("#tau-jet , p_{T} = %.0f TeV", ptVals.at(k) / 1000.);
[92237ca]2820
[341014c]2821 addResoGraph(mg_rectaujet_eff_eta, &gr_rectaujet_eff_eta[k], leg_rectaujet_eff_eta, markerStyles.at(k), colors.at(k), s_pt);
2822 }
[1408174]2823
[341014c]2824 TCanvas *c_rectaujet_eff_pt = new TCanvas("", "", 800, 600);
[1408174]2825
[341014c]2826 mg_rectaujet_eff_pt->Draw("APE");
2827 DrawAxis(mg_rectaujet_eff_pt, leg_rectaujet_eff_pt, ptMin, ptMax, 0.0, 100, "p_{T} [GeV]", "#tau - tag efficiency (%)", true, false);
2828 leg_rectaujet_eff_pt->Draw();
2829 pave->Draw();
[1408174]2830
[341014c]2831 c_rectaujet_eff_pt->Print(pdfOutput, "pdf");
2832 c_rectaujet_eff_pt->Print(figPath + "img_rectaujet_eff_pt.pdf", "pdf");
2833 c_rectaujet_eff_pt->Print(figPath + "img_rectaujet_eff_pt.png", "png");
[1408174]2834
[341014c]2835 TCanvas *c_rectaujet_eff_eta = new TCanvas("", "", 800, 600);
[1408174]2836
[341014c]2837 mg_rectaujet_eff_eta->Draw("APE");
2838 DrawAxis(mg_rectaujet_eff_eta, leg_rectaujet_eff_eta, etaMin, etaMax, 0.0, 100., " #eta ", "#tau - tag efficiency (%)", false, false);
2839 leg_rectaujet_eff_eta->Draw();
2840 pave->Draw();
[1408174]2841
[341014c]2842 c_rectaujet_eff_eta->Print(pdfOutput, "pdf");
2843 c_rectaujet_eff_eta->Print(figPath + "img_rectaujet_eff_eta.pdf", "pdf");
2844 c_rectaujet_eff_eta->Print(figPath + "img_rectaujet_eff_eta.png", "png");
[1408174]2845
[341014c]2846 //--------------- tau mistag rate ----------
[1408174]2847
[341014c]2848 TMultiGraph *mg_rectaujet_mis_pt = new TMultiGraph("", "");
2849 TMultiGraph *mg_rectaujet_mis_eta = new TMultiGraph("", "");
[1408174]2850
[341014c]2851 TLegend *leg_rectaujet_mis_pt = new TLegend(0.50, 0.64, 0.90, 0.90);
2852 TLegend *leg_rectaujet_mis_eta = new TLegend(0.50, 0.64, 0.90, 0.90);
[1408174]2853
[341014c]2854 TGraphErrors *gr_rectaujet_mis_pt = new TGraphErrors[n_etabins];
2855 TGraphErrors *gr_rectaujet_mis_eta = new TGraphErrors[n_ptbins];
2856 TH1D *h_rectaujet_mis_pt, *h_rectaujet_mis_eta;
[1408174]2857
[341014c]2858 // loop over eta bins
2859 for(k = 0; k < etaVals.size() - 1; k++)
2860 {
[1408174]2861
[341014c]2862 h_rectaujet_mis_pt = GetTauEffPt<Jet>(branchJet, branchParticleJet, "TauJet", 1, ptMin, ptMax, etaVals.at(k), etaVals.at(k + 1), treeReaderJet);
2863 gr_rectaujet_mis_pt[k] = TGraphErrors(h_rectaujet_mis_pt);
[1408174]2864
[341014c]2865 s_etaMin = Form("%.1f", etaVals.at(k));
2866 s_etaMax = Form("%.1f", etaVals.at(k + 1));
[1408174]2867
[341014c]2868 s_eta = "uds-jet , " + s_etaMin + " < | #eta | < " + s_etaMax;
[1408174]2869
[341014c]2870 gr_rectaujet_mis_pt[k].SetName("recEff_" + s_etaMin + "_" + s_etaMax);
[1408174]2871
[341014c]2872 addResoGraph(mg_rectaujet_mis_pt, &gr_rectaujet_mis_pt[k], leg_rectaujet_mis_pt, markerStyles.at(k), colors.at(k), s_eta);
2873 }
[1408174]2874
[341014c]2875 // loop over pt
2876 for(k = 0; k < ptVals.size(); k++)
2877 {
2878 h_rectaujet_mis_eta = GetTauEffEta<Jet>(branchJet, branchParticleJet, "TauJet", 1, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), etaMin, etaMax, treeReaderJet);
2879 gr_rectaujet_mis_eta[k] = TGraphErrors(h_rectaujet_mis_eta);
[1408174]2880
[341014c]2881 s_pt = Form("uds-jet , p_{T} = %.0f GeV", ptVals.at(k));
2882 if(ptVals.at(k) >= 1000.) s_pt = Form("uds-jet , p_{T} = %.0f TeV", ptVals.at(k) / 1000.);
[1408174]2883
[341014c]2884 addResoGraph(mg_rectaujet_mis_eta, &gr_rectaujet_mis_eta[k], leg_rectaujet_mis_eta, markerStyles.at(k), colors.at(k), s_pt);
2885 }
[1408174]2886
[341014c]2887 TCanvas *c_rectaujet_mis_pt = new TCanvas("", "", 800, 600);
[1408174]2888
[341014c]2889 mg_rectaujet_mis_pt->Draw("APE");
2890 DrawAxis(mg_rectaujet_mis_pt, leg_rectaujet_mis_pt, ptMin, ptMax, 0.0, 5., "p_{T} [GeV]", "#tau - mistag(%)", true, false);
2891 leg_rectaujet_mis_pt->Draw();
2892 pave->Draw();
[1408174]2893
[341014c]2894 c_rectaujet_mis_pt->Print(pdfOutput, "pdf");
2895 c_rectaujet_mis_pt->Print(figPath + "img_rectaujet_mis_pt.pdf", "pdf");
2896 c_rectaujet_mis_pt->Print(figPath + "img_rectaujet_mis_pt.png", "png");
[1408174]2897
[341014c]2898 TCanvas *c_rectaujet_mis_eta = new TCanvas("", "", 800, 600);
[1408174]2899
[341014c]2900 mg_rectaujet_mis_eta->Draw("APE");
2901 DrawAxis(mg_rectaujet_mis_eta, leg_rectaujet_mis_eta, etaMin, etaMax, 0.0, 5., " #eta ", "#tau - mistag (%)", false, false);
2902 leg_rectaujet_mis_eta->Draw();
2903 pave->Draw();
[1408174]2904
[341014c]2905 c_rectaujet_mis_eta->Print(pdfOutput + ")", "pdf");
2906 c_rectaujet_mis_eta->Print(figPath + "img_rectaujet_mis_eta.pdf", "pdf");
2907 c_rectaujet_mis_eta->Print(figPath + "img_rectaujet_mis_eta.png", "png");
[1408174]2908
[341014c]2909 // ---- store resolution histograms in the output (for leave efficiencies out) ---
[1408174]2910
[341014c]2911 TFile *fout = new TFile(outputFile, "recreate");
[1408174]2912
[341014c]2913 for(int bin = 0; bin < Nbins; bin++)
[1408174]2914 {
2915
[341014c]2916 for(k = 0; k < etaVals.size() - 1; k++)
[92237ca]2917 {
[341014c]2918 plots_trkpi_res_pt[k].at(bin).resolHist->Write();
2919 plots_trkele_res_pt[k].at(bin).resolHist->Write();
2920 plots_trkmu_res_pt[k].at(bin).resolHist->Write();
2921 plots_ecal_res_e[k].at(bin).resolHist->Write();
2922 plots_hcal_res_e[k].at(bin).resolHist->Write();
2923 plots_pfele_res_e[k].at(bin).resolHist->Write();
2924 plots_pfpi_res_e[k].at(bin).resolHist->Write();
2925 plots_pfjet_res_e[k].at(bin).resolHist->Write();
2926 plots_cajet_res_e[k].at(bin).resolHist->Write();
[92237ca]2927 }
2928
[341014c]2929 for(k = 0; k < ptVals.size(); k++)
[92237ca]2930 {
[341014c]2931 plots_trkpi_res_eta[k].at(bin).resolHist->Write();
2932 plots_trkele_res_eta[k].at(bin).resolHist->Write();
2933 plots_trkmu_res_eta[k].at(bin).resolHist->Write();
2934 plots_ecal_res_eta[k].at(bin).resolHist->Write();
2935 plots_hcal_res_eta[k].at(bin).resolHist->Write();
2936 plots_pfele_res_eta[k].at(bin).resolHist->Write();
2937 plots_pfpi_res_eta[k].at(bin).resolHist->Write();
2938 plots_pfjet_res_eta[k].at(bin).resolHist->Write();
2939 plots_cajet_res_eta[k].at(bin).resolHist->Write();
[92237ca]2940 }
[1408174]2941
[92237ca]2942 plots_pfmet.at(bin).resolHist->Write();
2943 plots_camet.at(bin).resolHist->Write();
[1408174]2944 }
2945
2946 fout->Write();
2947
2948 cout << "** Exiting..." << endl;
2949
2950 delete treeReaderElectron;
2951 delete treeReaderMuon;
2952 delete treeReaderPhoton;
2953 delete treeReaderJet;
2954 delete treeReaderBJet;
2955 delete treeReaderTauJet;
2956 delete chainElectron;
2957 delete chainMuon;
2958 delete chainPhoton;
2959 delete chainJet;
2960 delete chainBJet;
2961 delete chainTauJet;
2962}
[734b267]2963
2964//------------------------------------------------------------------------------
2965
2966int main(int argc, char *argv[])
2967{
[181c061]2968 char *appName = "DelphesValidation";
[734b267]2969
[85ad2b9]2970 if(argc != 12)
[734b267]2971 {
[181c061]2972 cout << " Usage: " << appName << " input_file_electron input_file_muon input_file_photon input_file_jet input_file_bjet input_file_taujet output_file version" << endl;
2973 cout << " input_file_pion - input file in ROOT format ('Delphes' tree)," << endl;
2974 cout << " input_file_electron - input file in ROOT format ('Delphes' tree)," << endl;
[1408174]2975 cout << " input_file_muon - input file in ROOT format ('Delphes' tree)," << endl;
2976 cout << " input_file_photon - input file in ROOT format ('Delphes' tree)," << endl;
[92237ca]2977 cout << " input_file_neutralhadron - input file in ROOT format ('Delphes' tree)," << endl;
[1408174]2978 cout << " input_file_jet - input file in ROOT format ('Delphes' tree)," << endl;
2979 cout << " input_file_bjet - input file in ROOT format ('Delphes' tree)," << endl;
[92237ca]2980 cout << " input_file_cjet - input file in ROOT format ('Delphes' tree)," << endl;
[1408174]2981 cout << " input_file_taujet - input file in ROOT format ('Delphes' tree)," << endl;
[181c061]2982 cout << " output_file - output file in ROOT format," << endl;
2983 cout << " version - Delphes version" << endl;
[85ad2b9]2984
[734b267]2985 return 1;
2986 }
2987
2988 gROOT->SetBatch();
2989
2990 int appargc = 1;
2991 char *appargv[] = {appName};
2992 TApplication app(appName, &appargc, appargv);
2993
[181c061]2994 DelphesValidation(argv[1], argv[2], argv[3], argv[4], argv[5], argv[6], argv[7], argv[8], argv[9], argv[10], argv[11]);
[734b267]2995}
Note: See TracBrowser for help on using the repository browser.