Fork me on GitHub

Changeset 92237ca in git for examples


Ignore:
Timestamp:
Sep 23, 2016, 1:33:07 PM (8 years ago)
Author:
Michele Selvaggi <michele.selvaggi@…>
Branches:
ImprovedOutputFile, Timing, dual_readout, llp, master
Children:
6387aea
Parents:
87ff5a0
Message:

include more plots, change style, etc...

File:
1 edited

Legend:

Unmodified
Added
Removed

  • examples/Validation.cpp

    r87ff5a0 r92237ca  
    5555//------------------------------------------------------------------------------
    5656
    57 double ptrangemin = 10;
    58 double ptrangemax = 10000;
    59 static const int Nbins = 20;
     57static const int Nbins = 50;
    6058
    6159int objStyle = 1;
    … …  
    8280double topLeftLegYmax = 0.85;
    8381
     82unsigned int k;
     83
    8484
    8585struct resolPlot
    8686{
    87   TH1 *cenResolHist;
    88   TH1 *fwdResolHist;
     87  TH1* resolHist;
    8988  double ptmin;
    9089  double ptmax;
     90  double etamin;
     91  double etamax;
    9192  double xmin;
    9293  double xmax;
    … …  
    9596  resolPlot();
    9697  resolPlot(double ptdown, double ptup, TString object);
     98  resolPlot(double etadown, double etaup, double ptdown, double ptup, TString object);
    9799  void set(double ptdown, double ptup, TString object, double xmin = 0, double xmax = 2);
     100  void set(double etadown, double etaup, double ptdown, double ptup, TString object, double xmin = 0, double xmax = 2);
    98101  void print(){std::cout << ptmin << std::endl;}
    99102};
    … …  
    107110{
    108111  this->set(ptdown,ptup,object);
     112}
     113
     114resolPlot::resolPlot(double etadown, double etaup, double ptdown, double ptup, TString object)
     115{
     116  this->set(etadown, etaup, ptdown, ptup, object);
    109117}
    110118
    … …  
    115123  obj = object;
    116124
    117   cenResolHist = new TH1D(obj+"_delta_pt_"+Form("%4.2f",ptmin)+"_"+Form("%4.2f",ptmax)+"_cen", obj+"_delta_pt_"+Form("%4.2f",ptmin)+"_"+Form("%4.2f",ptmax)+"_cen", 200,  xmin, xmax);
    118   fwdResolHist = new TH1D(obj+"_delta_pt_"+Form("%4.2f",ptmin)+"_"+Form("%4.2f",ptmax)+"_fwd", obj+"_delta_pt_"+Form("%4.2f",ptmin)+"_"+Form("%4.2f",ptmax)+"_fwd", 200,  xmin, xmax);
     125  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);
    119126}
     127
     128void resolPlot::set(double etadown, double etaup, double ptdown, double ptup, TString object, double xmin, double xmax)
     129{
     130  etamin = etadown;
     131  etamax = etaup;
     132  ptmin = ptdown;
     133  ptmax = ptup;
     134  obj = object;
     135
     136  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);
     137}
     138
     139
    120140
    121141void HistogramsCollection(std::vector<resolPlot> *histos, double ptmin, double ptmax, TString obj, double xmin = 0, double xmax = 2)
    … …  
    136156}
    137157
     158
     159void HistogramsCollectionVsEta(std::vector<resolPlot> *histos, double etamin, double etamax, double ptmin, double ptmax, TString obj, double xmin = 0, double xmax = 2)
     160{
     161  resolPlot ptemp;
     162  double width;
     163  double etadown;
     164  double etaup;
     165
     166  for (int i = 0; i < Nbins; i++)
     167  {
     168    width = (etamax - etamin) / Nbins;
     169    etadown = etamin + i * width;
     170    etaup = etamin + (i+1) * width;
     171
     172    ptemp.set(etadown, etaup, ptmin, ptmax, obj, xmin, xmax);
     173    histos->push_back(ptemp);
     174  }
     175}
     176
     177
    138178//------------------------------------------------------------------------------
    139179
    … …  
    165205
    166206template<typename T>
    167 std::pair<TH1D*, TH1D*> GetEff(TClonesArray *branchReco, TClonesArray *branchParticle, TString name, int pdgID, ExRootTreeReader *treeReader)
     207TH1D* GetEffPt(TClonesArray *branchReco, TClonesArray *branchParticle, TString name, int pdgID, double ptmin, double ptmax, double etamin, double etamax, ExRootTreeReader *treeReader)
    168208{
    169209
    … …  
    183223  Int_t i, j;
    184224
    185   TH1D *histGenPtcen = new TH1D(name+" gen spectra Pt",name+" gen spectra cen", Nbins, TMath::Log10(ptrangemin), TMath::Log10(ptrangemax));
    186   TH1D *histRecoPtcen = new TH1D(name+" reco spectra Pt",name+" reco spectra cen", Nbins, TMath::Log10(ptrangemin), TMath::Log10(ptrangemax));
    187   TH1D *histGenPtfwd  = new TH1D(name+" gen spectra Eta",name+" gen spectra fwd", Nbins, TMath::Log10(ptrangemin), TMath::Log10(ptrangemax));
    188   TH1D *histRecoPtfwd = new TH1D(name+" reco spectra Eta",name+" reco spectra fwd", Nbins, TMath::Log10(ptrangemin), TMath::Log10(ptrangemax));
    189 
    190   histGenPtcen->SetDirectory(0);
    191   histRecoPtcen->SetDirectory(0);
    192   histGenPtfwd->SetDirectory(0);
    193   histRecoPtfwd->SetDirectory(0);
    194 
    195   BinLogX(histGenPtcen);
    196   BinLogX(histRecoPtcen);
    197   BinLogX(histGenPtfwd);
    198   BinLogX(histRecoPtfwd);
     225  TH1D *histGenPt = new TH1D(name+" gen spectra Pt",name+" gen spectra cen", Nbins, TMath::Log10(ptmin), TMath::Log10(ptmax));
     226  TH1D *histRecoPt = new TH1D(name+" reco spectra Pt",name+" reco spectra cen", Nbins, TMath::Log10(ptmin), TMath::Log10(ptmax));
     227
     228  histGenPt->SetDirectory(0);
     229  histRecoPt->SetDirectory(0);
     230
     231  BinLogX(histGenPt);
     232  BinLogX(histRecoPt);
    199233
    200234  // Loop over all events
    … …  
    213247      deltaR = 999;
    214248
    215       if (particle->PID == pdgID && genMomentum.Pt() > ptrangemin && genMomentum.Pt() < ptrangemax )
     249      pt  = genMomentum.Pt();
     250      eta = TMath::Abs(genMomentum.Eta());
     251
     252      if(eta > etamax || eta < etamin ) continue;
     253
     254      if (particle->PID == pdgID && genMomentum.Pt() > ptmin && genMomentum.Pt() < ptmax )
    216255      {
    217 
    218256        // Loop over all reco object in event
    219257        for(j = 0; j < branchReco->GetEntriesFast(); ++j)
    … …  
    229267          {
    230268            Jet *jet = (Jet *)recoObj;
    231             if(jet->BTag != 1) continue;
     269            if( !(jet->BTag & (1 << 0)) ) continue;
     270
     271            //if(jet->BTag != ) continue;
    232272          }
     273
     274            if(TMath::Abs(pdgID) == 4)
     275          {
     276            Jet *jet = (Jet *)recoObj;
     277            if( !(jet->BTag & (1 << 0)) ) continue;
     278          }
     279
     280           if(TMath::Abs(pdgID) == 1)
     281          {
     282            Jet *jet = (Jet *)recoObj;
     283            if( !(jet->BTag & (1 << 0)) ) continue;
     284          }
     285
     286          if(TMath::Abs(pdgID) == 15)
     287          {
     288            Jet *jet = (Jet *)recoObj;
     289            if(jet->TauTag != 1) continue;
     290          }
     291
     292
     293          if(genMomentum.DeltaR(recoMomentum) < deltaR)
     294          {
     295            deltaR = genMomentum.DeltaR(recoMomentum);
     296            bestRecoMomentum = recoMomentum;
     297          }
     298        }
     299
     300        histGenPt->Fill(pt);
     301        if(deltaR < 0.3) { histRecoPt->Fill(pt); }
     302
     303      }
     304    }
     305  }
     306
     307  histRecoPt->Sumw2();
     308  histGenPt->Sumw2();
     309
     310  histRecoPt->Divide(histGenPt);
     311  histRecoPt->Scale(100.);
     312
     313  return histRecoPt;
     314}
     315
     316template<typename T>
     317TH1D* GetEffEta(TClonesArray *branchReco, TClonesArray *branchParticle, TString name, int pdgID, double ptmin, double ptmax, double etamin, double etamax, ExRootTreeReader *treeReader)
     318{
     319
     320  cout << "** Computing Efficiency of reconstructing "<< branchReco->GetName() << " induced by " << branchParticle->GetName() << " with PID " << pdgID << endl;
     321
     322  Long64_t allEntries = treeReader->GetEntries();
     323
     324  GenParticle *particle;
     325  T *recoObj;
     326
     327  TLorentzVector recoMomentum, genMomentum, bestRecoMomentum;
     328
     329  Float_t deltaR;
     330  Float_t pt, eta;
     331  Long64_t entry;
     332
     333  Int_t i, j;
     334
     335  TH1D *histGenEta = new TH1D(name+" gen spectra Eta",name+" gen spectra", Nbins, etamin, etamax);
     336  TH1D *histRecoEta = new TH1D(name+" reco spectra Eta",name+" reco spectra", Nbins, etamin, etamax);
     337
     338  histGenEta->SetDirectory(0);
     339  histRecoEta->SetDirectory(0);
     340
     341  // Loop over all events
     342  for(entry = 0; entry < allEntries; ++entry)
     343  {
     344    // Load selected branches with data from specified event
     345    treeReader->ReadEntry(entry);
     346
     347    // Loop over all generated particle in event
     348    for(i = 0; i < branchParticle->GetEntriesFast(); ++i)
     349    {
     350
     351      particle = (GenParticle*) branchParticle->At(i);
     352      genMomentum = particle->P4();
     353
     354      deltaR = 999;
     355
     356      pt  = genMomentum.Pt();
     357      eta = genMomentum.Eta();
     358
     359      if(pt > ptmax || pt < ptmin ) continue;
     360
     361      if (particle->PID == pdgID && genMomentum.Pt() > ptmin && genMomentum.Pt() < ptmax )
     362      {
     363        // Loop over all reco object in event
     364        for(j = 0; j < branchReco->GetEntriesFast(); ++j)
     365        {
     366          recoObj = (T*)branchReco->At(j);
     367          recoMomentum = recoObj->P4();
     368          // this is simply to avoid warnings from initial state particle
     369          // having infite rapidity ...
     370        //if(Momentum.Px() == 0 && genMomentum.Py() == 0) continue;
     371
     372          // take the closest parton candidate
     373          if(TMath::Abs(pdgID) == 5)
     374          {
     375            Jet *jet = (Jet *)recoObj;
     376            if( !(jet->BTag & (1 << 0)) ) continue;
     377
     378          }
     379
     380          if(TMath::Abs(pdgID) == 4)
     381          {
     382            Jet *jet = (Jet *)recoObj;
     383            if( !(jet->BTag & (1 << 0)) ) continue;
     384          }
     385
     386           if(TMath::Abs(pdgID) == 1)
     387          {
     388            Jet *jet = (Jet *)recoObj;
     389            if( !(jet->BTag & (1 << 0)) ) continue;
     390          }
     391
    233392          if(TMath::Abs(pdgID) == 15)
    234393          {
    … …  
    243402        }
    244403
    245         pt  = genMomentum.Pt();
    246         eta = genMomentum.Eta();
    247 
    248         if (TMath::Abs(eta) < 1.5)
     404        histGenEta->Fill(eta);
     405        if(deltaR < 0.3) { histRecoEta->Fill(eta); }
     406
     407      }
     408    }
     409  }
     410
     411  histRecoEta->Sumw2();
     412  histGenEta->Sumw2();
     413
     414  histRecoEta->Divide(histGenEta);
     415  histRecoEta->Scale(100.);
     416
     417  return histRecoEta;
     418}
     419
     420
     421
     422template<typename T>
     423TH1D* GetTauEffPt(TClonesArray *branchReco, TClonesArray *branchParticle, TString name, int pdgID, double ptmin, double ptmax, double etamin, double etamax, ExRootTreeReader *treeReader)
     424{
     425
     426  cout << "** Computing Efficiency of reconstructing "<< branchReco->GetName() << " induced by " << branchParticle->GetName() << " with PID " << pdgID << endl;
     427
     428  Long64_t allEntries = treeReader->GetEntries();
     429
     430  GenParticle *particle;
     431  T *recoObj;
     432
     433  TLorentzVector recoMomentum, genMomentum, bestRecoMomentum;
     434
     435  Float_t deltaR;
     436  Float_t pt, eta;
     437  Long64_t entry;
     438
     439  Int_t i, j;
     440
     441  TH1D *histGenPt = new TH1D(name+" gen spectra Pt",name+" gen spectra cen", Nbins, TMath::Log10(ptmin), TMath::Log10(ptmax));
     442  TH1D *histRecoPt = new TH1D(name+" reco spectra Pt",name+" reco spectra cen", Nbins, TMath::Log10(ptmin), TMath::Log10(ptmax));
     443
     444  histGenPt->SetDirectory(0);
     445  histRecoPt->SetDirectory(0);
     446
     447  BinLogX(histGenPt);
     448  BinLogX(histRecoPt);
     449
     450  // Loop over all events
     451  for(entry = 0; entry < allEntries; ++entry)
     452  {
     453    // Load selected branches with data from specified event
     454    treeReader->ReadEntry(entry);
     455
     456    // Loop over all generated particle in event
     457    for(i = 0; i < branchParticle->GetEntriesFast(); ++i)
     458    {
     459
     460      particle = (GenParticle*) branchParticle->At(i);
     461      genMomentum = particle->P4();
     462
     463      deltaR = 999;
     464
     465      pt  = genMomentum.Pt();
     466      eta = TMath::Abs(genMomentum.Eta());
     467
     468      if(eta > etamax || eta < etamin ) continue;
     469
     470      if (particle->PID == pdgID && genMomentum.Pt() > ptmin && genMomentum.Pt() < ptmax )
     471      {
     472        // Loop over all reco object in event
     473        for(j = 0; j < branchReco->GetEntriesFast(); ++j)
    249474        {
    250           histGenPtcen->Fill(pt);
    251           if(deltaR < 0.3) { histRecoPtcen->Fill(pt); }
     475          recoObj = (T*)branchReco->At(j);
     476          recoMomentum = recoObj->P4();
     477          // this is simply to avoid warnings from initial state particle
     478          // having infite rapidity ...
     479          //if(Momentum.Px() == 0 && genMomentum.Py() == 0) continue;
     480
     481          if(TMath::Abs(pdgID) == 1)
     482          {
     483            Jet *jet = (Jet *)recoObj;
     484            if( jet->TauTag != 1 ) continue;
     485          }
     486
     487          if(TMath::Abs(pdgID) == 15)
     488          {
     489            Jet *jet = (Jet *)recoObj;
     490            if(jet->TauTag != 1) continue;
     491          }
     492
     493          if(genMomentum.DeltaR(recoMomentum) < deltaR)
     494          {
     495            deltaR = genMomentum.DeltaR(recoMomentum);
     496            bestRecoMomentum = recoMomentum;
     497          }
    252498        }
    253         else if (TMath::Abs(eta) < 2.5)
     499
     500        histGenPt->Fill(pt);
     501        if(deltaR < 0.3) { histRecoPt->Fill(pt); }
     502
     503      }
     504    }
     505  }
     506
     507  histRecoPt->Sumw2();
     508  histGenPt->Sumw2();
     509
     510
     511  histRecoPt->Divide(histGenPt);
     512  histRecoPt->Scale(100.);
     513  if(TMath::Abs(pdgID) == 15) histRecoPt->Scale(1/0.648);
     514
     515  return histRecoPt;
     516}
     517
     518template<typename T>
     519TH1D* GetTauEffEta(TClonesArray *branchReco, TClonesArray *branchParticle, TString name, int pdgID, double ptmin, double ptmax, double etamin, double etamax, ExRootTreeReader *treeReader)
     520{
     521
     522  cout << "** Computing Efficiency of reconstructing "<< branchReco->GetName() << " induced by " << branchParticle->GetName() << " with PID " << pdgID << endl;
     523
     524  Long64_t allEntries = treeReader->GetEntries();
     525
     526  GenParticle *particle;
     527  T *recoObj;
     528
     529  TLorentzVector recoMomentum, genMomentum, bestRecoMomentum;
     530
     531  Float_t deltaR;
     532  Float_t pt, eta;
     533  Long64_t entry;
     534
     535  Int_t i, j;
     536
     537  TH1D *histGenEta = new TH1D(name+" gen spectra Eta",name+" gen spectra", Nbins, etamin, etamax);
     538  TH1D *histRecoEta = new TH1D(name+" reco spectra Eta",name+" reco spectra", Nbins, etamin, etamax);
     539
     540  histGenEta->SetDirectory(0);
     541  histRecoEta->SetDirectory(0);
     542
     543  // Loop over all events
     544  for(entry = 0; entry < allEntries; ++entry)
     545  {
     546    // Load selected branches with data from specified event
     547    treeReader->ReadEntry(entry);
     548
     549    // Loop over all generated particle in event
     550    for(i = 0; i < branchParticle->GetEntriesFast(); ++i)
     551    {
     552
     553      particle = (GenParticle*) branchParticle->At(i);
     554      genMomentum = particle->P4();
     555
     556      deltaR = 999;
     557
     558      pt  = genMomentum.Pt();
     559      eta = genMomentum.Eta();
     560
     561      if(pt > ptmax || pt < ptmin ) continue;
     562
     563      if (particle->PID == pdgID && genMomentum.Pt() > ptmin && genMomentum.Pt() < ptmax )
     564      {
     565        // Loop over all reco object in event
     566        for(j = 0; j < branchReco->GetEntriesFast(); ++j)
    254567        {
    255           histGenPtfwd->Fill(pt);
    256           if(deltaR < 0.3) { histRecoPtfwd->Fill(pt); }
    257 
     568          recoObj = (T*)branchReco->At(j);
     569          recoMomentum = recoObj->P4();
     570          // this is simply to avoid warnings from initial state particle
     571          // having infite rapidity ...
     572        //if(Momentum.Px() == 0 && genMomentum.Py() == 0) continue;
     573
     574         if(TMath::Abs(pdgID) == 1)
     575          {
     576            Jet *jet = (Jet *)recoObj;
     577            if( jet->TauTag != 1 ) continue;
     578          }
     579
     580          if(TMath::Abs(pdgID) == 15)
     581          {
     582            Jet *jet = (Jet *)recoObj;
     583            if(jet->TauTag != 1) continue;
     584          }
     585
     586          if(genMomentum.DeltaR(recoMomentum) < deltaR)
     587          {
     588            deltaR = genMomentum.DeltaR(recoMomentum);
     589            bestRecoMomentum = recoMomentum;
     590          }
    258591        }
     592
     593        histGenEta->Fill(eta);
     594        if(deltaR < 0.3) { histRecoEta->Fill(eta); }
     595
    259596      }
    260597    }
    261598  }
    262599
    263 
    264   std::pair<TH1D*,TH1D*> histos;
    265 
    266   histRecoPtcen->Divide(histGenPtcen);
    267   histRecoPtfwd->Divide(histGenPtfwd);
    268 
    269   histos.first = histRecoPtcen;
    270   histos.second = histRecoPtfwd;
    271 
    272   return histos;
     600  histRecoEta->Sumw2();
     601  histGenEta->Sumw2();
     602
     603  histRecoEta->Divide(histGenEta);
     604  histRecoEta->Scale(100.);
     605  if(TMath::Abs(pdgID) == 15) histRecoEta->Scale(1/0.648);
     606
     607  return histRecoEta;
    273608}
    274609
     610
    275611template<typename T>
    276 void GetEres(std::vector<resolPlot> *histos, TClonesArray *branchReco, TClonesArray *branchParticle, int pdgID, ExRootTreeReader *treeReader)
     612void GetPtres(std::vector<resolPlot> *histos, TClonesArray *branchReco, TClonesArray *branchParticle, int pdgID, Double_t etaMin, Double_t etaMax, ExRootTreeReader *treeReader)
    277613{
    278614  Long64_t allEntries = treeReader->GetEntries();
    279615
    280   cout << "** Computing resolution of " << branchReco->GetName() << " induced by " << branchParticle->GetName() << " with PID " << pdgID << endl;
     616  cout << "** Computing pt resolution of " << branchReco->GetName() << " induced by " << branchParticle->GetName() << " with PID " << pdgID << endl;
    281617
    282618  GenParticle *particle;
    … …  
    328664      if(deltaR < 0.3)
    329665      {
    330         pt  = bestGenMomentum.E();
     666        pt  = bestGenMomentum.Pt();
    331667        eta = TMath::Abs(bestGenMomentum.Eta());
    332668
    333669        for (bin = 0; bin < Nbins; bin++)
    334670        {
    335           if(pt > histos->at(bin).ptmin && pt < histos->at(bin).ptmax && eta < 2.5)
     671          if(pt > histos->at(bin).ptmin && pt < histos->at(bin).ptmax && eta > etaMin && eta < etaMax)
    336672          {
    337             if (eta < 1.5) {histos->at(bin).cenResolHist->Fill(recoMomentum.E()/bestGenMomentum.E());}
    338             else if (eta < 2.5) {histos->at(bin).fwdResolHist->Fill(recoMomentum.E()/bestGenMomentum.E());}
     673            histos->at(bin).resolHist->Fill(recoMomentum.Pt()/bestGenMomentum.Pt());
    339674          }
    340675        }
    … …  
    346681
    347682template<typename T>
    348 void GetPtres(std::vector<resolPlot> *histos, TClonesArray *branchReco, TClonesArray *branchParticle, int pdgID, ExRootTreeReader *treeReader)
     683void GetEres(std::vector<resolPlot> *histos, TClonesArray *branchReco, TClonesArray *branchParticle, int pdgID, Double_t etaMin, Double_t etaMax, ExRootTreeReader *treeReader)
    349684{
    350685  Long64_t allEntries = treeReader->GetEntries();
    351686
    352   cout << "** Computing pt resolution of " << branchReco->GetName() << " induced by " << branchParticle->GetName() << " with PID " << pdgID << endl;
     687  cout << "** Computing e resolution of " << branchReco->GetName() << " induced by " << branchParticle->GetName() << " with PID " << pdgID << endl;
    353688
    354689  GenParticle *particle;
    … …  
    358693
    359694  Float_t deltaR;
    360   Float_t pt, eta;
     695  Float_t e, eta;
    361696  Long64_t entry;
    362697
    … …  
    369704    treeReader->ReadEntry(entry);
    370705
    371     // Loop over all reconstructed jets in event
     706     // Loop over all reconstructed jets in event
    372707    for(i = 0; i < branchReco->GetEntriesFast(); ++i)
    373708    {
    … …  
    400735      if(deltaR < 0.3)
    401736      {
    402         pt  = bestGenMomentum.Pt();
     737        e  = bestGenMomentum.E();
    403738        eta = TMath::Abs(bestGenMomentum.Eta());
    404739
    405740        for (bin = 0; bin < Nbins; bin++)
    406741        {
    407           if(pt > histos->at(bin).ptmin && pt < histos->at(bin).ptmax && eta < 2.5)
     742          if(e > histos->at(bin).ptmin && e < histos->at(bin).ptmax && eta > etaMin && eta < etaMax)
    408743          {
    409             if (eta < 1.5) {histos->at(bin).cenResolHist->Fill(recoMomentum.Pt()/bestGenMomentum.Pt());}
    410             else if (eta < 2.5) {histos->at(bin).fwdResolHist->Fill(recoMomentum.Pt()/bestGenMomentum.Pt());}
     744            histos->at(bin).resolHist->Fill(recoMomentum.E()/bestGenMomentum.E());
    411745          }
    412746        }
    … …  
    417751
    418752
    419 void GetJetsEres(std::vector<resolPlot> *histos, TClonesArray *branchJet, TClonesArray *branchGenJet, ExRootTreeReader *treeReader)
     753template<typename T>
     754void GetPtresVsEta(std::vector<resolPlot> *histos, TClonesArray *branchReco, TClonesArray *branchParticle, int pdgID, Double_t ptMin, Double_t ptMax, ExRootTreeReader *treeReader)
     755{
     756  Long64_t allEntries = treeReader->GetEntries();
     757
     758  cout << "** Computing pt resolution of " << branchReco->GetName() << " induced by " << branchParticle->GetName() << " with PID " << pdgID << endl;
     759
     760  GenParticle *particle;
     761  T* recoObj;
     762
     763  TLorentzVector recoMomentum, genMomentum, bestGenMomentum;
     764
     765  Float_t deltaR;
     766  Float_t pt, eta;
     767  Long64_t entry;
     768
     769  Int_t i, j, bin;
     770
     771  // Loop over all events
     772  for(entry = 0; entry < allEntries; ++entry)
     773  {
     774    // Load selected branches with data from specified event
     775    treeReader->ReadEntry(entry);
     776
     777    // Loop over all reconstructed jets in event
     778    for(i = 0; i < branchReco->GetEntriesFast(); ++i)
     779    {
     780      recoObj = (T*) branchReco->At(i);
     781      recoMomentum = recoObj->P4();
     782
     783      deltaR = 999;
     784
     785     // Loop over all hard partons in event
     786     for(j = 0; j < branchParticle->GetEntriesFast(); ++j)
     787     {
     788        particle = (GenParticle*) branchParticle->At(j);
     789        if (particle->PID == pdgID && particle->Status == 1)
     790        {
     791          genMomentum = particle->P4();
     792
     793          // this is simply to avoid warnings from initial state particle
     794          // having infite rapidity ...
     795          if(genMomentum.Px() == 0 && genMomentum.Py() == 0) continue;
     796
     797          // take the closest parton candidate
     798          if(genMomentum.DeltaR(recoMomentum) < deltaR)
     799          {
     800             deltaR = genMomentum.DeltaR(recoMomentum);
     801             bestGenMomentum = genMomentum;
     802          }
     803        }
     804      }
     805
     806      if(deltaR < 0.3)
     807      {
     808        pt  = bestGenMomentum.Pt();
     809        eta = bestGenMomentum.Eta();
     810
     811        for (bin = 0; bin < Nbins; bin++)
     812        {
     813          if(eta > histos->at(bin).etamin && eta < histos->at(bin).etamax && pt > ptMin && pt < ptMax)
     814          {
     815             histos->at(bin).resolHist->Fill(recoMomentum.Pt()/bestGenMomentum.Pt());
     816          }
     817        }
     818      }
     819    }
     820  }
     821}
     822
     823template<typename T>
     824void GetEresVsEta(std::vector<resolPlot> *histos, TClonesArray *branchReco, TClonesArray *branchParticle, int pdgID, Double_t eMin, Double_t eMax, ExRootTreeReader *treeReader)
     825{
     826  Long64_t allEntries = treeReader->GetEntries();
     827
     828  cout << "** Computing E resolution of " << branchReco->GetName() << " induced by " << branchParticle->GetName() << " with PID " << pdgID << endl;
     829
     830  GenParticle *particle;
     831  T* recoObj;
     832
     833  TLorentzVector recoMomentum, genMomentum, bestGenMomentum;
     834
     835  Float_t deltaR;
     836  Float_t e, eta;
     837  Long64_t entry;
     838
     839  Int_t i, j, bin;
     840
     841  // Loop over all events
     842  for(entry = 0; entry < allEntries; ++entry)
     843  {
     844    // Load selected branches with data from specified event
     845    treeReader->ReadEntry(entry);
     846
     847    // Loop over all reconstructed jets in event
     848    for(i = 0; i < branchReco->GetEntriesFast(); ++i)
     849    {
     850      recoObj = (T*) branchReco->At(i);
     851      recoMomentum = recoObj->P4();
     852
     853      deltaR = 999;
     854
     855     // Loop over all hard partons in event
     856     for(j = 0; j < branchParticle->GetEntriesFast(); ++j)
     857     {
     858        particle = (GenParticle*) branchParticle->At(j);
     859        if (particle->PID == pdgID && particle->Status == 1)
     860        {
     861          genMomentum = particle->P4();
     862
     863          // this is simply to avoid warnings from initial state particle
     864          // having infite rapidity ...
     865          if(genMomentum.Px() == 0 && genMomentum.Py() == 0) continue;
     866
     867          // take the closest parton candidate
     868          if(genMomentum.DeltaR(recoMomentum) < deltaR)
     869          {
     870             deltaR = genMomentum.DeltaR(recoMomentum);
     871             bestGenMomentum = genMomentum;
     872          }
     873        }
     874      }
     875
     876      if(deltaR < 0.3)
     877      {
     878        e  = bestGenMomentum.E();
     879        eta = bestGenMomentum.Eta();
     880
     881
     882
     883        for (bin = 0; bin < Nbins; bin++)
     884        {
     885          if(eta > histos->at(bin).etamin && eta < histos->at(bin).etamax && e > eMin && e < eMax)
     886          {
     887             histos->at(bin).resolHist->Fill(recoMomentum.E()/bestGenMomentum.E());
     888          }
     889        }
     890      }
     891    }
     892  }
     893}
     894
     895
     896
     897void GetJetsEres(std::vector<resolPlot> *histos, TClonesArray *branchJet, TClonesArray *branchGenJet, ExRootTreeReader *treeReader, Double_t etaMin, Double_t etaMax)
    420898{
    421899
    … …  
    470948      }
    471949
    472       if(deltaR < 0.25)
     950      if(deltaR < 0.3)
    473951      {
    474952        pt  = genJetMomentum.E();
    475         eta = TMath::Abs(genJetMomentum.Eta());
     953        eta = genJetMomentum.Eta();
    476954
    477955        for (bin = 0; bin < Nbins; bin++)
    478956        {
    479             if(pt > histos->at(bin).ptmin && pt < histos->at(bin).ptmax && eta < 1.5)
    480             {
    481                 histos->at(bin).cenResolHist->Fill(jetMomentum.E()/bestGenJetMomentum.E());
    482             }
    483             else if(pt > histos->at(bin).ptmin && pt < histos->at(bin).ptmax && eta < 2.5)
    484             {
    485                 histos->at(bin).fwdResolHist->Fill(jetMomentum.E()/bestGenJetMomentum.E());
    486             }
     957          if(pt > histos->at(bin).ptmin && pt < histos->at(bin).ptmax && eta < etaMax && eta > etaMin)
     958          {
     959             histos->at(bin).resolHist->Fill(jetMomentum.E()/bestGenJetMomentum.E());
     960          }
    487961        }
    488962      }
    489963    }
    490964  }
     965}
     966
     967void GetJetsEresVsEta(std::vector<resolPlot> *histos, TClonesArray *branchJet, TClonesArray *branchGenJet, ExRootTreeReader *treeReader, Double_t eMin, Double_t eMax)
     968{
     969
     970  Long64_t allEntries = treeReader->GetEntries();
     971
     972  cout << "** Computing resolution of " << branchJet->GetName() << " induced by " << branchGenJet->GetName() << endl;
     973
     974  Jet *jet, *genjet;
     975
     976  TLorentzVector jetMomentum, genJetMomentum, bestGenJetMomentum;
     977
     978  Float_t deltaR;
     979  Float_t pt, eta;
     980  Long64_t entry;
     981
     982  Int_t i, j, bin;
     983
     984  // Loop over all events
     985  for(entry = 0; entry < allEntries; ++entry)
     986  {
     987    // Load selected branches with data from specified event
     988    treeReader->ReadEntry(entry);
     989
     990    if(entry%10000 == 0) cout << "Event number: "<< entry <<endl;
     991
     992    // Loop over all reconstructed jets in event
     993    for(i = 0; i < TMath::Min(2,branchJet->GetEntriesFast()); ++i) //branchJet->GetEntriesFast(); ++i)
     994    {
     995
     996      jet = (Jet*) branchJet->At(i);
     997      jetMomentum = jet->P4();
     998
     999      deltaR = 999;
     1000
     1001     // Loop over all hard partons in event
     1002     for(j = 0; j < TMath::Min(2,branchGenJet->GetEntriesFast()); ++j)
     1003     {
     1004        genjet = (Jet*) branchGenJet->At(j);
     1005
     1006        genJetMomentum = genjet->P4();
     1007
     1008        // this is simply to avoid warnings from initial state particle
     1009        // having infite rapidity ...
     1010        if(genJetMomentum.Px() == 0 && genJetMomentum.Py() == 0) continue;
     1011
     1012        // take the closest parton candidate
     1013        if(genJetMomentum.DeltaR(jetMomentum) < deltaR)
     1014        {
     1015           deltaR = genJetMomentum.DeltaR(jetMomentum);
     1016           bestGenJetMomentum = genJetMomentum;
     1017        }
     1018      }
     1019
     1020      if(deltaR < 0.3)
     1021      {
     1022
     1023        pt  = genJetMomentum.E();
     1024        eta = genJetMomentum.Eta();
     1025
     1026        for (bin = 0; bin < Nbins; bin++)
     1027        {
     1028          if(eta > histos->at(bin).etamin && eta < histos->at(bin).etamax && pt < eMax && pt > eMin)
     1029          {
     1030             histos->at(bin).resolHist->Fill(jetMomentum.E()/bestGenJetMomentum.E());
     1031          }
     1032        }
     1033      }
     1034    }
     1035  }
     1036}
     1037
     1038
     1039
     1040std::pair<Double_t, Double_t> GausFit(TH1* hist)
     1041{
     1042  TF1 *f1 = new TF1("f1", "gaus", hist->GetMean()-2*hist->GetRMS(), hist->GetMean()+2*hist->GetRMS());
     1043  hist->Fit("f1","RQ");
     1044
     1045  TF1 *f2 = new TF1("f2", "gaus", f1->GetParameter(1) - 2*f1->GetParameter(2), f1->GetParameter(1) + 2*f1->GetParameter(2));
     1046  hist->Fit("f2","RQ");
     1047
     1048  Double_t sig = f2->GetParameter(2);
     1049  Double_t sigErr = f2->GetParError(2);
     1050
     1051  delete f1;
     1052  delete f2;
     1053  return make_pair (sig, sigErr);
     1054}
     1055
     1056
     1057TGraphErrors EresGraph(std::vector<resolPlot> *histos, bool rms = false)
     1058{
     1059  Int_t bin;
     1060  Int_t count = 0;
     1061  TGraphErrors gr = TGraphErrors(Nbins/2);
     1062  double val, error;
     1063  for (bin = 0; bin < Nbins; bin++)
     1064  {
     1065     std::pair<Double_t, Double_t> sigvalues = GausFit(histos->at(bin).resolHist);
     1066     if (rms == true)
     1067     {
     1068       gr.SetPoint(count,(histos->at(bin).ptmin+histos->at(bin).ptmax)/2.0, 100*histos->at(bin).resolHist->GetRMS());
     1069       //gr.SetPointError(count,0, 100*sigvalues.second); // to correct
     1070       error = 100*histos->at(bin).resolHist->GetRMSError();
     1071       val = 100*histos->at(bin).resolHist->GetRMS();
     1072       if(error > 0.2*val) error = 0.2*val;
     1073       gr.SetPointError(count,0, error); // to correct
     1074     }
     1075     else
     1076     {
     1077
     1078       gr.SetPoint(count,(histos->at(bin).ptmin+histos->at(bin).ptmax)/2.0, 100*sigvalues.first);
     1079       error = 100*sigvalues.second;
     1080       val = 100*sigvalues.first;
     1081       if(error > 0.2*val) error = 0.2*val;
     1082       gr.SetPointError(count,0, error); // to correct
     1083       //gr.SetPointError(count,0, 100*sigvalues.second);
     1084     }
     1085     count++;
     1086  }
     1087
     1088  return gr;
     1089}
     1090
     1091TGraphErrors MetResGraph(std::vector<resolPlot> *histos, bool rms = false)
     1092{
     1093  Int_t bin;
     1094  Int_t count = 0;
     1095  TGraphErrors gr = TGraphErrors(Nbins/2);
     1096  double val, error;
     1097  for (bin = 0; bin < Nbins; bin++)
     1098  {
     1099     std::pair<Double_t, Double_t> sigvalues = GausFit(histos->at(bin).resolHist);
     1100     if (rms == true)
     1101     {
     1102       gr.SetPoint(count,(histos->at(bin).ptmin+histos->at(bin).ptmax)/2.0, histos->at(bin).resolHist->GetRMS());
     1103       error = histos->at(bin).resolHist->GetRMSError();
     1104       val = histos->at(bin).resolHist->GetRMS();
     1105       if(error > 0.2*val) error = 0.2*val;
     1106       gr.SetPointError(count,0, error); // to correct
     1107     }
     1108     else
     1109     {
     1110
     1111       gr.SetPoint(count,(histos->at(bin).ptmin+histos->at(bin).ptmax)/2.0, sigvalues.first);
     1112       val = sigvalues.first;
     1113       error = sigvalues.second;
     1114       if(error > 0.2*val) error = 0.2*val;
     1115       gr.SetPointError(count,0, error);
     1116       //gr.SetPointError(count,0, 100*sigvalues.second);
     1117     }
     1118     count++;
     1119  }
     1120
     1121  return gr;
     1122}
     1123
     1124
     1125
     1126TGraphErrors EresGraphVsEta(std::vector<resolPlot> *histos, bool rms = false)
     1127{
     1128  Int_t bin;
     1129  Int_t count = 0;
     1130  TGraphErrors gr = TGraphErrors(Nbins/2);
     1131  double val, error;
     1132  for (bin = 0; bin < Nbins; bin++)
     1133  {
     1134
     1135     std::pair<Double_t, Double_t> sigvalues = GausFit(histos->at(bin).resolHist);
     1136     if (rms == true)
     1137     {
     1138       gr.SetPoint(count,(histos->at(bin).etamin+histos->at(bin).etamax)/2.0, histos->at(bin).resolHist->GetRMS());
     1139       error = 100*histos->at(bin).resolHist->GetRMSError();
     1140       val = 100*histos->at(bin).resolHist->GetRMS();
     1141       if(error > 0.2*val) error = 0.2*val;
     1142       gr.SetPointError(count,0, error); // to correct
     1143     }
     1144     else
     1145     {
     1146       gr.SetPoint(count,(histos->at(bin).etamin+histos->at(bin).etamax)/2.0, 100*sigvalues.first);
     1147       val = 100*sigvalues.first;
     1148       error = 100*sigvalues.second;
     1149       if(error > 0.2*val) error = 0.2*val;
     1150       gr.SetPointError(count,0, error);
     1151       //gr.SetPointError(count,0, 100*sigvalues.second);
     1152     }
     1153     count++;
     1154  }
     1155
     1156  return gr;
    4911157}
    4921158
    … …  
    5361202        if(ht > histos->at(bin).ptmin && ht < histos->at(bin).ptmax )
    5371203        {
    538           histos->at(bin).cenResolHist->Fill(met->P4().Px());
    539           histos->at(bin).fwdResolHist->Fill(met->P4().Py());
     1204          histos->at(bin).resolHist->Fill(met->P4().Px());
    5401205        }
    5411206      }
    … …  
    5451210
    5461211
    547 std::pair<Double_t, Double_t> GausFit(TH1* hist)
     1212
     1213//------------------------------------------------------------------------------
     1214
     1215
     1216void addResoGraph(TMultiGraph *mg, TGraphErrors *gr, TLegend *leg, int style, Color_t color, TString text)
    5481217{
    549   TF1 *f1 = new TF1("f1", "gaus", hist->GetMean()-2*hist->GetRMS(), hist->GetMean()+2*hist->GetRMS());
    550   hist->Fit("f1","RQ");
    551 
    552   TF1 *f2 = new TF1("f2", "gaus", f1->GetParameter(1) - 2*f1->GetParameter(2), f1->GetParameter(1) + 2*f1->GetParameter(2));
    553   hist->Fit("f2","RQ");
    554 
    555   Double_t sig = f2->GetParameter(2);
    556   Double_t sigErr = f2->GetParError(2);
    557 
    558   delete f1;
    559   delete f2;
    560   return make_pair (sig, sigErr);
     1218
     1219  gr->SetLineWidth(2);
     1220  gr->SetLineColor(color);
     1221  gr->SetMarkerStyle(style);
     1222  gr->SetMarkerColor(color);
     1223  gr->SetMarkerSize(1.5);
     1224
     1225  std::cout << "Adding " << gr->GetName() << std::endl;
     1226  mg->Add(gr);
     1227  leg->AddEntry(gr,text,"p");
     1228
    5611229}
    5621230
    5631231
    564 TGraphErrors EresGraph(std::vector<resolPlot> *histos, bool central, bool rms = false)
     1232void DrawAxis(TMultiGraph *mg, TLegend *leg, double xmin, double xmax, double ymin, double ymax, TString tx, TString ty, bool logx = 0, bool logy = 0)
    5651233{
    566   Int_t bin;
    567   Int_t count = 0;
    568   TGraphErrors gr = TGraphErrors(Nbins/2);
    569   Double_t sig = 0;
    570   Double_t sigErr = 0;
    571   for (bin = 0; bin < Nbins; bin++)
    572   {
    573     if (central == true && histos->at(bin).cenResolHist->GetEntries() > 100)
    574     {
    575       std::pair<Double_t, Double_t> sigvalues = GausFit(histos->at(bin).cenResolHist);
    576       if (rms == true)
    577       {
    578         gr.SetPoint(count,(histos->at(bin).ptmin+histos->at(bin).ptmax)/2.0, sigvalues.second);
    579         gr.SetPointError(count,0, sigvalues.second); // to correct
    580       }
    581       else
    582       {
    583         gr.SetPoint(count,(histos->at(bin).ptmin+histos->at(bin).ptmax)/2.0, sigvalues.first);
    584         gr.SetPointError(count,0, sigvalues.second);
    585       }
    586       count++;
    587     }
    588 
    589     else if (central == false && histos->at(bin).fwdResolHist->GetEntries() > 10)
    590     {
    591       std::pair<Double_t, Double_t> sigvalues = GausFit(histos->at(bin).fwdResolHist);
    592       if (rms == true)
    593       {
    594         gr.SetPoint(count,(histos->at(bin).ptmin+histos->at(bin).ptmax)/2.0, sigvalues.second);
    595         gr.SetPointError(count,0, sigvalues.second); // to correct
    596       }
    597       else
    598       {
    599         gr.SetPoint(count,(histos->at(bin).ptmin+histos->at(bin).ptmax)/2.0, sigvalues.first);
    600         gr.SetPointError(count,0, sigvalues.second);
    601       }
    602       count++;
    603     }
    604 
    605   }
    606   return gr;
    607 }
    608 
    609 
    610 //------------------------------------------------------------------------------
    611 
    612 
    613 // type 1 : object, 2 : track, 3 : tower
    614 
    615 void addGraph(TMultiGraph *mg, TGraphErrors *gr, TLegend *leg, int type)
    616 {
    617 
    618   gr->SetLineWidth(2);
    619 
    620   switch ( type )
    621   {
    622     case 1:
    623       gr->SetLineColor(objColor);
    624       gr->SetLineStyle(objStyle);
    625       std::cout << "Adding " << gr->GetName() << std::endl;
    626       mg->Add(gr);
    627       leg->AddEntry(gr,"Reco","l");
    628       break;
    629 
    630     case 2:
    631       gr->SetLineColor(trackColor);
    632       gr->SetLineStyle(trackStyle);
    633       mg->Add(gr);
    634       leg->AddEntry(gr,"Track","l");
    635       break;
    636 
    637     case 3:
    638       gr->SetLineColor(towerColor);
    639       gr->SetLineStyle(towerStyle);
    640       mg->Add(gr);
    641       leg->AddEntry(gr,"Tower","l");
    642       break;
    643 
    644     case 0:
    645       gr->SetLineColor(objColor);
    646       gr->SetLineStyle(objStyle);
    647       mg->Add(gr);
    648       break;
    649 
    650     default:
    651       std::cout << "wrong type, possibles choices are Object, Track and Tower" << std::endl;
    652       break;
    653   }
    654 }
    655 
    656 void addHist(TH1D *h, TLegend *leg, int type)
    657 {
    658   h->SetLineWidth(2);
    659 
    660   switch ( type )
    661   {
    662     case 1:
    663       h->SetLineColor(objColor);
    664       h->SetLineStyle(objStyle);
    665       leg->AddEntry(h,"Reco","l");
    666       break;
    667 
    668     case 2:
    669       h->SetLineColor(trackColor);
    670       h->SetLineStyle(trackStyle);
    671       leg->AddEntry(h,"Track","l");
    672       break;
    673 
    674     case 3:
    675       h->SetLineColor(towerColor);
    676       h->SetLineStyle(towerStyle);
    677       leg->AddEntry(h,"Tower","l");
    678       break;
    679 
    680     case 0:
    681       h->SetLineColor(objColor);
    682       h->SetLineStyle(objStyle);
    683       break;
    684 
    685     default:
    686       std::cout << "wrong type, possibles choices are Object, Track and Tower" << std::endl;
    687       break;
    688   }
    689 }
    690 
    691 void DrawAxis(TMultiGraph *mg, TLegend *leg, double max, int type = 0)
    692 {
    693   mg->SetMinimum(0.);
    694   mg->SetMaximum(max);
    695   mg->GetXaxis()->SetLimits(ptrangemin,ptrangemax);
    696   mg->GetYaxis()->SetTitle("resolution");
    697   if (type == 0) mg->GetXaxis()->SetTitle("E [GeV]");
    698   else mg->GetXaxis()->SetTitle("p_{T} [GeV]");
     1234  mg->SetMinimum(ymin);
     1235  mg->SetMaximum(ymax);
     1236  mg->GetXaxis()->SetLimits(xmin,xmax);
     1237
     1238  mg->GetXaxis()->SetTitle(tx);
     1239  mg->GetYaxis()->SetTitle(ty);
     1240
    6991241  mg->GetYaxis()->SetTitleSize(0.07);
    7001242  mg->GetXaxis()->SetTitleSize(0.07);
    … …  
    7051247  mg->GetXaxis()->SetTitleOffset(1.4);
    7061248
     1249  mg->GetXaxis()->SetTitleFont(132);
     1250  mg->GetYaxis()->SetTitleFont(132);
     1251  mg->GetXaxis()->SetLabelFont(132);
     1252  mg->GetYaxis()->SetLabelFont(132);
     1253
    7071254  mg->GetYaxis()->SetNdivisions(505);
    7081255
     1256  leg->SetTextFont(132);
    7091257  leg->SetBorderSize(0);
    7101258  leg->SetShadowColor(0);
    … …  
    7131261
    7141262  gStyle->SetOptTitle(0);
    715   gPad->SetLogx();
     1263
     1264  if(logx) gPad->SetLogx();
     1265  if(logy) gPad->SetLogy();
     1266
     1267  //gPad->SetGridx();
     1268  //gPad->SetGridy();
    7161269  gPad->SetBottomMargin(0.2);
    7171270  gPad->SetLeftMargin(0.2);
    … …  
    7211274}
    7221275
    723 void DrawAxis(TH1D *h, TLegend *leg, int type = 0)
     1276
     1277void Validation(const char *inputFilePion, const char *inputFileElectron, const char *inputFileMuon, const char *inputFilePhoton, const char *inputFileNeutralHadron, const char *inputFileJet, const char *inputFileBJet, const char *inputFileCJet, const char *inputFileTauJet, const char *outputFile)
    7241278{
    7251279
    726   h->GetYaxis()->SetRangeUser(0,1.0);
    727   if (type == 0) h->GetXaxis()->SetTitle("E [GeV]");
    728   else h->GetXaxis()->SetTitle("p_{T} [GeV]");
    729   h->GetYaxis()->SetTitle("efficiency");
    730   h->GetYaxis()->SetTitleSize(0.07);
    731   h->GetXaxis()->SetTitleSize(0.07);
    732   h->GetYaxis()->SetLabelSize(0.06);
    733   h->GetXaxis()->SetLabelSize(0.06);
    734   h->GetYaxis()->SetLabelOffset(0.03);
    735   h->GetYaxis()->SetTitleOffset(1.3);
    736   h->GetXaxis()->SetTitleOffset(1.4);
    737 
    738   h->GetYaxis()->SetNdivisions(505);
    739 
    740   leg->SetBorderSize(0);
    741   leg->SetShadowColor(0);
    742   leg->SetFillColor(0);
    743   leg->SetFillStyle(0);
    744 
    745   gStyle->SetOptTitle(0);
    746   gStyle->SetOptStat(0);
    747   gPad->SetBottomMargin(0.2);
    748   gPad->SetLeftMargin(0.2);
    749 
    750   gPad->Modified();
    751   gPad->Update();
    752 
    753 }
    754 
    755 
    756 void Validation(const char *inputFileElectron, const char *inputFileMuon, const char *inputFilePhoton, const char *inputFileJet, const char *inputFileBJet, const char *inputFileTauJet, const char *outputFile)
    757 {
     1280  TChain *chainPion = new TChain("Delphes");
     1281  chainPion->Add(inputFilePion);
     1282  ExRootTreeReader *treeReaderPion = new ExRootTreeReader(chainPion);
     1283
    7581284  TChain *chainElectron = new TChain("Delphes");
    7591285  chainElectron->Add(inputFileElectron);
    … …  
    7681294  ExRootTreeReader *treeReaderPhoton = new ExRootTreeReader(chainPhoton);
    7691295
     1296  TChain *chainNeutralHadron = new TChain("Delphes");
     1297  chainNeutralHadron->Add(inputFileNeutralHadron);
     1298  ExRootTreeReader *treeReaderNeutralHadron = new ExRootTreeReader(chainNeutralHadron);
     1299
    7701300  TChain *chainJet = new TChain("Delphes");
    7711301  chainJet->Add(inputFileJet);
    … …  
    7761306  ExRootTreeReader *treeReaderBJet = new ExRootTreeReader(chainBJet);
    7771307
     1308  TChain *chainCJet = new TChain("Delphes");
     1309  chainCJet->Add(inputFileCJet);
     1310  ExRootTreeReader *treeReaderCJet = new ExRootTreeReader(chainCJet);
     1311
    7781312  TChain *chainTauJet = new TChain("Delphes");
    7791313  chainTauJet->Add(inputFileTauJet);
    … …  
    7821316  TClonesArray *branchParticleElectron = treeReaderElectron->UseBranch("Particle");
    7831317  TClonesArray *branchTrackElectron = treeReaderElectron->UseBranch("Track");
    784   TClonesArray *branchTowerElectron = treeReaderElectron->UseBranch("Tower");
    7851318  TClonesArray *branchElectron = treeReaderElectron->UseBranch("Electron");
    7861319
    … …  
    7891322  TClonesArray *branchMuon = treeReaderMuon->UseBranch("Muon");
    7901323
     1324  TClonesArray *branchParticlePion = treeReaderPion->UseBranch("Particle");
     1325  TClonesArray *branchTrackPion = treeReaderPion->UseBranch("Track");
     1326  TClonesArray *branchPion = treeReaderPion->UseBranch("Pion");
     1327
    7911328  TClonesArray *branchParticlePhoton = treeReaderPhoton->UseBranch("Particle");
    7921329  TClonesArray *branchTowerPhoton = treeReaderPhoton->UseBranch("Tower");
    7931330  TClonesArray *branchPhoton = treeReaderPhoton->UseBranch("Photon");
    7941331
     1332  TClonesArray *branchParticleNeutralHadron = treeReaderNeutralHadron->UseBranch("Particle");
     1333  TClonesArray *branchTowerNeutralHadron = treeReaderNeutralHadron->UseBranch("Tower");
     1334
    7951335  TClonesArray *branchGenJet = treeReaderJet->UseBranch("GenJet");
     1336  TClonesArray *branchParticleJet = treeReaderJet->UseBranch("Particle");
    7961337  TClonesArray *branchPFJet = treeReaderJet->UseBranch("Jet");
    7971338  TClonesArray *branchCaloJet = treeReaderJet->UseBranch("CaloJet");
    … …  
    8001341  TClonesArray *branchPFBJet = treeReaderBJet->UseBranch("Jet");
    8011342
     1343  TClonesArray *branchParticleCJet = treeReaderCJet->UseBranch("Particle");
     1344  TClonesArray *branchPFCJet = treeReaderCJet->UseBranch("Jet");
     1345
    8021346  TClonesArray *branchParticleTauJet = treeReaderTauJet->UseBranch("Particle");
    8031347  TClonesArray *branchPFTauJet = treeReaderTauJet->UseBranch("Jet");
    8041348
    805   TClonesArray *branchScalarHT = treeReaderJet->UseBranch("ScalarHT");
     1349  TClonesArray *branchGenScalarHT = treeReaderJet->UseBranch("GenScalarHT");
    8061350  TClonesArray *branchMet = treeReaderJet->UseBranch("MissingET");
    807 
    808   ///////////////
    809   // Electrons //
    810   ///////////////
    811 
    812   // Reconstruction efficiency
    813   TString elecs = "Electron";
    814   int elID = 11;
    815   std::pair<TH1D*,TH1D*> histos_el = GetEff<Electron>(branchElectron, branchParticleElectron, "Electron", elID, treeReaderElectron);
    816 
    817   // tracking reconstruction efficiency
    818   std::pair <TH1D*,TH1D*> histos_eltrack = GetEff<Track>(branchTrackElectron, branchParticleElectron, "electronTrack", elID, treeReaderElectron);
    819 
    820   // Tower reconstruction efficiency
    821   std::pair <TH1D*,TH1D*> histos_eltower = GetEff<Tower>(branchTowerElectron, branchParticleElectron, "electronTower", elID, treeReaderElectron);
    822 
    823   // Electron Energy Resolution
    824   std::vector<resolPlot> plots_el;
    825   HistogramsCollection(&plots_el, TMath::Log10(ptrangemin), TMath::Log10(ptrangemax), "electrons");
    826   GetEres<Electron>(&plots_el, branchElectron, branchParticleElectron, elID, treeReaderElectron);
    827   TGraphErrors gr_el = EresGraph(&plots_el, true);
    828   TGraphErrors gr_elFwd = EresGraph(&plots_el, false);
    829   gr_el.SetName("Electron");
    830   gr_elFwd.SetName("ElectronFwd");
    831 
    832   // Electron Track Energy Resolution
    833   std::vector<resolPlot> plots_eltrack;
    834   HistogramsCollection(&plots_eltrack, TMath::Log10(ptrangemin), TMath::Log10(ptrangemax), "electronsTracks");
    835   GetEres<Track>(&plots_eltrack, branchTrackElectron, branchParticleElectron, elID, treeReaderElectron);
    836   TGraphErrors gr_eltrack = EresGraph(&plots_eltrack, true);
    837   TGraphErrors gr_eltrackFwd = EresGraph(&plots_eltrack, false);
    838   gr_eltrack.SetName("ElectronTracks");
    839   gr_eltrackFwd.SetName("ElectronTracksFwd");
    840 
    841   // Electron Tower Energy Resolution
    842   std::vector<resolPlot> plots_eltower;
    843   HistogramsCollection(&plots_eltower, TMath::Log10(ptrangemin), TMath::Log10(ptrangemax), "electronsTower");
    844   GetEres<Tower>(&plots_eltower, branchTowerElectron, branchParticleElectron, elID, treeReaderElectron);
    845   TGraphErrors gr_eltower = EresGraph(&plots_eltower, true);
    846   TGraphErrors gr_eltowerFwd = EresGraph(&plots_eltower, false);
    847   gr_eltower.SetName("ElectronTower");
    848   gr_eltrackFwd.SetName("ElectronTracksFwd");
    849 
    850   // Canvases
    851   TString elEff = "electronEff";
    852   TCanvas *C_el1 = new TCanvas(elEff,elEff, 1600, 600);
    853   C_el1->Divide(2);
    854   C_el1->cd(1);
    855   TLegend *leg_el1 = new TLegend(effLegXmin,effLegYmin,effLegXmax,effLegYmax);
    856   leg_el1->SetHeader("#splitline{electrons}{|#eta| < 1.5}");
    857   leg_el1->AddEntry("","","");
    858 
    859   gPad->SetLogx();
    860   histos_eltrack.first->Draw("][");
    861   addHist(histos_eltrack.first, leg_el1, 2);
    862   histos_el.first->Draw("same ][");
    863   addHist(histos_el.first, leg_el1, 1);
    864   DrawAxis(histos_eltrack.first, leg_el1,1);
    865 
    866   leg_el1->Draw();
    867 
    868   C_el1->cd(2);
    869   TLegend *leg_el2 = new TLegend(effLegXmin,effLegYmin,effLegXmax,effLegYmax);
    870   leg_el2->SetHeader("#splitline{electrons}{1.5 < |#eta| < 2.5}");
    871   leg_el2->AddEntry("","","");
    872 
    873   gPad->SetLogx();
    874   histos_eltrack.second->Draw("][");
    875   addHist(histos_eltrack.second, leg_el2, 2);
    876   histos_el.second->Draw("same ][");
    877   addHist(histos_el.second, leg_el2, 1);
    878 
    879   DrawAxis(histos_eltrack.second, leg_el2, 1);
    880   leg_el2->Draw();
    881 
    882   TString elRes = "electronERes";
    883   TString elResFwd = "electronEResForward";
    884   TCanvas *C_el2 = new TCanvas(elRes,elRes, 1600, 600);
    885   C_el2->Divide(2);
    886   C_el2->cd(1);
    887   TMultiGraph *mg_el = new TMultiGraph(elRes,elRes);
    888   TLegend *leg_el = new TLegend(resLegXmin,resLegYmin,resLegXmax,resLegYmax);
    889   leg_el->SetHeader("#splitline{electrons}{|#eta| < 1.5}");
    890   leg_el->AddEntry("","","");
    891 
    892   addGraph(mg_el, &gr_eltower, leg_el, 3);
    893   addGraph(mg_el, &gr_eltrack, leg_el, 2);
    894   addGraph(mg_el, &gr_el, leg_el, 1);
    895 
    896   mg_el->Draw("ACX");
    897   leg_el->Draw();
    898 
    899   DrawAxis(mg_el, leg_el, 0.1);
    900 
    901   C_el2->cd(2);
    902   TMultiGraph *mg_elFwd = new TMultiGraph(elResFwd,elResFwd);
    903   TLegend *leg_elFwd = new TLegend(resLegXmin,resLegYmin,resLegXmax,resLegYmax);
    904   leg_elFwd->SetHeader("#splitline{electrons}{1.5 < |#eta| < 2.5}");
    905   leg_elFwd->AddEntry("","","");
    906 
    907   addGraph(mg_elFwd, &gr_eltowerFwd, leg_elFwd, 3);
    908   addGraph(mg_elFwd, &gr_eltrackFwd, leg_elFwd, 2);
    909   addGraph(mg_elFwd, &gr_elFwd, leg_elFwd, 1);
    910 
    911   mg_elFwd->Draw("ACX");
    912   leg_elFwd->Draw();
    913 
    914   DrawAxis(mg_elFwd, leg_elFwd, 0.2);
     1351  TClonesArray *branchCaloMet = treeReaderJet->UseBranch("CaloMissingET");
     1352
     1353  std::vector<Color_t> colors;
     1354
     1355  colors.push_back(kBlack);
     1356  colors.push_back(kBlue);
     1357  colors.push_back(kRed);
     1358  colors.push_back(kGreen+1);
     1359  colors.push_back(kMagenta+1);
     1360  colors.push_back(kOrange);
     1361
     1362  std::vector<Int_t> markerStyles;
     1363
     1364  markerStyles.push_back(20);
     1365  markerStyles.push_back(21);
     1366  markerStyles.push_back(22);
     1367  markerStyles.push_back(23);
     1368  markerStyles.push_back(33);
     1369  markerStyles.push_back(34);
    9151370
    9161371  TString pdfOutput(outputFile);
    9171372  pdfOutput.ReplaceAll(".root", ".pdf");
    9181373
    919   C_el1->Print(pdfOutput+"(","pdf");
    920   C_el2->Print(pdfOutput,"pdf");
    921 
    922   gDirectory->cd(0);
    923 
    924   ///////////
    925   // Muons //
    926   ///////////
    927 
    928   // Reconstruction efficiency
    929   int muID = 13;
    930   std::pair<TH1D*,TH1D*> histos_mu = GetEff<Muon>(branchMuon, branchParticleMuon,"Muon", muID, treeReaderMuon);
    931 
    932   // muon tracking reconstruction efficiency
    933   std::pair <TH1D*,TH1D*> histos_mutrack = GetEff<Track>(branchTrackMuon, branchParticleMuon, "muonTrack", muID, treeReaderMuon);
    934 
    935   // Muon Pt Resolution
    936   std::vector<resolPlot> plots_mu;
    937   HistogramsCollection(&plots_mu, TMath::Log10(ptrangemin), TMath::Log10(ptrangemax), "muons");
    938   GetPtres<Muon>(&plots_mu, branchMuon, branchParticleMuon, muID, treeReaderMuon);
    939   TGraphErrors gr_mu = EresGraph(&plots_mu, true);
    940   TGraphErrors gr_muFwd = EresGraph(&plots_mu, false);
    941   gr_mu.SetName("Muon");
    942   gr_muFwd.SetName("MuonFwd");
    943 
    944   // Muon Track Energy Resolution
    945   std::vector<resolPlot> plots_mutrack;
    946   HistogramsCollection(&plots_mutrack, TMath::Log10(ptrangemin), TMath::Log10(ptrangemax), "muonsTracks");
    947   GetPtres<Track>(&plots_mutrack, branchTrackMuon, branchParticleMuon, muID, treeReaderMuon);
    948   TGraphErrors gr_mutrack = EresGraph(&plots_mutrack, true);
    949   TGraphErrors gr_mutrackFwd = EresGraph(&plots_mutrack, false);
    950   gr_mutrackFwd.SetName("MuonTracksFwd");
    951 
    952   // Canvas
    953 
    954   TString muEff = "muonEff";
    955   TCanvas *C_mu1 = new TCanvas(muEff,muEff, 1600, 600);
    956   C_mu1->Divide(2);
    957   C_mu1->cd(1);
    958   TLegend *leg_mu1 = new TLegend(effLegXmin,effLegYmin,effLegXmax,effLegYmax);
    959   leg_mu1->SetHeader("#splitline{muons}{|#eta| < 1.5}");
    960   leg_mu1->AddEntry("","","");
    961 
    962 
    963   gPad->SetLogx();
    964   histos_mutrack.first->Draw("][");
    965   addHist(histos_mutrack.first, leg_mu1, 2);
    966   histos_mu.first->Draw("same ][");
    967   addHist(histos_mu.first, leg_mu1, 1);
    968 
    969   DrawAxis(histos_mutrack.first, leg_mu1, 1);
    970 
    971   leg_mu1->Draw();
    972 
    973   C_mu1->cd(2);
    974   TLegend *leg_mu2 = new TLegend(effLegXmin,effLegYmin,effLegXmax,effLegYmax);
    975   leg_mu2->SetHeader("#splitline{muons}{1.5 < |#eta| < 2.5}");
    976   leg_mu2->AddEntry("","","");
    977 
    978   gPad->SetLogx();
    979   histos_mutrack.second->Draw("][");
    980   addHist(histos_mutrack.second, leg_mu2, 2);
    981   histos_mu.second->Draw("same ][");
    982   addHist(histos_mu.second, leg_mu2, 1);
    983 
    984   DrawAxis(histos_mutrack.second, leg_mu2, 1);
    985   leg_mu2->Draw();
    986 
    987   TString muRes = "muonERes";
    988   TString muResFwd = "muonEResFwd";
    989 
    990   TCanvas *C_mu = new TCanvas(muRes,muRes, 1600, 600);
    991   C_mu->Divide(2);
    992   C_mu->cd(1);
    993   TMultiGraph *mg_mu = new TMultiGraph(muRes,muRes);
    994   TLegend *leg_mu = new TLegend(topLeftLegXmin,topLeftLegYmin,topLeftLegXmax,topLeftLegYmax);
    995   leg_mu->SetHeader("#splitline{muons}{|#eta| < 1.5}");
    996   leg_mu->AddEntry("","","");
    997 
    998   addGraph(mg_mu, &gr_mutrack, leg_mu, 2);
    999   addGraph(mg_mu, &gr_mu, leg_mu, 1);
    1000 
    1001   mg_mu->Draw("ACX");
    1002   leg_mu->Draw();
    1003 
    1004   DrawAxis(mg_mu, leg_mu, 0.3, 1);
    1005 
    1006   C_mu->cd(2);
    1007   TMultiGraph *mg_muFwd = new TMultiGraph(muResFwd,muResFwd);
    1008   TLegend *leg_muFwd = new TLegend(topLeftLegXmin,topLeftLegYmin,topLeftLegXmax,topLeftLegYmax);
    1009   leg_muFwd->SetHeader("#splitline{muons}{1.5 < |#eta| < 2.5}");
    1010   leg_muFwd->AddEntry("","","");
    1011 
    1012   addGraph(mg_muFwd, &gr_mutrackFwd, leg_muFwd, 2);
    1013   addGraph(mg_muFwd, &gr_muFwd, leg_muFwd, 1);
    1014 
    1015   mg_muFwd->Draw("ACX");
    1016   leg_muFwd->Draw();
    1017 
    1018   DrawAxis(mg_muFwd, leg_muFwd, 0.3, 1);
    1019 
    1020   C_mu1->Print(pdfOutput,"pdf");
    1021   C_mu->Print(pdfOutput,"pdf");
    1022 
    1023   gDirectory->cd(0);
    1024 
    1025   /////////////
    1026   // Photons //
    1027   /////////////
    1028 
    1029   // Reconstruction efficiency
    1030   int phID = 22;
    1031   std::pair<TH1D*,TH1D*> histos_ph = GetEff<Electron>(branchPhoton, branchParticlePhoton, "Photon", phID, treeReaderPhoton);
    1032   std::pair<TH1D*,TH1D*> histos_phtower = GetEff<Electron>(branchTowerPhoton, branchParticlePhoton, "Photon", phID, treeReaderPhoton);
    1033 
    1034   // Photon Energy Resolution
    1035   std::vector<resolPlot> plots_ph;
    1036   HistogramsCollection(&plots_ph, TMath::Log10(ptrangemin), TMath::Log10(ptrangemax), "photons");
    1037   GetEres<Photon>(&plots_ph, branchPhoton, branchParticlePhoton, phID, treeReaderPhoton);
    1038   TGraphErrors gr_ph = EresGraph(&plots_ph, true);
    1039   TGraphErrors gr_phFwd = EresGraph(&plots_ph, false);
    1040   gr_ph.SetName("Photon");
    1041   gr_phFwd.SetName("PhotonFwd");
    1042 
    1043 
    1044   // Photon Tower Energy Resolution
    1045   std::vector<resolPlot> plots_phtower;
    1046   HistogramsCollection(&plots_phtower, TMath::Log10(ptrangemin), TMath::Log10(ptrangemax), "photonsTower");
    1047   GetEres<Tower>(&plots_phtower, branchTowerPhoton, branchParticlePhoton, phID, treeReaderPhoton);
    1048   TGraphErrors gr_phtower = EresGraph(&plots_phtower, true);
    1049   TGraphErrors gr_phtowerFwd = EresGraph(&plots_phtower, false);
    1050   gr_phtower.SetName("PhotonTower");
    1051   gr_phtowerFwd.SetName("PhotonTowerFwd");
    1052 
    1053   // Canvas
    1054 
    1055   TString phEff = "photonEff";
    1056   TCanvas *C_ph1 = new TCanvas(phEff,phEff, 1600, 600);
    1057   C_ph1->Divide(2);
    1058   C_ph1->cd(1);
    1059   TLegend *leg_ph1 = new TLegend(effLegXmin,effLegYmin,effLegXmax,effLegYmax);
    1060   leg_ph1->SetHeader("#splitline{photons}{|#eta| < 1.5}");
    1061   leg_ph1->AddEntry("","","");
    1062 
    1063 
    1064   gPad->SetLogx();
    1065   histos_phtower.first->Draw("][");
    1066   addHist(histos_phtower.first, leg_ph1, 3);
    1067   histos_ph.first->Draw("same ][");
    1068   addHist(histos_ph.first, leg_ph1, 1);
    1069 
    1070   DrawAxis(histos_phtower.first, leg_ph1, 1);
    1071   leg_ph1->Draw();
    1072 
    1073   C_ph1->cd(2);
    1074   TLegend *leg_ph2 = new TLegend(effLegXmin,effLegYmin,effLegXmax,effLegYmax);
    1075   leg_ph2->SetHeader("#splitline{photons}{1.5 < |#eta| < 2.5}");
    1076   leg_ph2->AddEntry("","","");
    1077 
    1078 
    1079   gPad->SetLogx();
    1080   histos_phtower.second->Draw("][");
    1081   addHist(histos_phtower.second, leg_ph2, 3);
    1082   histos_ph.second->Draw("same ][");
    1083   addHist(histos_ph.second, leg_ph2, 1);
    1084 
    1085   DrawAxis(histos_phtower.second, leg_ph2, 1);
    1086   leg_ph2->Draw();
    1087 
    1088   TString phRes = "phERes";
    1089   TString phResFwd = "phEResFwd";
    1090 
    1091   TCanvas *C_ph = new TCanvas(phRes,phRes, 1600, 600);
    1092   C_ph->Divide(2);
    1093   C_ph->cd(1);
    1094   TMultiGraph *mg_ph = new TMultiGraph(phRes,phRes);
    1095   TLegend *leg_ph = new TLegend(resLegXmin,resLegYmin,resLegXmax,resLegYmax);
    1096   leg_ph->SetHeader("#splitline{photons}{|#eta| < 1.5}");
    1097   leg_ph->AddEntry("","","");
    1098 
    1099   addGraph(mg_ph, &gr_phtower, leg_ph, 3);
    1100   addGraph(mg_ph, &gr_ph, leg_ph, 1);
    1101 
    1102   mg_ph->Draw("ACX");
    1103   leg_ph->Draw();
    1104 
    1105   DrawAxis(mg_ph, leg_ph, 0.1);
    1106 
    1107   C_ph->cd(2);
    1108   TMultiGraph *mg_phFwd = new TMultiGraph(phResFwd,phResFwd);
    1109   TLegend *leg_phFwd = new TLegend(resLegXmin,resLegYmin,resLegXmax,resLegYmax);
    1110   leg_phFwd->SetHeader("#splitline{photons}{1.5 < |#eta| < 2.5}");
    1111   leg_phFwd->AddEntry("","","");
    1112 
    1113   addGraph(mg_phFwd, &gr_phtowerFwd, leg_phFwd, 3);
    1114   addGraph(mg_phFwd, &gr_phFwd, leg_phFwd, 1);
    1115 
    1116   mg_phFwd->Draw("ACX");
    1117   leg_phFwd->Draw();
    1118 
    1119   DrawAxis(mg_phFwd, leg_phFwd, 0.1);
    1120 
    1121   C_ph1->Print(pdfOutput,"pdf");
    1122   C_ph->Print(pdfOutput,"pdf");
    1123 
    1124   gDirectory->cd(0);
    1125 
    1126   //////////
    1127   // Jets //
    1128   //////////
    1129 
    1130   // BJets Reconstruction efficiency
    1131   int bID = 5;
    1132   std::pair<TH1D*,TH1D*> histos_btag = GetEff<Jet>(branchPFBJet, branchParticleBJet,"BTag", bID, treeReaderBJet);
    1133 
    1134   // TauJets Reconstruction efficiency
    1135   int tauID = 15;
    1136   std::pair<TH1D*,TH1D*> histos_tautag = GetEff<Jet>(branchPFTauJet, branchParticleTauJet,"TauTag", tauID, treeReaderTauJet);
    1137 
    1138   // PFJets Energy Resolution
    1139   std::vector<resolPlot> plots_pfjets;
    1140   HistogramsCollection(&plots_pfjets, TMath::Log10(ptrangemin), TMath::Log10(ptrangemax), "PFJet");
    1141   GetJetsEres(&plots_pfjets, branchPFJet, branchGenJet, treeReaderJet);
    1142   TGraphErrors gr_pfjets = EresGraph(&plots_pfjets, true);
    1143   TGraphErrors gr_pfjetsFwd = EresGraph(&plots_pfjets, false);
    1144   gr_pfjets.SetName("pfJet");
    1145   gr_pfjetsFwd.SetName("pfJetFwd");
    1146 
    1147   // CaloJets Energy Resolution
    1148   std::vector<resolPlot> plots_calojets;
    1149   HistogramsCollection(&plots_calojets, TMath::Log10(ptrangemin), TMath::Log10(ptrangemax), "CaloJet");
    1150   GetJetsEres(&plots_calojets, branchCaloJet, branchGenJet, treeReaderJet);
    1151   TGraphErrors gr_calojets = EresGraph(&plots_calojets, true);
    1152   TGraphErrors gr_calojetsFwd = EresGraph(&plots_calojets, false);
    1153   gr_calojets.SetName("caloJet");
    1154   gr_calojetsFwd.SetName("caloJetFwd");
    1155 
    1156   // MET Resolution vs HT
    1157   std::vector<resolPlot> plots_met;
    1158   HistogramsCollection(&plots_met, TMath::Log10(ptrangemin), TMath::Log10(ptrangemax), "MET", -500, 500);
    1159   GetMetres(&plots_met, branchScalarHT, branchMet, branchPFJet, treeReaderJet);
    1160   TGraphErrors gr_met = EresGraph(&plots_met, true);
    1161   gr_calojets.SetName("MET");
    1162 
    1163   // Canvas
    1164   TString btagEff = "btagEff";
    1165   TCanvas *C_btag1 = new TCanvas(btagEff,btagEff, 1600, 600);
    1166   C_btag1->Divide(2);
    1167   C_btag1->cd(1);
    1168   TLegend *leg_btag1 = new TLegend(resLegXmin,resLegYmin,resLegXmax,resLegYmax);
    1169   leg_btag1->SetHeader("#splitline{B-tagging}{|#eta| < 1.5}");
    1170   leg_btag1->AddEntry("","","");
    1171 
    1172   gPad->SetLogx();
    1173   histos_btag.first->Draw();
    1174   addHist(histos_btag.first, leg_btag1, 0);
    1175 
    1176   DrawAxis(histos_btag.first, leg_btag1, 1);
    1177   leg_btag1->Draw();
    1178 
    1179   C_btag1->cd(2);
    1180   TLegend *leg_btag2 = new TLegend(resLegXmin,resLegYmin,resLegXmax+0.05,resLegYmax);
    1181   leg_btag2->SetHeader("#splitline{B-tagging}{1.5 < |#eta| < 2.5}");
    1182   leg_btag2->AddEntry("","","");
    1183 
    1184   gPad->SetLogx();
    1185   histos_btag.second->Draw();
    1186   addHist(histos_btag.second, leg_btag2, 0);
    1187 
    1188   DrawAxis(histos_btag.second, leg_btag2, 1);
    1189   leg_btag2->Draw();
    1190   C_btag1->cd(0);
    1191 
    1192   TString tautagEff = "tautagEff";
    1193   TCanvas *C_tautag1 = new TCanvas(tautagEff,tautagEff, 1600, 600);
    1194   C_tautag1->Divide(2);
    1195   C_tautag1->cd(1);
    1196   TLegend *leg_tautag1 = new TLegend(resLegXmin,resLegYmin,resLegXmax,resLegYmax);
    1197   leg_tautag1->SetHeader("#splitline{#tau-tagging}{|#eta| < 1.5}");
    1198   leg_tautag1->AddEntry("","","");
    1199 
    1200   gPad->SetLogx();
    1201   histos_tautag.first->Draw();
    1202   addHist(histos_tautag.first, leg_tautag1, 0);
    1203 
    1204   DrawAxis(histos_tautag.first, leg_tautag1, 1);
    1205   leg_tautag1->Draw();
    1206 
    1207   C_tautag1->cd(2);
    1208   TLegend *leg_tautag2 = new TLegend(resLegXmin,resLegYmin,resLegXmax+0.05,resLegYmax);
    1209   leg_tautag2->SetHeader("#splitline{#tau-tagging}{1.5 < |#eta| < 2.5}");
    1210   leg_tautag2->AddEntry("","","");
    1211 
    1212   gPad->SetLogx();
    1213   histos_tautag.second->Draw();
    1214   addHist(histos_tautag.second, leg_tautag2, 0);
    1215 
    1216   DrawAxis(histos_tautag.second, leg_tautag2, 1);
    1217   leg_tautag2->Draw();
    1218   C_tautag1->cd(0);
    1219 
    1220   TString jetRes = "jetERes";
    1221   TString jetResFwd = "jetEResFwd";
    1222   TCanvas *C_jet = new TCanvas(jetRes,jetRes, 1600, 600);
    1223   C_jet->Divide(2);
    1224 
    1225   C_jet->cd(1);
    1226   TMultiGraph *mg_jet = new TMultiGraph(jetRes,jetRes);
    1227   TLegend *leg_jet = new TLegend(resLegXmin,resLegYmin,resLegXmax,resLegYmax);
    1228   leg_jet->SetHeader("#splitline{jets}{|#eta| < 1.5}");
    1229   leg_jet->AddEntry("","","");
    1230 
    1231   addGraph(mg_jet, &gr_calojets, leg_jet, 3);
    1232   addGraph(mg_jet, &gr_pfjets, leg_jet, 1);
    1233 
    1234   mg_jet->Draw("ALX");
    1235   leg_jet->Draw();
    1236 
    1237   DrawAxis(mg_jet, leg_jet, 0.5);
    1238 
    1239   C_jet->cd(2);
    1240   TMultiGraph *mg_jetFwd = new TMultiGraph(jetResFwd,jetResFwd);
    1241   TLegend *leg_jetFwd = new TLegend(resLegXmin,resLegYmin,resLegXmax,resLegYmax);
    1242   leg_jetFwd->SetHeader("#splitline{jets}{1.5 < |#eta| < 2.5}");
    1243   leg_jetFwd->AddEntry("","","");
    1244 
    1245   addGraph(mg_jetFwd, &gr_calojetsFwd, leg_jetFwd, 3);
    1246   addGraph(mg_jetFwd, &gr_pfjetsFwd, leg_jetFwd, 1);
    1247 
    1248   mg_jetFwd->Draw("ALX");
    1249   leg_jetFwd->Draw();
    1250 
    1251   DrawAxis(mg_jetFwd, leg_jetFwd, 0.5);
    1252 
    1253   TString metRes = "MetRes";
    1254   TCanvas *C_met = new TCanvas(metRes,metRes, 800, 600);
    1255 
    1256   TMultiGraph *mg_met = new TMultiGraph(metRes,metRes);
    1257   TLegend *leg_met = new TLegend(topLeftLegXmin,topLeftLegYmin+0.2,topLeftLegXmax-0.2,topLeftLegYmax);
    1258   leg_met->SetHeader("E_{T}^{miss}");
    1259   leg_met->AddEntry("","","");
    1260 
    1261 
    1262   addGraph(mg_met, &gr_met, leg_met, 0);
    1263 
    1264   mg_met->Draw("ACX");
    1265   leg_met->Draw();
    1266 
    1267   DrawAxis(mg_met, leg_met, 300);
    1268 
    1269   mg_met->GetXaxis()->SetTitle("H_{T} [GeV]");
    1270   mg_met->GetYaxis()->SetTitle("#sigma(ME_{x}) [GeV]");
    1271 
    1272   C_jet->Print(pdfOutput,"pdf");
    1273   C_btag1->Print(pdfOutput,"pdf");
    1274   C_tautag1->Print(pdfOutput,"pdf");
    1275   C_met->Print(pdfOutput+")","pdf");
     1374  TString figPath = inputFilePion;
     1375  figPath.ReplaceAll(".root", "");
     1376  figPath.ReplaceAll("root", "www/fig");
     1377  Int_t lastSlash = figPath.Last('/');
     1378  Int_t sizePath = figPath.Length();
     1379  figPath.Remove(lastSlash+1,sizePath);
     1380
     1381  TString s_etaMin, s_etaMax, s_eta, s_pt, s_e;
     1382
     1383  Double_t ptMin = 1.;
     1384  Double_t ptMax = 50000.;
     1385
     1386  Double_t etaMin = -6.;
     1387  Double_t etaMax = 6.;
     1388
     1389  std::vector<Double_t> ptVals;
     1390
     1391  ptVals.push_back(1.);
     1392  ptVals.push_back(10.);
     1393  ptVals.push_back(100.);
     1394  ptVals.push_back(1000.);
     1395  ptVals.push_back(10000.);
     1396
     1397  std::vector<Double_t> etaVals;
     1398
     1399  etaVals.push_back(0.);
     1400  etaVals.push_back(1.5);
     1401  etaVals.push_back(2.5);
     1402  etaVals.push_back(4.0);
     1403  etaVals.push_back(6.0);
     1404
     1405  const int n_etabins = etaVals.size()-1;
     1406  const int n_ptbins = ptVals.size();
     1407
     1408
     1409  //////////////////////////
     1410  // Tracking performance //
     1411  //////////////////////////
     1412
     1413  // --------- Pion Tracks  --------- //
     1414
     1415  TMultiGraph *mg_trkpi_res_pt  = new TMultiGraph("","");
     1416  TMultiGraph *mg_trkpi_eff_pt  = new TMultiGraph("","");
     1417  TMultiGraph *mg_trkpi_res_eta = new TMultiGraph("","");
     1418  TMultiGraph *mg_trkpi_eff_eta = new TMultiGraph("","");
     1419
     1420  TLegend *leg_trkpi_res_pt = new TLegend(0.55,0.22,0.90,0.48);
     1421  TLegend *leg_trkpi_eff_pt = (TLegend*)leg_trkpi_res_pt->Clone();
     1422  TLegend *leg_trkpi_res_eta = (TLegend*)leg_trkpi_res_pt->Clone();
     1423  TLegend *leg_trkpi_eff_eta = (TLegend*)leg_trkpi_res_eta->Clone();
     1424
     1425
     1426  TGraphErrors gr_trkpi_res_pt[n_etabins], gr_trkpi_eff_pt[n_etabins], gr_trkpi_res_eta[n_ptbins], gr_trkpi_eff_eta[n_ptbins];
     1427  TH1D* h_trkpi_eff_pt, *h_trkpi_eff_eta;
     1428
     1429  std::vector<resolPlot> plots_trkpi_res_pt[n_etabins], plots_trkpi_res_eta[n_ptbins];
     1430
     1431  // loop over eta bins
     1432  for (k = 0; k < etaVals.size()-1; k++)
     1433  {
     1434     HistogramsCollection(&plots_trkpi_res_pt[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "trkpi");
     1435     GetPtres<Track>(&plots_trkpi_res_pt[k], branchTrackPion, branchParticlePion, 211, etaVals.at(k), etaVals.at(k+1), treeReaderPion);
     1436     gr_trkpi_res_pt[k] = EresGraph(&plots_trkpi_res_pt[k]);
     1437
     1438     h_trkpi_eff_pt = GetEffPt<Track>(branchTrackPion, branchParticlePion, "Pion", 211, ptMin, ptMax, etaVals.at(k), etaVals.at(k+1), treeReaderPion);
     1439     gr_trkpi_eff_pt[k] = TGraphErrors(h_trkpi_eff_pt);
     1440
     1441     s_etaMin = Form("%.1f",etaVals.at(k));
     1442     s_etaMax = Form("%.1f",etaVals.at(k+1));
     1443
     1444     s_eta = "#pi^{ #pm} , " + s_etaMin + " < | #eta | < "+s_etaMax;
     1445
     1446     gr_trkpi_res_pt[k].SetName("trkRes_"+s_etaMin+"_"+s_etaMax);
     1447     gr_trkpi_eff_pt[k].SetName("trkEff_"+s_etaMin+"_"+s_etaMax);
     1448
     1449     addResoGraph(mg_trkpi_res_pt, &gr_trkpi_res_pt[k], leg_trkpi_res_pt, markerStyles.at(k), colors.at(k), s_eta);
     1450     addResoGraph(mg_trkpi_eff_pt, &gr_trkpi_eff_pt[k], leg_trkpi_eff_pt, markerStyles.at(k), colors.at(k), s_eta);
     1451  }
     1452
     1453  // loop over pt
     1454  for (k = 0; k < ptVals.size(); k++)
     1455  {
     1456     HistogramsCollectionVsEta(&plots_trkpi_res_eta[k], etaMin, etaMax, 0.5*ptVals.at(k), 2.0*ptVals.at(k), "trkpi", 0.0, 2.0);
     1457     GetPtresVsEta<Track>(&plots_trkpi_res_eta[k], branchTrackPion, branchParticlePion, 211, 0.5*ptVals.at(k), 2.0*ptVals.at(k), treeReaderPion);
     1458     gr_trkpi_res_eta[k] = EresGraphVsEta(&plots_trkpi_res_eta[k]);
     1459
     1460     h_trkpi_eff_eta = GetEffEta<Track>(branchTrackPion, branchParticlePion, "Pion", 211, 0.5*ptVals.at(k), 2.0*ptVals.at(k) ,etaMin, etaMax , treeReaderPion);
     1461     gr_trkpi_eff_eta[k] = TGraphErrors(h_trkpi_eff_eta);
     1462
     1463     s_pt = Form("#pi^{ #pm} , p_{T} = %.0f GeV",ptVals.at(k));
     1464     if(ptVals.at(k) >= 1000.) s_pt = Form("#pi^{ #pm} , p_{T} = %.0f TeV",ptVals.at(k)/1000.);
     1465
     1466     addResoGraph(mg_trkpi_res_eta, &gr_trkpi_res_eta[k], leg_trkpi_res_eta, markerStyles.at(k), colors.at(k), s_pt );
     1467     addResoGraph(mg_trkpi_eff_eta, &gr_trkpi_eff_eta[k], leg_trkpi_eff_eta, markerStyles.at(k), colors.at(k), s_pt );
     1468
     1469  }
     1470
     1471  TCanvas *c_trkpi_res_pt = new TCanvas("","", 800, 600);
     1472
     1473  mg_trkpi_res_pt->Draw("APE");
     1474  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);
     1475  leg_trkpi_res_pt->Draw();
     1476
     1477  c_trkpi_res_pt->Print(pdfOutput+"(","pdf");
     1478  c_trkpi_res_pt->Print(figPath+"img_trkpi_res_pt.pdf","pdf");
     1479  c_trkpi_res_pt->Print(figPath+"img_trkpi_res_pt.png","png");
     1480
     1481  TCanvas *c_trkpi_res_eta = new TCanvas("","", 800, 600);
     1482
     1483  mg_trkpi_res_eta->Draw("APE");
     1484  DrawAxis(mg_trkpi_res_eta, leg_trkpi_res_eta, etaMin, etaMax, 0.01, 100, " #eta ", "(track resolution in p_{T})/p_{T} (%)", false, true);
     1485  leg_trkpi_res_eta->Draw();
     1486
     1487  c_trkpi_res_eta->Print(pdfOutput,"pdf");
     1488  c_trkpi_res_eta->Print(figPath+"img_trkpi_res_eta.pdf","pdf");
     1489  c_trkpi_res_eta->Print(figPath+"img_trkpi_res_eta.png","png");
     1490
     1491  TCanvas *c_trkpi_eff_pt = new TCanvas("","", 800, 600);
     1492
     1493  mg_trkpi_eff_pt->Draw("APE");
     1494  DrawAxis(mg_trkpi_eff_pt, leg_trkpi_eff_pt, ptMin, ptMax, 0.0, 100, "p_{T} [GeV]", "tracking efficiency (%)", true, false);
     1495  leg_trkpi_eff_pt->Draw();
     1496
     1497  c_trkpi_eff_pt->Print(pdfOutput,"pdf");
     1498  c_trkpi_eff_pt->Print(figPath+"img_trkpi_eff_pt.pdf","pdf");
     1499  c_trkpi_eff_pt->Print(figPath+"img_trkpi_eff_pt.png","png");
     1500
     1501  TCanvas *c_trkpi_eff_eta = new TCanvas("","", 800, 600);
     1502
     1503  mg_trkpi_eff_eta->Draw("APE");
     1504  DrawAxis(mg_trkpi_eff_eta, leg_trkpi_eff_eta, etaMin, etaMax, 0.0, 100, " #eta ", "tracking efficiency (%)", false, false);
     1505  leg_trkpi_eff_eta->Draw();
     1506
     1507  c_trkpi_eff_eta->Print(pdfOutput,"pdf");
     1508  c_trkpi_eff_eta->Print(figPath+"img_trkpi_eff_eta.pdf","pdf");
     1509  c_trkpi_eff_eta->Print(figPath+"img_trkpi_eff_eta.png","png");
     1510
     1511  // --------- Electron Tracks  --------- //
     1512
     1513  TMultiGraph *mg_trkele_res_pt  = new TMultiGraph("","");
     1514  TMultiGraph *mg_trkele_eff_pt  = new TMultiGraph("","");
     1515  TMultiGraph *mg_trkele_res_eta = new TMultiGraph("","");
     1516  TMultiGraph *mg_trkele_eff_eta = new TMultiGraph("","");
     1517
     1518  TLegend *leg_trkele_res_pt = new TLegend(0.55,0.22,0.90,0.48);
     1519  TLegend *leg_trkele_eff_pt = (TLegend*)leg_trkele_res_pt->Clone();
     1520  TLegend *leg_trkele_res_eta = (TLegend*)leg_trkele_res_pt->Clone();
     1521  TLegend *leg_trkele_eff_eta = (TLegend*)leg_trkele_res_eta->Clone();
     1522
     1523  TGraphErrors gr_trkele_res_pt[n_etabins], gr_trkele_eff_pt[n_etabins], gr_trkele_res_eta[n_ptbins], gr_trkele_eff_eta[n_ptbins];
     1524  TH1D* h_trkele_eff_pt, *h_trkele_eff_eta;
     1525
     1526  std::vector<resolPlot> plots_trkele_res_pt[n_etabins], plots_trkele_res_eta[n_ptbins];
     1527
     1528  // loop over eta bins
     1529  for (k = 0; k < etaVals.size()-1; k++)
     1530  {
     1531     HistogramsCollection(&plots_trkele_res_pt[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "trkele");
     1532     GetPtres<Track>(&plots_trkele_res_pt[k], branchTrackElectron, branchParticleElectron, 11, etaVals.at(k), etaVals.at(k+1), treeReaderElectron);
     1533     gr_trkele_res_pt[k] = EresGraph(&plots_trkele_res_pt[k]);
     1534
     1535     h_trkele_eff_pt = GetEffPt<Track>(branchTrackElectron, branchParticleElectron, "Electron", 11, ptMin, ptMax, etaVals.at(k), etaVals.at(k+1), treeReaderElectron);
     1536     gr_trkele_eff_pt[k] = TGraphErrors(h_trkele_eff_pt);
     1537
     1538     s_etaMin = Form("%.1f",etaVals.at(k));
     1539     s_etaMax = Form("%.1f",etaVals.at(k+1));
     1540
     1541     s_eta = "e^{ #pm} , " + s_etaMin + " < | #eta | < "+s_etaMax;
     1542
     1543     gr_trkele_res_pt[k].SetName("trkRes_"+s_etaMin+"_"+s_etaMax);
     1544     gr_trkele_eff_pt[k].SetName("trkEff_"+s_etaMin+"_"+s_etaMax);
     1545
     1546     addResoGraph(mg_trkele_res_pt, &gr_trkele_res_pt[k], leg_trkele_res_pt, markerStyles.at(k), colors.at(k), s_eta);
     1547     addResoGraph(mg_trkele_eff_pt, &gr_trkele_eff_pt[k], leg_trkele_eff_pt, markerStyles.at(k), colors.at(k), s_eta);
     1548  }
     1549
     1550  // loop over pt
     1551  for (k = 0; k < ptVals.size(); k++)
     1552  {
     1553     HistogramsCollectionVsEta(&plots_trkele_res_eta[k], etaMin, etaMax, 0.5*ptVals.at(k), 2.0*ptVals.at(k), "trkele", 0.0, 2.0);
     1554     GetPtresVsEta<Track>(&plots_trkele_res_eta[k], branchTrackElectron, branchParticleElectron, 11, 0.5*ptVals.at(k), 2.0*ptVals.at(k), treeReaderElectron);
     1555     gr_trkele_res_eta[k] = EresGraphVsEta(&plots_trkele_res_eta[k]);
     1556
     1557     h_trkele_eff_eta = GetEffEta<Track>(branchTrackElectron, branchParticleElectron, "Electron", 11, 0.5*ptVals.at(k), 2.0*ptVals.at(k) ,etaMin, etaMax , treeReaderElectron);
     1558     gr_trkele_eff_eta[k] = TGraphErrors(h_trkele_eff_eta);
     1559
     1560     s_pt = Form("e^{ #pm} , p_{T} = %.0f GeV",ptVals.at(k));
     1561     if(ptVals.at(k) >= 1000.) s_pt = Form("e^{ #pm} , p_{T} = %.0f TeV",ptVals.at(k)/1000.);
     1562
     1563     addResoGraph(mg_trkele_res_eta, &gr_trkele_res_eta[k], leg_trkele_res_eta, markerStyles.at(k), colors.at(k), s_pt );
     1564     addResoGraph(mg_trkele_eff_eta, &gr_trkele_eff_eta[k], leg_trkele_eff_eta, markerStyles.at(k), colors.at(k), s_pt );
     1565
     1566  }
     1567
     1568  TCanvas *c_trkele_res_pt = new TCanvas("","", 800, 600);
     1569
     1570  mg_trkele_res_pt->Draw("APE");
     1571  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);
     1572  leg_trkele_res_pt->Draw();
     1573
     1574  c_trkele_res_pt->Print(pdfOutput,"pdf");
     1575  c_trkele_res_pt->Print(figPath+"img_trkele_res_pt.pdf","pdf");
     1576  c_trkele_res_pt->Print(figPath+"img_trkele_res_pt.png","png");
     1577
     1578  TCanvas *c_trkele_res_eta = new TCanvas("","", 800, 600);
     1579
     1580  mg_trkele_res_eta->Draw("APE");
     1581  DrawAxis(mg_trkele_res_eta, leg_trkele_res_eta, etaMin, etaMax, 0.01, 100, " #eta ", "(track resolution in p_{T})/p_{T} (%)", false, true);
     1582  leg_trkele_res_eta->Draw();
     1583
     1584  c_trkele_res_eta->Print(pdfOutput,"pdf");
     1585  c_trkele_res_eta->Print(figPath+"img_trkele_res_eta.pdf","pdf");
     1586  c_trkele_res_eta->Print(figPath+"img_trkele_res_eta.png","png");
     1587
     1588  TCanvas *c_trkele_eff_pt = new TCanvas("","", 800, 600);
     1589
     1590  mg_trkele_eff_pt->Draw("APE");
     1591  DrawAxis(mg_trkele_eff_pt, leg_trkele_eff_pt, ptMin, ptMax, 0.0, 100, "p_{T} [GeV]", "tracking efficiency (%)", true, false);
     1592  leg_trkele_eff_pt->Draw();
     1593
     1594  c_trkele_eff_pt->Print(pdfOutput,"pdf");
     1595  c_trkele_eff_pt->Print(figPath+"img_trkele_eff_pt.pdf","pdf");
     1596  c_trkele_eff_pt->Print(figPath+"img_trkele_eff_pt.png","png");
     1597
     1598  TCanvas *c_trkele_eff_eta = new TCanvas("","", 800, 600);
     1599
     1600  mg_trkele_eff_eta->Draw("APE");
     1601  DrawAxis(mg_trkele_eff_eta, leg_trkele_eff_eta, etaMin, etaMax, 0.0, 100, " #eta ", "tracking efficiency (%)", false, false);
     1602  leg_trkele_eff_eta->Draw();
     1603
     1604  c_trkele_eff_eta->Print(pdfOutput,"pdf");
     1605  c_trkele_eff_eta->Print(figPath+"img_trkele_eff_eta.pdf","pdf");
     1606  c_trkele_eff_eta->Print(figPath+"img_trkele_eff_eta.png","png");
     1607
     1608
     1609  // --------- Muon Tracks  --------- //
     1610
     1611  TMultiGraph *mg_trkmu_res_pt  = new TMultiGraph("","");
     1612  TMultiGraph *mg_trkmu_eff_pt  = new TMultiGraph("","");
     1613  TMultiGraph *mg_trkmu_res_eta = new TMultiGraph("","");
     1614  TMultiGraph *mg_trkmu_eff_eta = new TMultiGraph("","");
     1615
     1616  TLegend *leg_trkmu_res_pt = new TLegend(0.55,0.22,0.90,0.48);
     1617  TLegend *leg_trkmu_eff_pt = (TLegend*)leg_trkmu_res_pt->Clone();
     1618
     1619  TLegend *leg_trkmu_res_eta = (TLegend*)leg_trkmu_res_pt->Clone();
     1620  TLegend *leg_trkmu_eff_eta = (TLegend*)leg_trkmu_res_eta->Clone();
     1621
     1622
     1623  TGraphErrors gr_trkmu_res_pt[n_etabins], gr_trkmu_eff_pt[n_etabins], gr_trkmu_res_eta[n_ptbins], gr_trkmu_eff_eta[n_ptbins];
     1624  TH1D* h_trkmu_eff_pt, *h_trkmu_eff_eta;
     1625
     1626  std::vector<resolPlot> plots_trkmu_res_pt[n_etabins], plots_trkmu_res_eta[n_ptbins];
     1627
     1628  // loop over eta bins
     1629  for (k = 0; k < etaVals.size()-1; k++)
     1630  {
     1631     HistogramsCollection(&plots_trkmu_res_pt[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "trkmu");
     1632     GetPtres<Track>(&plots_trkmu_res_pt[k], branchTrackMuon, branchParticleMuon, 13, etaVals.at(k), etaVals.at(k+1), treeReaderMuon);
     1633     gr_trkmu_res_pt[k] = EresGraph(&plots_trkmu_res_pt[k]);
     1634
     1635     h_trkmu_eff_pt = GetEffPt<Track>(branchTrackMuon, branchParticleMuon, "Muon", 13, ptMin, ptMax, etaVals.at(k), etaVals.at(k+1), treeReaderMuon);
     1636     gr_trkmu_eff_pt[k] = TGraphErrors(h_trkmu_eff_pt);
     1637
     1638     s_etaMin = Form("%.1f",etaVals.at(k));
     1639     s_etaMax = Form("%.1f",etaVals.at(k+1));
     1640
     1641     s_eta = "#mu^{ #pm} , " + s_etaMin + " < | #eta | < "+s_etaMax;
     1642
     1643     gr_trkmu_res_pt[k].SetName("trkRes_"+s_etaMin+"_"+s_etaMax);
     1644     gr_trkmu_eff_pt[k].SetName("trkEff_"+s_etaMin+"_"+s_etaMax);
     1645
     1646     addResoGraph(mg_trkmu_res_pt, &gr_trkmu_res_pt[k], leg_trkmu_res_pt, markerStyles.at(k), colors.at(k), s_eta);
     1647     addResoGraph(mg_trkmu_eff_pt, &gr_trkmu_eff_pt[k], leg_trkmu_eff_pt, markerStyles.at(k), colors.at(k), s_eta);
     1648  }
     1649
     1650  // loop over pt
     1651  for (k = 0; k < ptVals.size(); k++)
     1652  {
     1653     HistogramsCollectionVsEta(&plots_trkmu_res_eta[k], etaMin, etaMax, 0.5*ptVals.at(k), 2.0*ptVals.at(k), "trkmu", 0.0, 2.0);
     1654     GetPtresVsEta<Track>(&plots_trkmu_res_eta[k], branchTrackMuon, branchParticleMuon, 13, 0.5*ptVals.at(k), 2.0*ptVals.at(k), treeReaderMuon);
     1655     gr_trkmu_res_eta[k] = EresGraphVsEta(&plots_trkmu_res_eta[k]);
     1656
     1657     h_trkmu_eff_eta = GetEffEta<Track>(branchTrackMuon, branchParticleMuon, "Muon", 13, 0.5*ptVals.at(k), 2.0*ptVals.at(k) ,etaMin, etaMax , treeReaderMuon);
     1658     gr_trkmu_eff_eta[k] = TGraphErrors(h_trkmu_eff_eta);
     1659
     1660     s_pt = Form("#mu^{ #pm} , p_{T} = %.0f GeV",ptVals.at(k));
     1661     if(ptVals.at(k) >= 1000.) s_pt = Form("#mu^{ #pm} , p_{T} = %.0f TeV",ptVals.at(k)/1000.);
     1662
     1663     addResoGraph(mg_trkmu_res_eta, &gr_trkmu_res_eta[k], leg_trkmu_res_eta, markerStyles.at(k), colors.at(k), s_pt );
     1664     addResoGraph(mg_trkmu_eff_eta, &gr_trkmu_eff_eta[k], leg_trkmu_eff_eta, markerStyles.at(k), colors.at(k), s_pt );
     1665
     1666  }
     1667
     1668  TCanvas *c_trkmu_res_pt = new TCanvas("","", 800, 600);
     1669
     1670  mg_trkmu_res_pt->Draw("APE");
     1671  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);
     1672  leg_trkmu_res_pt->Draw();
     1673
     1674  c_trkmu_res_pt->Print(pdfOutput,"pdf");
     1675  c_trkmu_res_pt->Print(figPath+"img_trkmu_res_pt.pdf","pdf");
     1676  c_trkmu_res_pt->Print(figPath+"img_trkmu_res_pt.png","png");
     1677
     1678  TCanvas *c_trkmu_res_eta = new TCanvas("","", 800, 600);
     1679
     1680  mg_trkmu_res_eta->Draw("APE");
     1681  DrawAxis(mg_trkmu_res_eta, leg_trkmu_res_eta, etaMin, etaMax, 0.01, 100, " #eta ", "(track resolution in p_{T})/p_{T} (%)", false, true);
     1682  leg_trkmu_res_eta->Draw();
     1683
     1684  c_trkmu_res_eta->Print(pdfOutput,"pdf");
     1685  c_trkmu_res_eta->Print(figPath+"img_trkmu_res_eta.pdf","pdf");
     1686  c_trkmu_res_eta->Print(figPath+"img_trkmu_res_eta.png","png");
     1687
     1688  TCanvas *c_trkmu_eff_pt = new TCanvas("","", 800, 600);
     1689
     1690  mg_trkmu_eff_pt->Draw("APE");
     1691  DrawAxis(mg_trkmu_eff_pt, leg_trkmu_eff_pt, ptMin, ptMax, 0.0, 100, "p_{T} [GeV]", "tracking efficiency (%)", true, false);
     1692  leg_trkmu_eff_pt->Draw();
     1693
     1694  c_trkmu_eff_pt->Print(pdfOutput,"pdf");
     1695  c_trkmu_eff_pt->Print(figPath+"img_trkmu_eff_pt.pdf","pdf");
     1696  c_trkmu_eff_pt->Print(figPath+"img_trkmu_eff_pt.png","png");
     1697
     1698  TCanvas *c_trkmu_eff_eta = new TCanvas("","", 800, 600);
     1699
     1700  mg_trkmu_eff_eta->Draw("APE");
     1701  DrawAxis(mg_trkmu_eff_eta, leg_trkmu_eff_eta, etaMin, etaMax, 0.0, 100, " #eta ", "tracking efficiency (%)", false, false);
     1702  leg_trkmu_eff_eta->Draw();
     1703
     1704  c_trkmu_eff_eta->Print(pdfOutput,"pdf");
     1705  c_trkmu_eff_eta->Print(figPath+"img_trkmu_eff_eta.pdf","pdf");
     1706  c_trkmu_eff_eta->Print(figPath+"img_trkmu_eff_eta.png","png");
     1707
     1708
     1709  //////////////////////
     1710  // Ecal performance //
     1711  //////////////////////
     1712
     1713
     1714  TMultiGraph *mg_ecal_res_e  = new TMultiGraph("","");
     1715  TMultiGraph *mg_ecal_res_eta = new TMultiGraph("","");
     1716
     1717  TLegend *leg_ecal_res_e = new TLegend(0.55,0.64,0.90,0.90);
     1718  TLegend *leg_ecal_res_eta = new TLegend(0.60,0.59,0.95,0.90);
     1719
     1720  TGraphErrors gr_ecal_res_e[n_etabins], gr_ecal_res_eta[n_ptbins];
     1721
     1722  std::vector<resolPlot> plots_ecal_res_e[n_etabins], plots_ecal_res_eta[n_ptbins];
     1723
     1724  // loop over eta bins
     1725  for (k = 0; k < etaVals.size()-1; k++)
     1726  {
     1727     HistogramsCollection(&plots_ecal_res_e[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "ecal");
     1728     GetEres<Tower>(&plots_ecal_res_e[k], branchTowerPhoton, branchParticlePhoton, 22, etaVals.at(k), etaVals.at(k+1), treeReaderPhoton);
     1729     gr_ecal_res_e[k] = EresGraph(&plots_ecal_res_e[k]);
     1730
     1731     s_etaMin = Form("%.1f",etaVals.at(k));
     1732     s_etaMax = Form("%.1f",etaVals.at(k+1));
     1733
     1734     s_eta = "#gamma , " + s_etaMin + " < | #eta | < "+s_etaMax;
     1735
     1736     gr_ecal_res_e[k].SetName("trkRes_"+s_etaMin+"_"+s_etaMax);
     1737
     1738     addResoGraph(mg_ecal_res_e, &gr_ecal_res_e[k], leg_ecal_res_e, markerStyles.at(k), colors.at(k), s_eta);
     1739  }
     1740
     1741  // loop over pt
     1742  for (k = 0; k < ptVals.size(); k++)
     1743  {
     1744     HistogramsCollectionVsEta(&plots_ecal_res_eta[k], etaMin, etaMax, 0.5*ptVals.at(k), 2.0*ptVals.at(k), "ecal", 0.0, 2.0);
     1745     GetEresVsEta<Tower>(&plots_ecal_res_eta[k], branchTowerPhoton, branchParticlePhoton, 22, 0.5*ptVals.at(k), 2.0*ptVals.at(k), treeReaderPhoton);
     1746     gr_ecal_res_eta[k] = EresGraphVsEta(&plots_ecal_res_eta[k]);
     1747
     1748     s_e = Form("#gamma , E = %.0f GeV",ptVals.at(k));
     1749     if(ptVals.at(k) >= 1000.) s_e = Form("#gamma , E = %.0f TeV",ptVals.at(k)/1000.);
     1750
     1751     addResoGraph(mg_ecal_res_eta, &gr_ecal_res_eta[k], leg_ecal_res_eta, markerStyles.at(k), colors.at(k), s_e );
     1752  }
     1753
     1754  TCanvas *c_ecal_res_e = new TCanvas("","", 800, 600);
     1755
     1756  mg_ecal_res_e->Draw("APE");
     1757 // DrawAxis(mg_ecal_res_e, leg_ecal_res_e, ptMin, ptMax, 0.5, 100, "E [GeV]", "(ECAL resolution in E)/E (%)", true, true);
     1758  DrawAxis(mg_ecal_res_e, leg_ecal_res_e, ptMin, ptMax, 0.0, 20, "E [GeV]", "(ECAL resolution in E)/E (%)", true, false);
     1759  leg_ecal_res_e->Draw();
     1760
     1761  c_ecal_res_e->Print(pdfOutput,"pdf");
     1762  c_ecal_res_e->Print(figPath+"img_ecal_res_e.pdf","pdf");
     1763  c_ecal_res_e->Print(figPath+"img_ecal_res_e.png","png");
     1764
     1765  TCanvas *c_ecal_res_eta = new TCanvas("","", 800, 600);
     1766
     1767  mg_ecal_res_eta->Draw("APE");
     1768  //DrawAxis(mg_ecal_res_eta, leg_ecal_res_eta, etaMin, etaMax, 0.5, 100, " #eta ", "(ECAL resolution in E)/E (%)", false, true);
     1769  DrawAxis(mg_ecal_res_eta, leg_ecal_res_eta, etaMin, etaMax, 0.0, 20, " #eta ", "(ECAL resolution in E)/E (%)", false, false);
     1770  leg_ecal_res_eta->Draw();
     1771
     1772  c_ecal_res_eta->Print(pdfOutput,"pdf");
     1773  c_ecal_res_eta->Print(figPath+"img_ecal_res_eta.pdf","pdf");
     1774  c_ecal_res_eta->Print(figPath+"img_ecal_res_eta.png","png");
     1775
     1776  //////////////////////
     1777  // Hcal performance //
     1778  //////////////////////
     1779
     1780
     1781  TMultiGraph *mg_hcal_res_e  = new TMultiGraph("","");
     1782  TMultiGraph *mg_hcal_res_eta = new TMultiGraph("","");
     1783
     1784  TLegend *leg_hcal_res_e = new TLegend(0.55,0.64,0.90,0.90);
     1785  TLegend *leg_hcal_res_eta = new TLegend(0.60,0.59,0.95,0.90);
     1786
     1787  TGraphErrors gr_hcal_res_e[n_etabins], gr_hcal_res_eta[n_ptbins];
     1788
     1789  std::vector<resolPlot> plots_hcal_res_e[n_etabins], plots_hcal_res_eta[n_ptbins];
     1790
     1791  // loop over eta bins
     1792  for (k = 0; k < etaVals.size()-1; k++)
     1793  {
     1794     HistogramsCollection(&plots_hcal_res_e[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "hcal");
     1795     GetEres<Tower>(&plots_hcal_res_e[k], branchTowerNeutralHadron, branchParticleNeutralHadron, 2112, etaVals.at(k), etaVals.at(k+1), treeReaderNeutralHadron);
     1796
     1797     gr_hcal_res_e[k] = EresGraph(&plots_hcal_res_e[k]);
     1798
     1799     s_etaMin = Form("%.1f",etaVals.at(k));
     1800     s_etaMax = Form("%.1f",etaVals.at(k+1));
     1801
     1802     s_eta = "n , " + s_etaMin + " < | #eta | < "+s_etaMax;
     1803
     1804     gr_hcal_res_e[k].SetName("trkRes_"+s_etaMin+"_"+s_etaMax);
     1805
     1806     addResoGraph(mg_hcal_res_e, &gr_hcal_res_e[k], leg_hcal_res_e, markerStyles.at(k), colors.at(k), s_eta);
     1807  }
     1808
     1809  // loop over pt
     1810  for (k = 0; k < ptVals.size(); k++)
     1811  {
     1812     HistogramsCollectionVsEta(&plots_hcal_res_eta[k], etaMin, etaMax, 0.5*ptVals.at(k), 2.0*ptVals.at(k), "hcal", 0.0, 2.0);
     1813     GetEresVsEta<Tower>(&plots_hcal_res_eta[k], branchTowerNeutralHadron, branchParticleNeutralHadron, 2112, 0.5*ptVals.at(k), 2.0*ptVals.at(k), treeReaderNeutralHadron);
     1814     gr_hcal_res_eta[k] = EresGraphVsEta(&plots_hcal_res_eta[k]);
     1815
     1816     s_e = Form("n , E = %.0f GeV",ptVals.at(k));
     1817     if(ptVals.at(k) >= 1000.) s_e = Form("n , E = %.0f TeV",ptVals.at(k)/1000.);
     1818
     1819     addResoGraph(mg_hcal_res_eta, &gr_hcal_res_eta[k], leg_hcal_res_eta, markerStyles.at(k), colors.at(k), s_e );
     1820  }
     1821
     1822
     1823  TCanvas *c_hcal_res_e = new TCanvas("","", 800, 600);
     1824
     1825  mg_hcal_res_e->Draw("APE");
     1826  //DrawAxis(mg_hcal_res_e, leg_hcal_res_e, ptMin, ptMax, 1, 100, "E [GeV]", "(HCAL resolution in E)/E (%)", true, true);
     1827  DrawAxis(mg_hcal_res_e, leg_hcal_res_e, ptMin, ptMax, 0.0, 50, "E [GeV]", "(HCAL resolution in E)/E (%)", true, false);
     1828  leg_hcal_res_e->Draw();
     1829
     1830  c_hcal_res_e->Print(pdfOutput,"pdf");
     1831  c_hcal_res_e->Print(figPath+"img_hcal_res_e.pdf","pdf");
     1832  c_hcal_res_e->Print(figPath+"img_hcal_res_e.png","png");
     1833
     1834  TCanvas *c_hcal_res_eta = new TCanvas("","", 800, 600);
     1835
     1836  mg_hcal_res_eta->Draw("APE");
     1837  //DrawAxis(mg_hcal_res_eta, leg_hcal_res_eta, etaMin, etaMax, 1, 100, " #eta ", "(HCAL resolution in E)/E (%)", false, true);
     1838  DrawAxis(mg_hcal_res_eta, leg_hcal_res_eta, etaMin, etaMax, 0.0, 50, " #eta ", "(HCAL resolution in E)/E (%)", false, false);
     1839  leg_hcal_res_eta->Draw();
     1840
     1841  c_hcal_res_eta->Print(pdfOutput,"pdf");
     1842  c_hcal_res_eta->Print(figPath+"img_hcal_res_eta.pdf","pdf");
     1843  c_hcal_res_eta->Print(figPath+"img_hcal_res_eta.png","png");
     1844
     1845  ////////////////////
     1846  // PF - electrons //
     1847  ////////////////////
     1848
     1849  TMultiGraph *mg_pfele_res_e[n_etabins];
     1850  TMultiGraph *mg_pfele_res_eta[n_ptbins];
     1851
     1852  TLegend *leg_pfele_res_e[n_etabins];
     1853  TLegend *leg_pfele_res_eta[n_ptbins];
     1854
     1855  TGraphErrors gr_pfele_res_e[n_etabins];
     1856  TGraphErrors gr_pfele_res_eta[n_ptbins];
     1857
     1858  std::vector<resolPlot> plots_pfele_res_e[n_etabins], plots_pfele_res_eta[n_ptbins];
     1859
     1860  TCanvas *c_pfele_res_e[n_etabins];
     1861  TCanvas *c_pfele_res_eta[n_ptbins];
     1862
     1863
     1864  // loop over eta bins
     1865  for (k = 0; k < etaVals.size()-1; k++)
     1866  {
     1867     mg_pfele_res_e[k] = new TMultiGraph("","");
     1868     leg_pfele_res_e[k] = new TLegend(0.40,0.60,0.75,0.90);
     1869
     1870     HistogramsCollection(&plots_pfele_res_e[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "pfele");
     1871     GetEres<Electron>(&plots_pfele_res_e[k], branchElectron, branchParticleElectron, 11, etaVals.at(k), etaVals.at(k+1), treeReaderElectron);
     1872     gr_pfele_res_e[k] = EresGraph(&plots_pfele_res_e[k]);
     1873
     1874     s_etaMin = Form("%.1f",etaVals.at(k));
     1875     s_etaMax = Form("%.1f",etaVals.at(k+1));
     1876     s_eta    = "e^{ #pm}, "+ s_etaMin + " < | #eta | < " + s_etaMax;
     1877
     1878     leg_pfele_res_e[k]->SetTextFont(132);
     1879     leg_pfele_res_e[k]->SetHeader(s_eta);
     1880
     1881     addResoGraph(mg_pfele_res_e[k], &gr_ecal_res_e[k], leg_pfele_res_e[k], markerStyles.at(0), colors.at(0), "ecal");
     1882     addResoGraph(mg_pfele_res_e[k], &gr_trkele_res_pt[k], leg_pfele_res_e[k], markerStyles.at(1), colors.at(1), "track");
     1883     addResoGraph(mg_pfele_res_e[k], &gr_pfele_res_e[k], leg_pfele_res_e[k], markerStyles.at(2), colors.at(2), "p-flow");
     1884
     1885     c_pfele_res_e[k] = new TCanvas("","", 800, 600);
     1886
     1887     mg_pfele_res_e[k]->Draw("APE");
     1888     //DrawAxis(mg_pfele_res_e[k], leg_pfele_res_e[k], ptMin, ptMax, 0.1, 100, "E [GeV]", "(resolution in E)/E (%)", true, true);
     1889     DrawAxis(mg_pfele_res_e[k], leg_pfele_res_e[k], ptMin, ptMax, 0.0, 20, "E [GeV]", "(resolution in E)/E (%)", true, false);
     1890     leg_pfele_res_e[k]->Draw();
     1891
     1892     TString s_etarange = "eta_"+s_etaMin+"_"+s_etaMax+"_";
     1893
     1894     c_pfele_res_e[k]->Print(pdfOutput,"pdf");
     1895     c_pfele_res_e[k]->Print(figPath+"img_pfele_res_"+s_etarange+"e.pdf","pdf");
     1896     c_pfele_res_e[k]->Print(figPath+"img_pfele_res_"+s_etarange+"e.png","png");
     1897
     1898  }
     1899
     1900
     1901    // loop over eta bins
     1902  for (k = 0; k < ptVals.size(); k++)
     1903  {
     1904
     1905     mg_pfele_res_eta[k] = new TMultiGraph("","");
     1906     leg_pfele_res_eta[k] = new TLegend(0.40,0.60,0.75,0.90);
     1907
     1908     HistogramsCollectionVsEta(&plots_pfele_res_eta[k], etaMin, etaMax, 0.5*ptVals.at(k), 2.0*ptVals.at(k), "pfele", 0.0, 2.0);
     1909     GetEresVsEta<Electron>(&plots_pfele_res_eta[k], branchElectron, branchParticleElectron, 11, 0.5*ptVals.at(k), 2.0*ptVals.at(k), treeReaderElectron);
     1910     gr_pfele_res_eta[k] = EresGraphVsEta(&plots_pfele_res_eta[k]);
     1911
     1912     s_e = Form("e^{ #pm}, E = %.0f GeV",ptVals.at(k));
     1913     if(ptVals.at(k) >= 1000.) s_e = Form("e^{ #pm}, E = %.0f TeV",ptVals.at(k)/1000.);
     1914
     1915
     1916     leg_pfele_res_eta[k]->SetTextFont(132);
     1917     leg_pfele_res_eta[k]->SetHeader(s_e);
     1918
     1919     addResoGraph(mg_pfele_res_eta[k], &gr_ecal_res_eta[k], leg_pfele_res_eta[k], markerStyles.at(0), colors.at(0), "ecal");
     1920     addResoGraph(mg_pfele_res_eta[k], &gr_trkele_res_eta[k], leg_pfele_res_eta[k], markerStyles.at(1), colors.at(1), "track");
     1921     addResoGraph(mg_pfele_res_eta[k], &gr_pfele_res_eta[k], leg_pfele_res_eta[k], markerStyles.at(2), colors.at(2), "p-flow");
     1922
     1923     c_pfele_res_eta[k] = new TCanvas("","", 800, 600);
     1924
     1925     mg_pfele_res_eta[k]->Draw("APE");
     1926     //DrawAxis(mg_pfele_res_eta[k], leg_pfele_res_eta[k], etaMin, etaMax, 0.1, 1000, "#eta", "(resolution in E)/E (%)", false, true);
     1927     DrawAxis(mg_pfele_res_eta[k], leg_pfele_res_eta[k], etaMin, etaMax, 0.0, 50, "#eta", "(resolution in E)/E (%)", false, false);
     1928     leg_pfele_res_eta[k]->Draw();
     1929
     1930     TString s_ptrange = Form("pt_%.0f_",ptVals.at(k));
     1931
     1932     c_pfele_res_eta[k]->Print(pdfOutput,"pdf");
     1933     c_pfele_res_eta[k]->Print(figPath+"img_pfele_res_"+s_ptrange+"eta.pdf","pdf");
     1934     c_pfele_res_eta[k]->Print(figPath+"img_pfele_res_"+s_ptrange+"eta.png","png");
     1935
     1936  }
     1937
     1938  /////////////////
     1939  // PF - Pions  //
     1940  /////////////////
     1941
     1942  TMultiGraph *mg_pfpi_res_e[n_etabins];
     1943  TMultiGraph *mg_pfpi_res_eta[n_ptbins];
     1944
     1945  TLegend *leg_pfpi_res_e[n_etabins];
     1946  TLegend *leg_pfpi_res_eta[n_ptbins];
     1947
     1948  TGraphErrors gr_pfpi_res_e[n_etabins];
     1949  TGraphErrors gr_pfpi_res_eta[n_ptbins];
     1950
     1951  std::vector<resolPlot> plots_pfpi_res_e[n_etabins], plots_pfpi_res_eta[n_ptbins];
     1952
     1953  TCanvas *c_pfpi_res_e[n_etabins];
     1954  TCanvas *c_pfpi_res_eta[n_ptbins];
     1955
     1956
     1957  // loop over eta bins
     1958  for (k = 0; k < etaVals.size()-1; k++)
     1959  {
     1960     mg_pfpi_res_e[k] = new TMultiGraph("","");
     1961     leg_pfpi_res_e[k] = new TLegend(0.40,0.60,0.75,0.90);
     1962
     1963     HistogramsCollection(&plots_pfpi_res_e[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "pfpi");
     1964     GetEres<Track>(&plots_pfpi_res_e[k], branchPion, branchParticlePion, 211, etaVals.at(k), etaVals.at(k+1), treeReaderPion);
     1965     gr_pfpi_res_e[k] = EresGraph(&plots_pfpi_res_e[k]);
     1966
     1967     s_etaMin = Form("%.1f",etaVals.at(k));
     1968     s_etaMax = Form("%.1f",etaVals.at(k+1));
     1969     s_eta    = "#pi^{ #pm}, "+ s_etaMin + " < | #eta | < " + s_etaMax;
     1970
     1971     leg_pfpi_res_e[k]->SetTextFont(132);
     1972     leg_pfpi_res_e[k]->SetHeader(s_eta);
     1973
     1974     addResoGraph(mg_pfpi_res_e[k], &gr_hcal_res_e[k], leg_pfpi_res_e[k], markerStyles.at(0), colors.at(0), "hcal");
     1975     addResoGraph(mg_pfpi_res_e[k], &gr_trkpi_res_pt[k], leg_pfpi_res_e[k], markerStyles.at(1), colors.at(1), "track");
     1976     addResoGraph(mg_pfpi_res_e[k], &gr_pfpi_res_e[k], leg_pfpi_res_e[k], markerStyles.at(2), colors.at(2), "p-flow");
     1977
     1978     c_pfpi_res_e[k] = new TCanvas("","", 800, 600);
     1979
     1980     mg_pfpi_res_e[k]->Draw("APE");
     1981     //DrawAxis(mg_pfpi_res_e[k], leg_pfpi_res_e[k], ptMin, ptMax, 0.1, 100, "E [GeV]", "(resolution in E)/E (%)", true, true);
     1982     DrawAxis(mg_pfpi_res_e[k], leg_pfpi_res_e[k], ptMin, ptMax, 0.1, 50, "E [GeV]", "(resolution in E)/E (%)", true, false);
     1983     leg_pfpi_res_e[k]->Draw();
     1984
     1985     TString s_etarange = "eta_"+s_etaMin+"_"+s_etaMax+"_";
     1986
     1987     c_pfpi_res_e[k]->Print(pdfOutput,"pdf");
     1988     c_pfpi_res_e[k]->Print(figPath+"img_pfpi_res_"+s_etarange+"e.pdf","pdf");
     1989     c_pfpi_res_e[k]->Print(figPath+"img_pfpi_res_"+s_etarange+"e.png","png");
     1990
     1991  }
     1992
     1993
     1994    // loop over eta bins
     1995  for (k = 0; k < ptVals.size(); k++)
     1996  {
     1997
     1998     mg_pfpi_res_eta[k] = new TMultiGraph("","");
     1999     leg_pfpi_res_eta[k] = new TLegend(0.40,0.60,0.75,0.90);
     2000
     2001     HistogramsCollectionVsEta(&plots_pfpi_res_eta[k], etaMin, etaMax, 0.5*ptVals.at(k), 2.0*ptVals.at(k), "pfpi", 0.0, 2.0);
     2002     GetEresVsEta<Track>(&plots_pfpi_res_eta[k], branchPion, branchParticlePion, 211, 0.5*ptVals.at(k), 2.0*ptVals.at(k), treeReaderPion);
     2003     gr_pfpi_res_eta[k] = EresGraphVsEta(&plots_pfpi_res_eta[k]);
     2004
     2005     s_e = Form("#pi^{ #pm}, E = %.0f GeV",ptVals.at(k));
     2006     if(ptVals.at(k) >= 1000.) s_e = Form("#pi^{ #pm}, E = %.0f TeV",ptVals.at(k)/1000.);
     2007
     2008     leg_pfpi_res_eta[k]->SetTextFont(132);
     2009     leg_pfpi_res_eta[k]->SetHeader(s_e);
     2010
     2011     addResoGraph(mg_pfpi_res_eta[k], &gr_hcal_res_eta[k], leg_pfpi_res_eta[k], markerStyles.at(0), colors.at(0), "hcal");
     2012     addResoGraph(mg_pfpi_res_eta[k], &gr_trkpi_res_eta[k], leg_pfpi_res_eta[k], markerStyles.at(1), colors.at(1), "track");
     2013     addResoGraph(mg_pfpi_res_eta[k], &gr_pfpi_res_eta[k], leg_pfpi_res_eta[k], markerStyles.at(2), colors.at(2), "p-flow");
     2014
     2015     c_pfpi_res_eta[k] = new TCanvas("","", 800, 600);
     2016
     2017     mg_pfpi_res_eta[k]->Draw("APE");
     2018     //DrawAxis(mg_pfpi_res_eta[k], leg_pfpi_res_eta[k], etaMin, etaMax, 0.1, 1000, "#eta", "(resolution in E)/E (%)", false, true);
     2019     DrawAxis(mg_pfpi_res_eta[k], leg_pfpi_res_eta[k], etaMin, etaMax, 0.0, 50, "#eta", "(resolution in E)/E (%)", false, false);
     2020     leg_pfpi_res_eta[k]->Draw();
     2021
     2022     TString s_ptrange = Form("pt_%.0f_",ptVals.at(k));
     2023
     2024     c_pfpi_res_eta[k]->Print(pdfOutput,"pdf");
     2025     c_pfpi_res_eta[k]->Print(figPath+"img_pfpi_res_"+s_ptrange+"eta.pdf","pdf");
     2026     c_pfpi_res_eta[k]->Print(figPath+"img_pfpi_res_"+s_ptrange+"eta.png","png");
     2027
     2028  }
     2029
     2030
     2031  /////////////////
     2032  // PF - jets   //
     2033  /////////////////
     2034
     2035  TMultiGraph *mg_pfjet_res_e[n_etabins];
     2036  TMultiGraph *mg_pfjet_res_eta[n_ptbins];
     2037
     2038  TLegend *leg_pfjet_res_e[n_etabins];
     2039  TLegend *leg_pfjet_res_eta[n_ptbins];
     2040
     2041  TGraphErrors gr_pfjet_res_e[n_etabins];
     2042  TGraphErrors gr_pfjet_res_eta[n_ptbins];
     2043
     2044  TGraphErrors gr_cajet_res_e[n_etabins];
     2045  TGraphErrors gr_cajet_res_eta[n_ptbins];
     2046
     2047  std::vector<resolPlot> plots_pfjet_res_e[n_etabins], plots_pfjet_res_eta[n_ptbins];
     2048  std::vector<resolPlot> plots_cajet_res_e[n_etabins], plots_cajet_res_eta[n_ptbins];
     2049
     2050  TCanvas *c_pfjet_res_e[n_etabins];
     2051  TCanvas *c_pfjet_res_eta[n_ptbins];
     2052
     2053
     2054  // loop over eta bins
     2055  for (k = 0; k < etaVals.size()-1; k++)
     2056  {
     2057
     2058     mg_pfjet_res_e[k] = new TMultiGraph("","");
     2059     leg_pfjet_res_e[k] = new TLegend(0.60,0.60,0.90,0.90);
     2060
     2061     HistogramsCollection(&plots_pfjet_res_e[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "pfjet");
     2062     HistogramsCollection(&plots_cajet_res_e[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "cajet");
     2063
     2064     GetJetsEres(&plots_pfjet_res_e[k], branchPFJet, branchGenJet, treeReaderJet, etaVals.at(k), etaVals.at(k+1));
     2065     GetJetsEres(&plots_cajet_res_e[k], branchCaloJet, branchGenJet, treeReaderJet, etaVals.at(k), etaVals.at(k+1));
     2066
     2067     gr_pfjet_res_e[k] = EresGraph(&plots_pfjet_res_e[k]);
     2068     gr_cajet_res_e[k] = EresGraph(&plots_cajet_res_e[k]);
     2069
     2070     s_etaMin = Form("%.1f",etaVals.at(k));
     2071     s_etaMax = Form("%.1f",etaVals.at(k+1));
     2072     s_eta    = "jets, "+ s_etaMin + " < | #eta | < " + s_etaMax;
     2073
     2074     leg_pfjet_res_e[k]->SetTextFont(132);
     2075     leg_pfjet_res_e[k]->SetHeader(s_eta);
     2076
     2077     addResoGraph(mg_pfjet_res_e[k], &gr_cajet_res_e[k], leg_pfjet_res_e[k], markerStyles.at(0), colors.at(0), "calo");
     2078     addResoGraph(mg_pfjet_res_e[k], &gr_pfjet_res_e[k], leg_pfjet_res_e[k], markerStyles.at(1), colors.at(1), "p-flow");
     2079
     2080     c_pfjet_res_e[k] = new TCanvas("","", 800, 600);
     2081
     2082     mg_pfjet_res_e[k]->Draw("APE");
     2083     //DrawAxis(mg_pfjet_res_e[k], leg_pfjet_res_e[k], 10, ptMax, 0.5, 100, "E [GeV]", "(resolution in E)/E (%)", true, true);
     2084     DrawAxis(mg_pfjet_res_e[k], leg_pfjet_res_e[k], 10, ptMax, 0.0, 30, "E [GeV]", "(resolution in E)/E (%)", true, false);
     2085     leg_pfjet_res_e[k]->Draw();
     2086
     2087     TString s_etarange = "eta_"+s_etaMin+"_"+s_etaMax+"_";
     2088
     2089     c_pfjet_res_e[k]->Print(pdfOutput,"pdf");
     2090     c_pfjet_res_e[k]->Print(figPath+"img_pfjet_res_"+s_etarange+"e.pdf","pdf");
     2091     c_pfjet_res_e[k]->Print(figPath+"img_pfjet_res_"+s_etarange+"e.png","png");
     2092
     2093  }
     2094
     2095
     2096    // loop over eta bins
     2097  for (k = 0; k < ptVals.size(); k++)
     2098  {
     2099
     2100     mg_pfjet_res_eta[k] = new TMultiGraph("","");
     2101     leg_pfjet_res_eta[k] = new TLegend(0.40,0.60,0.75,0.90);
     2102
     2103     HistogramsCollectionVsEta(&plots_pfjet_res_eta[k], etaMin, etaMax, 0.5*ptVals.at(k), 2.0*ptVals.at(k), "pfjet", 0.0, 2.0);
     2104     HistogramsCollectionVsEta(&plots_cajet_res_eta[k], etaMin, etaMax, 0.5*ptVals.at(k), 2.0*ptVals.at(k), "cajet", 0.0, 2.0);
     2105
     2106     GetJetsEresVsEta(&plots_pfjet_res_eta[k], branchPFJet, branchGenJet, treeReaderJet, 0.5*ptVals.at(k), 2.0*ptVals.at(k));
     2107     GetJetsEresVsEta(&plots_cajet_res_eta[k], branchCaloJet, branchGenJet, treeReaderJet, 0.5*ptVals.at(k), 2.0*ptVals.at(k));
     2108
     2109     gr_pfjet_res_eta[k] = EresGraphVsEta(&plots_pfjet_res_eta[k]);
     2110     gr_cajet_res_eta[k] = EresGraphVsEta(&plots_cajet_res_eta[k]);
     2111
     2112     s_e = Form("jets, E = %.0f GeV",ptVals.at(k));
     2113     if(ptVals.at(k) >= 1000.) s_e = Form("jets, E = %.0f TeV",ptVals.at(k)/1000.);
     2114
     2115     leg_pfjet_res_eta[k]->SetTextFont(132);
     2116     leg_pfjet_res_eta[k]->SetHeader(s_e);
     2117
     2118     addResoGraph(mg_pfjet_res_eta[k], &gr_cajet_res_eta[k], leg_pfjet_res_eta[k], markerStyles.at(0), colors.at(0), "calo");
     2119     addResoGraph(mg_pfjet_res_eta[k], &gr_pfjet_res_eta[k], leg_pfjet_res_eta[k], markerStyles.at(1), colors.at(1), "p-flow");
     2120
     2121     c_pfjet_res_eta[k] = new TCanvas("","", 800, 600);
     2122
     2123     mg_pfjet_res_eta[k]->Draw("APE");
     2124     //DrawAxis(mg_pfjet_res_eta[k], leg_pfjet_res_eta[k], etaMin, etaMax, 0.1, 1000, "#eta", "(resolution in E)/E (%)", false, true);
     2125     DrawAxis(mg_pfjet_res_eta[k], leg_pfjet_res_eta[k], etaMin, etaMax, 0.0, 50, "#eta", "(resolution in E)/E (%)", false, false);
     2126     leg_pfjet_res_eta[k]->Draw();
     2127
     2128     TString s_ptrange = Form("pt_%.0f_",ptVals.at(k));
     2129
     2130     c_pfjet_res_eta[k]->Print(pdfOutput,"pdf");
     2131     c_pfjet_res_eta[k]->Print(figPath+"img_pfjet_res_"+s_ptrange+"eta.pdf","pdf");
     2132     c_pfjet_res_eta[k]->Print(figPath+"img_pfjet_res_"+s_ptrange+"eta.png","png");
     2133
     2134  }
     2135
     2136
     2137    /////////////////////
     2138    // PF Missing ET  ///
     2139    /////////////////////
     2140
     2141    TMultiGraph *mg_met_res_ht = new TMultiGraph("","");
     2142    TLegend *leg_met_res_ht    = new TLegend(0.60,0.22,0.90,0.42);
     2143
     2144    std::vector<resolPlot> plots_pfmet, plots_camet;
     2145
     2146    HistogramsCollection(&plots_pfmet, TMath::Log10(ptMin), TMath::Log10(ptMax), "pfMET", -500, 500);
     2147    HistogramsCollection(&plots_camet, TMath::Log10(ptMin), TMath::Log10(ptMax), "caMET", -500, 500);
     2148
     2149    GetMetres(&plots_pfmet, branchGenScalarHT, branchMet, branchPFJet, treeReaderJet);
     2150    GetMetres(&plots_camet, branchGenScalarHT, branchCaloMet, branchCaloJet, treeReaderJet);
     2151
     2152    TGraphErrors gr_pfmet_res_ht = MetResGraph(&plots_pfmet, true);
     2153    TGraphErrors gr_camet_res_ht = MetResGraph(&plots_camet, true);
     2154
     2155    addResoGraph(mg_met_res_ht, &gr_camet_res_ht, leg_met_res_ht, markerStyles.at(0), colors.at(0), "calo");
     2156    addResoGraph(mg_met_res_ht, &gr_pfmet_res_ht, leg_met_res_ht, markerStyles.at(1), colors.at(1), "p-flow");
     2157
     2158    TCanvas *c_met_res_ht = new TCanvas("","", 800, 600);
     2159
     2160    mg_met_res_ht->Draw("APE");
     2161    DrawAxis(mg_met_res_ht, leg_met_res_ht, 10, 10000, 0.1, 1000, " #sum p_{T} [GeV]", "resolution in E_{x,y} [GeV]", true, true);
     2162
     2163    leg_met_res_ht->Draw();
     2164
     2165    c_met_res_ht->Print(pdfOutput,"pdf");
     2166    c_met_res_ht->Print(figPath+"img_met_res_ht.pdf","pdf");
     2167    c_met_res_ht->Print(figPath+"img_met_res_ht.png","png");
     2168
     2169
     2170    /////////////////////////////////////////
     2171    // Electron Reconstruction Efficiency ///
     2172    /////////////////////////////////////////
     2173
     2174    TMultiGraph *mg_recele_eff_pt  = new TMultiGraph("","");
     2175    TMultiGraph *mg_recele_eff_eta = new TMultiGraph("","");
     2176
     2177    TLegend *leg_recele_eff_pt  = new TLegend(0.55,0.22,0.90,0.48);
     2178    TLegend *leg_recele_eff_eta = new TLegend(0.55,0.22,0.90,0.48);
     2179
     2180    TGraphErrors gr_recele_eff_pt[n_etabins], gr_recele_eff_eta[n_ptbins];
     2181    TH1D* h_recele_eff_pt, *h_recele_eff_eta;
     2182
     2183    // loop over eta bins
     2184    for (k = 0; k < etaVals.size()-1; k++)
     2185    {
     2186
     2187       h_recele_eff_pt = GetEffPt<Electron>(branchElectron, branchParticleElectron, "Electron", 11, ptMin, ptMax, etaVals.at(k), etaVals.at(k+1), treeReaderElectron);
     2188       gr_recele_eff_pt[k] = TGraphErrors(h_recele_eff_pt);
     2189
     2190       s_etaMin = Form("%.1f",etaVals.at(k));
     2191       s_etaMax = Form("%.1f",etaVals.at(k+1));
     2192
     2193       s_eta = "e^{ #pm} , " + s_etaMin + " < | #eta | < "+s_etaMax;
     2194
     2195       gr_recele_eff_pt[k].SetName("recEff_"+s_etaMin+"_"+s_etaMax);
     2196
     2197       addResoGraph(mg_recele_eff_pt, &gr_recele_eff_pt[k], leg_recele_eff_pt, markerStyles.at(k), colors.at(k), s_eta);
     2198    }
     2199
     2200    // loop over pt
     2201    for (k = 0; k < ptVals.size(); k++)
     2202    {
     2203       h_recele_eff_eta = GetEffEta<Electron>(branchElectron, branchParticleElectron, "Electron", 11, 0.5*ptVals.at(k), 2.0*ptVals.at(k) ,etaMin, etaMax , treeReaderElectron);
     2204       gr_recele_eff_eta[k] = TGraphErrors(h_recele_eff_eta);
     2205
     2206       s_pt = Form("e^{ #pm} , p_{T} = %.0f GeV",ptVals.at(k));
     2207       if(ptVals.at(k) >= 1000.) s_pt = Form("e^{ #pm} , p_{T} = %.0f TeV",ptVals.at(k)/1000.);
     2208
     2209       addResoGraph(mg_recele_eff_eta, &gr_recele_eff_eta[k], leg_recele_eff_eta, markerStyles.at(k), colors.at(k), s_pt );
     2210    }
     2211
     2212    TCanvas *c_recele_eff_pt = new TCanvas("","", 800, 600);
     2213
     2214    mg_recele_eff_pt->Draw("APE");
     2215    DrawAxis(mg_recele_eff_pt, leg_recele_eff_pt, ptMin, ptMax, 0.0, 100, "p_{T} [GeV]", "reconstruction efficiency (%)", true, false);
     2216    leg_recele_eff_pt->Draw();
     2217
     2218    c_recele_eff_pt->Print(pdfOutput,"pdf");
     2219    c_recele_eff_pt->Print(figPath+"img_recele_eff_pt.pdf","pdf");
     2220    c_recele_eff_pt->Print(figPath+"img_recele_eff_pt.png","png");
     2221
     2222    TCanvas *c_recele_eff_eta = new TCanvas("","", 800, 600);
     2223
     2224    mg_recele_eff_eta->Draw("APE");
     2225    DrawAxis(mg_recele_eff_eta, leg_recele_eff_eta, etaMin, etaMax, 0.0, 100, " #eta ", "reconstruction efficiency (%)", false, false);
     2226    leg_recele_eff_eta->Draw();
     2227
     2228    c_recele_eff_eta->Print(pdfOutput,"pdf");
     2229    c_recele_eff_eta->Print(figPath+"img_recele_eff_eta.pdf","pdf");
     2230    c_recele_eff_eta->Print(figPath+"img_recele_eff_eta.png","png");
     2231
     2232
     2233    /////////////////////////////////////////
     2234    // Muon Reconstruction Efficiency ///
     2235    /////////////////////////////////////////
     2236
     2237    TMultiGraph *mg_recmu_eff_pt  = new TMultiGraph("","");
     2238    TMultiGraph *mg_recmu_eff_eta = new TMultiGraph("","");
     2239
     2240    TLegend *leg_recmu_eff_pt  = new TLegend(0.55,0.22,0.90,0.48);
     2241    TLegend *leg_recmu_eff_eta = new TLegend(0.55,0.22,0.90,0.48);
     2242
     2243    TGraphErrors gr_recmu_eff_pt[n_etabins], gr_recmu_eff_eta[n_ptbins];
     2244    TH1D* h_recmu_eff_pt, *h_recmu_eff_eta;
     2245
     2246    // loop over eta bins
     2247    for (k = 0; k < etaVals.size()-1; k++)
     2248    {
     2249
     2250       h_recmu_eff_pt = GetEffPt<Muon>(branchMuon, branchParticleMuon, "muon", 13, ptMin, ptMax, etaVals.at(k), etaVals.at(k+1), treeReaderMuon);
     2251       gr_recmu_eff_pt[k] = TGraphErrors(h_recmu_eff_pt);
     2252
     2253       s_etaMin = Form("%.1f",etaVals.at(k));
     2254       s_etaMax = Form("%.1f",etaVals.at(k+1));
     2255
     2256       s_eta = "#mu^{ #pm} , " + s_etaMin + " < | #eta | < "+s_etaMax;
     2257
     2258       gr_recmu_eff_pt[k].SetName("recEff_"+s_etaMin+"_"+s_etaMax);
     2259
     2260       addResoGraph(mg_recmu_eff_pt, &gr_recmu_eff_pt[k], leg_recmu_eff_pt, markerStyles.at(k), colors.at(k), s_eta);
     2261    }
     2262
     2263    // loop over pt
     2264    for (k = 0; k < ptVals.size(); k++)
     2265    {
     2266       h_recmu_eff_eta = GetEffEta<Muon>(branchMuon, branchParticleMuon, "muon", 13, 0.5*ptVals.at(k), 2.0*ptVals.at(k) ,etaMin, etaMax , treeReaderMuon);
     2267       gr_recmu_eff_eta[k] = TGraphErrors(h_recmu_eff_eta);
     2268
     2269       s_pt = Form("#mu^{ #pm} , p_{T} = %.0f GeV",ptVals.at(k));
     2270       if(ptVals.at(k) >= 1000.) s_pt = Form("#mu^{ #pm} , p_{T} = %.0f TeV",ptVals.at(k)/1000.);
     2271
     2272       addResoGraph(mg_recmu_eff_eta, &gr_recmu_eff_eta[k], leg_recmu_eff_eta, markerStyles.at(k), colors.at(k), s_pt );
     2273    }
     2274
     2275    TCanvas *c_recmu_eff_pt = new TCanvas("","", 800, 600);
     2276
     2277    mg_recmu_eff_pt->Draw("APE");
     2278    DrawAxis(mg_recmu_eff_pt, leg_recmu_eff_pt, ptMin, ptMax, 0.0, 100, "p_{T} [GeV]", "reconstruction efficiency (%)", true, false);
     2279    leg_recmu_eff_pt->Draw();
     2280
     2281    c_recmu_eff_pt->Print(pdfOutput,"pdf");
     2282    c_recmu_eff_pt->Print(figPath+"img_recmu_eff_pt.pdf","pdf");
     2283    c_recmu_eff_pt->Print(figPath+"img_recmu_eff_pt.png","png");
     2284
     2285    TCanvas *c_recmu_eff_eta = new TCanvas("","", 800, 600);
     2286
     2287    mg_recmu_eff_eta->Draw("APE");
     2288    DrawAxis(mg_recmu_eff_eta, leg_recmu_eff_eta, etaMin, etaMax, 0.0, 100, " #eta ", "reconstruction efficiency (%)", false, false);
     2289    leg_recmu_eff_eta->Draw();
     2290
     2291    c_recmu_eff_eta->Print(pdfOutput,"pdf");
     2292    c_recmu_eff_eta->Print(figPath+"img_recmu_eff_eta.pdf","pdf");
     2293    c_recmu_eff_eta->Print(figPath+"img_recmu_eff_eta.png","png");
     2294
     2295
     2296    /////////////////////////////////////////
     2297    // Photon Reconstruction Efficiency   ///
     2298    /////////////////////////////////////////
     2299
     2300    TMultiGraph *mg_recpho_eff_pt  = new TMultiGraph("","");
     2301    TMultiGraph *mg_recpho_eff_eta = new TMultiGraph("","");
     2302
     2303    TLegend *leg_recpho_eff_pt  = new TLegend(0.55,0.22,0.90,0.48);
     2304    TLegend *leg_recpho_eff_eta = new TLegend(0.55,0.22,0.90,0.48);
     2305
     2306    TGraphErrors gr_recpho_eff_pt[n_etabins], gr_recpho_eff_eta[n_ptbins];
     2307    TH1D* h_recpho_eff_pt, *h_recpho_eff_eta;
     2308
     2309    // loop over eta bins
     2310    for (k = 0; k < etaVals.size()-1; k++)
     2311    {
     2312
     2313       h_recpho_eff_pt = GetEffPt<Photon>(branchPhoton, branchParticlePhoton, "Photon", 22, ptMin, ptMax, etaVals.at(k), etaVals.at(k+1), treeReaderPhoton);
     2314       gr_recpho_eff_pt[k] = TGraphErrors(h_recpho_eff_pt);
     2315
     2316       s_etaMin = Form("%.1f",etaVals.at(k));
     2317       s_etaMax = Form("%.1f",etaVals.at(k+1));
     2318
     2319       s_eta = "#gamma , " + s_etaMin + " < | #eta | < "+s_etaMax;
     2320
     2321       gr_recpho_eff_pt[k].SetName("recEff_"+s_etaMin+"_"+s_etaMax);
     2322
     2323       addResoGraph(mg_recpho_eff_pt, &gr_recpho_eff_pt[k], leg_recpho_eff_pt, markerStyles.at(k), colors.at(k), s_eta);
     2324    }
     2325
     2326    // loop over pt
     2327    for (k = 0; k < ptVals.size(); k++)
     2328    {
     2329       h_recpho_eff_eta = GetEffEta<Photon>(branchPhoton, branchParticlePhoton, "Photon", 22, 0.5*ptVals.at(k), 2.0*ptVals.at(k) ,etaMin, etaMax , treeReaderPhoton);
     2330       gr_recpho_eff_eta[k] = TGraphErrors(h_recpho_eff_eta);
     2331
     2332       s_pt = Form("#gamma , p_{T} = %.0f GeV",ptVals.at(k));
     2333       if(ptVals.at(k) >= 1000.) s_pt = Form("#gamma , p_{T} = %.0f TeV",ptVals.at(k)/1000.);
     2334
     2335       addResoGraph(mg_recpho_eff_eta, &gr_recpho_eff_eta[k], leg_recpho_eff_eta, markerStyles.at(k), colors.at(k), s_pt );
     2336    }
     2337
     2338    TCanvas *c_recpho_eff_pt = new TCanvas("","", 800, 600);
     2339
     2340    mg_recpho_eff_pt->Draw("APE");
     2341    DrawAxis(mg_recpho_eff_pt, leg_recpho_eff_pt, ptMin, ptMax, 0.0, 100, "p_{T} [GeV]", "reconstruction efficiency (%)", true, false);
     2342    leg_recpho_eff_pt->Draw();
     2343
     2344    c_recpho_eff_pt->Print(pdfOutput,"pdf");
     2345    c_recpho_eff_pt->Print(figPath+"img_recpho_eff_pt.pdf","pdf");
     2346    c_recpho_eff_pt->Print(figPath+"img_recpho_eff_pt.png","png");
     2347
     2348    TCanvas *c_recpho_eff_eta = new TCanvas("","", 800, 600);
     2349
     2350    mg_recpho_eff_eta->Draw("APE");
     2351    DrawAxis(mg_recpho_eff_eta, leg_recpho_eff_eta, etaMin, etaMax, 0.0, 100, " #eta ", "reconstruction efficiency (%)", false, false);
     2352    leg_recpho_eff_eta->Draw();
     2353
     2354    c_recpho_eff_eta->Print(pdfOutput,"pdf");
     2355    c_recpho_eff_eta->Print(figPath+"img_recpho_eff_eta.pdf","pdf");
     2356    c_recpho_eff_eta->Print(figPath+"img_recpho_eff_eta.png","png");
     2357
     2358    /////////////////////////////////////////
     2359    // B-jets  Efficiency/ mistag rates   ///
     2360    /////////////////////////////////////////
     2361
     2362    TMultiGraph *mg_recbjet_eff_pt  = new TMultiGraph("","");
     2363    TMultiGraph *mg_recbjet_eff_eta = new TMultiGraph("","");
     2364
     2365    TLegend *leg_recbjet_eff_pt  = new TLegend(0.50,0.22,0.90,0.48);
     2366    TLegend *leg_recbjet_eff_eta = new TLegend(0.50,0.22,0.90,0.48);
     2367
     2368    TGraphErrors gr_recbjet_eff_pt[n_etabins], gr_recbjet_eff_eta[n_ptbins];
     2369    TH1D* h_recbjet_eff_pt, *h_recbjet_eff_eta;
     2370
     2371    // loop over eta bins
     2372    for (k = 0; k < etaVals.size()-1; k++)
     2373    {
     2374
     2375       h_recbjet_eff_pt = GetEffPt<Jet>(branchPFBJet, branchParticleBJet, "BJet", 5, ptMin, ptMax, etaVals.at(k), etaVals.at(k+1), treeReaderBJet);
     2376       gr_recbjet_eff_pt[k] = TGraphErrors(h_recbjet_eff_pt);
     2377
     2378       s_etaMin = Form("%.1f",etaVals.at(k));
     2379       s_etaMax = Form("%.1f",etaVals.at(k+1));
     2380
     2381       s_eta = "b-jet , " + s_etaMin + " < | #eta | < "+s_etaMax;
     2382
     2383       gr_recbjet_eff_pt[k].SetName("recEff_"+s_etaMin+"_"+s_etaMax);
     2384
     2385       addResoGraph(mg_recbjet_eff_pt, &gr_recbjet_eff_pt[k], leg_recbjet_eff_pt, markerStyles.at(k), colors.at(k), s_eta);
     2386    }
     2387
     2388    // loop over pt
     2389    for (k = 0; k < ptVals.size(); k++)
     2390    {
     2391       h_recbjet_eff_eta = GetEffEta<Jet>(branchPFBJet, branchParticleBJet, "BJet", 5, 0.5*ptVals.at(k), 2.0*ptVals.at(k) ,etaMin, etaMax , treeReaderBJet);
     2392       gr_recbjet_eff_eta[k] = TGraphErrors(h_recbjet_eff_eta);
     2393
     2394       s_pt = Form("b-jet , p_{T} = %.0f GeV",ptVals.at(k));
     2395       if(ptVals.at(k) >= 1000.) s_pt = Form("b-jet , p_{T} = %.0f TeV",ptVals.at(k)/1000.);
     2396
     2397       addResoGraph(mg_recbjet_eff_eta, &gr_recbjet_eff_eta[k], leg_recbjet_eff_eta, markerStyles.at(k), colors.at(k), s_pt );
     2398    }
     2399
     2400    TCanvas *c_recbjet_eff_pt = new TCanvas("","", 800, 600);
     2401
     2402    mg_recbjet_eff_pt->Draw("APE");
     2403    DrawAxis(mg_recbjet_eff_pt, leg_recbjet_eff_pt, ptMin, ptMax, 0.0, 100, "p_{T} [GeV]", "b - tag efficiency (%)", true, false);
     2404    leg_recbjet_eff_pt->Draw();
     2405
     2406    c_recbjet_eff_pt->Print(pdfOutput,"pdf");
     2407    c_recbjet_eff_pt->Print(figPath+"img_recbjet_eff_pt.pdf","pdf");
     2408    c_recbjet_eff_pt->Print(figPath+"img_recbjet_eff_pt.png","png");
     2409
     2410    TCanvas *c_recbjet_eff_eta = new TCanvas("","", 800, 600);
     2411
     2412    mg_recbjet_eff_eta->Draw("APE");
     2413    DrawAxis(mg_recbjet_eff_eta, leg_recbjet_eff_eta, etaMin, etaMax, 0.0, 100, " #eta ", "b - tag efficiency (%)", false, false);
     2414    leg_recbjet_eff_eta->Draw();
     2415
     2416    c_recbjet_eff_eta->Print(pdfOutput,"pdf");
     2417    c_recbjet_eff_eta->Print(figPath+"img_recbjet_eff_eta.pdf","pdf");
     2418    c_recbjet_eff_eta->Print(figPath+"img_recbjet_eff_eta.png","png");
     2419
     2420
     2421
     2422    // ------ c - mistag  ------
     2423
     2424    TMultiGraph *mg_recbjet_cmis_pt  = new TMultiGraph("","");
     2425    TMultiGraph *mg_recbjet_cmis_eta = new TMultiGraph("","");
     2426
     2427    TLegend *leg_recbjet_cmis_pt  = new TLegend(0.50,0.64,0.90,0.90);
     2428    TLegend *leg_recbjet_cmis_eta = new TLegend(0.50,0.64,0.90,0.90);
     2429
     2430    TGraphErrors gr_recbjet_cmis_pt[n_etabins], gr_recbjet_cmis_eta[n_ptbins];
     2431    TH1D* h_recbjet_cmis_pt, *h_recbjet_cmis_eta;
     2432
     2433    // loop over eta bins
     2434    for (k = 0; k < etaVals.size()-1; k++)
     2435    {
     2436
     2437       h_recbjet_cmis_pt = GetEffPt<Jet>(branchPFCJet, branchParticleCJet, "CJet", 4, ptMin, ptMax, etaVals.at(k), etaVals.at(k+1), treeReaderCJet);
     2438       gr_recbjet_cmis_pt[k] = TGraphErrors(h_recbjet_cmis_pt);
     2439
     2440       s_etaMin = Form("%.1f",etaVals.at(k));
     2441       s_etaMax = Form("%.1f",etaVals.at(k+1));
     2442
     2443       s_eta = "c-jet , " + s_etaMin + " < | #eta | < "+s_etaMax;
     2444
     2445       gr_recbjet_cmis_pt[k].SetName("recEff_"+s_etaMin+"_"+s_etaMax);
     2446
     2447       addResoGraph(mg_recbjet_cmis_pt, &gr_recbjet_cmis_pt[k], leg_recbjet_cmis_pt, markerStyles.at(k), colors.at(k), s_eta);
     2448    }
     2449
     2450    // loop over pt
     2451    for (k = 0; k < ptVals.size(); k++)
     2452    {
     2453       h_recbjet_cmis_eta = GetEffEta<Jet>(branchPFCJet, branchParticleCJet, "CJet", 4, 0.5*ptVals.at(k), 2.0*ptVals.at(k) ,etaMin, etaMax , treeReaderCJet);
     2454       gr_recbjet_cmis_eta[k] = TGraphErrors(h_recbjet_cmis_eta);
     2455
     2456       s_pt = Form("c-jet , p_{T} = %.0f GeV",ptVals.at(k));
     2457       if(ptVals.at(k) >= 1000.) s_pt = Form("c-jet , p_{T} = %.0f TeV",ptVals.at(k)/1000.);
     2458
     2459       addResoGraph(mg_recbjet_cmis_eta, &gr_recbjet_cmis_eta[k], leg_recbjet_cmis_eta, markerStyles.at(k), colors.at(k), s_pt );
     2460    }
     2461
     2462    TCanvas *c_recbjet_cmis_pt = new TCanvas("","", 800, 600);
     2463
     2464    mg_recbjet_cmis_pt->Draw("APE");
     2465    DrawAxis(mg_recbjet_cmis_pt, leg_recbjet_cmis_pt, ptMin, ptMax, 0.0, 20, "p_{T} [GeV]", "c - mistag rate (%)", true, false);
     2466    leg_recbjet_cmis_pt->Draw();
     2467
     2468    c_recbjet_cmis_pt->Print(pdfOutput,"pdf");
     2469    c_recbjet_cmis_pt->Print(figPath+"img_recbjet_cmis_pt.pdf","pdf");
     2470    c_recbjet_cmis_pt->Print(figPath+"img_recbjet_cmis_pt.png","png");
     2471
     2472    TCanvas *c_recbjet_cmis_eta = new TCanvas("","", 800, 600);
     2473
     2474    mg_recbjet_cmis_eta->Draw("APE");
     2475    DrawAxis(mg_recbjet_cmis_eta, leg_recbjet_cmis_eta, etaMin, etaMax, 0.0, 20, " #eta ", "c - mistag rate (%)", false, false);
     2476    leg_recbjet_cmis_eta->Draw();
     2477
     2478    c_recbjet_cmis_eta->Print(pdfOutput,"pdf");
     2479    c_recbjet_cmis_eta->Print(figPath+"img_recbjet_cmis_eta.pdf","pdf");
     2480    c_recbjet_cmis_eta->Print(figPath+"img_recbjet_cmis_eta.png","png");
     2481
     2482    // ------ light - mistag  ------
     2483
     2484    TMultiGraph *mg_recbjet_lmis_pt  = new TMultiGraph("","");
     2485    TMultiGraph *mg_recbjet_lmis_eta = new TMultiGraph("","");
     2486
     2487    TLegend *leg_recbjet_lmis_pt  = new TLegend(0.50,0.64,0.90,0.90);
     2488    TLegend *leg_recbjet_lmis_eta = new TLegend(0.50,0.64,0.90,0.90);
     2489
     2490    TGraphErrors gr_recbjet_lmis_pt[n_etabins], gr_recbjet_lmis_eta[n_ptbins];
     2491    TH1D* h_recbjet_lmis_pt, *h_recbjet_lmis_eta;
     2492
     2493    // loop over eta bins
     2494    for (k = 0; k < etaVals.size()-1; k++)
     2495    {
     2496
     2497       h_recbjet_lmis_pt = GetEffPt<Jet>(branchPFJet, branchParticleJet, "Jet", 1, ptMin, ptMax, etaVals.at(k), etaVals.at(k+1), treeReaderJet);
     2498       gr_recbjet_lmis_pt[k] = TGraphErrors(h_recbjet_lmis_pt);
     2499
     2500       s_etaMin = Form("%.1f",etaVals.at(k));
     2501       s_etaMax = Form("%.1f",etaVals.at(k+1));
     2502
     2503       s_eta = "uds-jet , " + s_etaMin + " < | #eta | < "+s_etaMax;
     2504
     2505       gr_recbjet_lmis_pt[k].SetName("recEff_"+s_etaMin+"_"+s_etaMax);
     2506
     2507       addResoGraph(mg_recbjet_lmis_pt, &gr_recbjet_lmis_pt[k], leg_recbjet_lmis_pt, markerStyles.at(k), colors.at(k), s_eta);
     2508    }
     2509
     2510    // loop over pt
     2511    for (k = 0; k < ptVals.size(); k++)
     2512    {
     2513       h_recbjet_lmis_eta = GetEffEta<Jet>(branchPFJet, branchParticleJet, "Jet", 1, 0.5*ptVals.at(k), 2.0*ptVals.at(k) ,etaMin, etaMax , treeReaderJet);
     2514       gr_recbjet_lmis_eta[k] = TGraphErrors(h_recbjet_lmis_eta);
     2515
     2516       s_pt = Form("uds-jet , p_{T} = %.0f GeV",ptVals.at(k));
     2517       if(ptVals.at(k) >= 1000.) s_pt = Form("uds-jet , p_{T} = %.0f TeV",ptVals.at(k)/1000.);
     2518
     2519       addResoGraph(mg_recbjet_lmis_eta, &gr_recbjet_lmis_eta[k], leg_recbjet_lmis_eta, markerStyles.at(k), colors.at(k), s_pt );
     2520    }
     2521
     2522    TCanvas *c_recbjet_lmis_pt = new TCanvas("","", 800, 600);
     2523
     2524    mg_recbjet_lmis_pt->Draw("APE");
     2525    DrawAxis(mg_recbjet_lmis_pt, leg_recbjet_lmis_pt, ptMin, ptMax, 0.0, 0.1, "p_{T} [GeV]", "light - mistag rate (%)", true, false);
     2526    leg_recbjet_lmis_pt->Draw();
     2527
     2528    c_recbjet_lmis_pt->Print(pdfOutput,"pdf");
     2529    c_recbjet_lmis_pt->Print(figPath+"img_recbjet_lmis_pt.pdf","pdf");
     2530    c_recbjet_lmis_pt->Print(figPath+"img_recbjet_lmis_pt.png","png");
     2531
     2532    TCanvas *c_recbjet_lmis_eta = new TCanvas("","", 800, 600);
     2533
     2534    mg_recbjet_lmis_eta->Draw("APE");
     2535    DrawAxis(mg_recbjet_lmis_eta, leg_recbjet_lmis_eta, etaMin, etaMax, 0.0, 0.1, " #eta ", "light - mistag rate (%)", false, false);
     2536    leg_recbjet_lmis_eta->Draw();
     2537
     2538    c_recbjet_lmis_eta->Print(pdfOutput,"pdf");
     2539    c_recbjet_lmis_eta->Print(figPath+"img_recbjet_lmis_eta.pdf","pdf");
     2540    c_recbjet_lmis_eta->Print(figPath+"img_recbjet_lmis_eta.png","png");
     2541
     2542
     2543    ///////////////////////////////////////////
     2544    // tau-jets  Efficiency/ mistag rates   ///
     2545    ///////////////////////////////////////////
     2546
     2547    TMultiGraph *mg_rectaujet_eff_pt  = new TMultiGraph("","");
     2548    TMultiGraph *mg_rectaujet_eff_eta = new TMultiGraph("","");
     2549
     2550    TLegend *leg_rectaujet_eff_pt  = new TLegend(0.50,0.22,0.90,0.48);
     2551    TLegend *leg_rectaujet_eff_eta = new TLegend(0.50,0.22,0.90,0.48);
     2552
     2553    TGraphErrors gr_rectaujet_eff_pt[n_etabins], gr_rectaujet_eff_eta[n_ptbins];
     2554    TH1D* h_rectaujet_eff_pt, *h_rectaujet_eff_eta;
     2555
     2556    // loop over eta bins
     2557    for (k = 0; k < etaVals.size()-1; k++)
     2558    {
     2559
     2560       h_rectaujet_eff_pt = GetTauEffPt<Jet>(branchPFTauJet, branchParticleTauJet, "TauJet", 15, ptMin, ptMax, etaVals.at(k), etaVals.at(k+1), treeReaderTauJet);
     2561       gr_rectaujet_eff_pt[k] = TGraphErrors(h_rectaujet_eff_pt);
     2562
     2563       s_etaMin = Form("%.1f",etaVals.at(k));
     2564       s_etaMax = Form("%.1f",etaVals.at(k+1));
     2565
     2566       s_eta = "#tau-jet , " + s_etaMin + " < | #eta | < "+s_etaMax;
     2567
     2568       gr_rectaujet_eff_pt[k].SetName("recEff_"+s_etaMin+"_"+s_etaMax);
     2569
     2570       addResoGraph(mg_rectaujet_eff_pt, &gr_rectaujet_eff_pt[k], leg_rectaujet_eff_pt, markerStyles.at(k), colors.at(k), s_eta);
     2571    }
     2572
     2573    // loop over pt
     2574    for (k = 0; k < ptVals.size(); k++)
     2575    {
     2576       h_rectaujet_eff_eta = GetTauEffEta<Jet>(branchPFTauJet, branchParticleTauJet, "TauJet", 15, 0.5*ptVals.at(k), 2.0*ptVals.at(k) ,etaMin, etaMax , treeReaderTauJet);
     2577       gr_rectaujet_eff_eta[k] = TGraphErrors(h_rectaujet_eff_eta);
     2578
     2579       s_pt = Form("#tau-jet , p_{T} = %.0f GeV",ptVals.at(k));
     2580       if(ptVals.at(k) >= 1000.) s_pt = Form("#tau-jet , p_{T} = %.0f TeV",ptVals.at(k)/1000.);
     2581
     2582       addResoGraph(mg_rectaujet_eff_eta, &gr_rectaujet_eff_eta[k], leg_rectaujet_eff_eta, markerStyles.at(k), colors.at(k), s_pt );
     2583    }
     2584
     2585
     2586    TCanvas *c_rectaujet_eff_pt = new TCanvas("","", 800, 600);
     2587
     2588    mg_rectaujet_eff_pt->Draw("APE");
     2589    DrawAxis(mg_rectaujet_eff_pt, leg_rectaujet_eff_pt, ptMin, ptMax, 0.0, 100, "p_{T} [GeV]", "#tau - tag efficiency (%)", true, false);
     2590    leg_rectaujet_eff_pt->Draw();
     2591
     2592    c_rectaujet_eff_pt->Print(pdfOutput,"pdf");
     2593    c_rectaujet_eff_pt->Print(figPath+"img_rectaujet_eff_pt.pdf","pdf");
     2594    c_rectaujet_eff_pt->Print(figPath+"img_rectaujet_eff_pt.png","png");
     2595
     2596    TCanvas *c_rectaujet_eff_eta = new TCanvas("","", 800, 600);
     2597
     2598    mg_rectaujet_eff_eta->Draw("APE");
     2599    DrawAxis(mg_rectaujet_eff_eta, leg_rectaujet_eff_eta, etaMin, etaMax, 0.0000001, 100, " #eta ", "#tau - tag efficiency (%)", false, true);
     2600    leg_rectaujet_eff_eta->Draw();
     2601
     2602    c_rectaujet_eff_eta->Print(pdfOutput,"pdf");
     2603    c_rectaujet_eff_eta->Print(figPath+"img_rectaujet_eff_eta.pdf","pdf");
     2604    c_rectaujet_eff_eta->Print(figPath+"img_rectaujet_eff_eta.png","png");
     2605
     2606
     2607    //--------------- tau mistag rate ----------
     2608
     2609    TMultiGraph *mg_rectaujet_mis_pt  = new TMultiGraph("","");
     2610    TMultiGraph *mg_rectaujet_mis_eta = new TMultiGraph("","");
     2611
     2612    TLegend *leg_rectaujet_mis_pt  = new TLegend(0.50,0.64,0.90,0.90);
     2613    TLegend *leg_rectaujet_mis_eta = new TLegend(0.50,0.64,0.90,0.90);
     2614
     2615    TGraphErrors gr_rectaujet_mis_pt[n_etabins], gr_rectaujet_mis_eta[n_ptbins];
     2616    TH1D* h_rectaujet_mis_pt, *h_rectaujet_mis_eta;
     2617
     2618    // loop over eta bins
     2619    for (k = 0; k < etaVals.size()-1; k++)
     2620    {
     2621
     2622       h_rectaujet_mis_pt = GetTauEffPt<Jet>(branchPFJet, branchParticleJet, "TauJet", 1, ptMin, ptMax, etaVals.at(k), etaVals.at(k+1), treeReaderJet);
     2623       gr_rectaujet_mis_pt[k] = TGraphErrors(h_rectaujet_mis_pt);
     2624
     2625       s_etaMin = Form("%.1f",etaVals.at(k));
     2626       s_etaMax = Form("%.1f",etaVals.at(k+1));
     2627
     2628       s_eta = "uds-jet , " + s_etaMin + " < | #eta | < "+s_etaMax;
     2629
     2630       gr_rectaujet_mis_pt[k].SetName("recEff_"+s_etaMin+"_"+s_etaMax);
     2631
     2632       addResoGraph(mg_rectaujet_mis_pt, &gr_rectaujet_mis_pt[k], leg_rectaujet_mis_pt, markerStyles.at(k), colors.at(k), s_eta);
     2633    }
     2634
     2635    // loop over pt
     2636    for (k = 0; k < ptVals.size(); k++)
     2637    {
     2638       h_rectaujet_mis_eta = GetTauEffEta<Jet>(branchPFJet, branchParticleJet, "TauJet", 1, 0.5*ptVals.at(k), 2.0*ptVals.at(k) ,etaMin, etaMax , treeReaderJet);
     2639       gr_rectaujet_mis_eta[k] = TGraphErrors(h_rectaujet_mis_eta);
     2640
     2641       s_pt = Form("uds-jet , p_{T} = %.0f GeV",ptVals.at(k));
     2642       if(ptVals.at(k) >= 1000.) s_pt = Form("uds-jet , p_{T} = %.0f TeV",ptVals.at(k)/1000.);
     2643
     2644       addResoGraph(mg_rectaujet_mis_eta, &gr_rectaujet_mis_eta[k], leg_rectaujet_mis_eta, markerStyles.at(k), colors.at(k), s_pt );
     2645    }
     2646
     2647    TCanvas *c_rectaujet_mis_pt = new TCanvas("","", 800, 600);
     2648
     2649    mg_rectaujet_mis_pt->Draw("APE");
     2650    DrawAxis(mg_rectaujet_mis_pt, leg_rectaujet_mis_pt, ptMin, ptMax, 0.0, 10., "p_{T} [GeV]", "#tau - mistag(%)", true, false);
     2651    leg_rectaujet_mis_pt->Draw();
     2652
     2653    c_rectaujet_mis_pt->Print(pdfOutput,"pdf");
     2654    c_rectaujet_mis_pt->Print(figPath+"img_rectaujet_mis_pt.pdf","pdf");
     2655    c_rectaujet_mis_pt->Print(figPath+"img_rectaujet_mis_pt.png","png");
     2656
     2657    TCanvas *c_rectaujet_mis_eta = new TCanvas("","", 800, 600);
     2658
     2659    mg_rectaujet_mis_eta->Draw("APE");
     2660    DrawAxis(mg_rectaujet_mis_eta, leg_rectaujet_mis_eta, etaMin, etaMax, 0.0, 5., " #eta ", "#tau - mistag (%)", false, false);
     2661    leg_rectaujet_mis_eta->Draw();
     2662
     2663    c_rectaujet_mis_eta->Print(pdfOutput+")","pdf");
     2664    c_rectaujet_mis_eta->Print(figPath+"img_rectaujet_mis_eta.pdf","pdf");
     2665    c_rectaujet_mis_eta->Print(figPath+"img_rectaujet_mis_eta.png","png");
     2666
     2667
     2668//   ----   store resolution histograms in the output (for leave efficiencies out) ---
     2669
    12762670
    12772671  TFile *fout = new TFile(outputFile,"recreate");
    … …  
    12792673  for (int bin = 0; bin < Nbins; bin++)
    12802674  {
    1281     plots_el.at(bin).cenResolHist->Write();
    1282     plots_eltrack.at(bin).cenResolHist->Write();
    1283     plots_eltower.at(bin).cenResolHist->Write();
    1284     plots_el.at(bin).fwdResolHist->Write();
    1285     plots_eltrack.at(bin).fwdResolHist->Write();
    1286     plots_eltower.at(bin).fwdResolHist->Write();
    1287   }
    1288 
    1289   histos_el.first->Write();
    1290   histos_el.second->Write();
    1291   histos_eltrack.first->Write();
    1292   histos_eltrack.second->Write();
    1293   histos_eltower.first->Write();
    1294   histos_eltower.second->Write();
    1295   C_el1->Write();
    1296   C_el2->Write();
    1297 
    1298   for (int bin = 0; bin < Nbins; bin++)
    1299   {
    1300     plots_mu.at(bin).cenResolHist->Write();
    1301     plots_mutrack.at(bin).cenResolHist->Write();
    1302     plots_mu.at(bin).fwdResolHist->Write();
    1303     plots_mutrack.at(bin).fwdResolHist->Write();
    1304   }
    1305 
    1306   histos_mu.first->Write();
    1307   histos_mu.second->Write();
    1308   histos_mutrack.first->Write();
    1309   histos_mutrack.second->Write();
    1310   C_mu1->Write();
    1311   C_mu->Write();
    1312 
    1313   histos_ph.first->Write();
    1314   histos_ph.second->Write();
    1315   C_ph1->Write();
    1316   C_ph->Write();
    1317 
    1318   for (int bin = 0; bin < Nbins; bin++)
    1319   {
    1320     plots_pfjets.at(bin).cenResolHist->Write();
    1321     plots_pfjets.at(bin).fwdResolHist->Write();
    1322     plots_calojets.at(bin).cenResolHist->Write();
    1323     plots_calojets.at(bin).fwdResolHist->Write();
    1324     plots_met.at(bin).cenResolHist->Write();
    1325   }
    1326   histos_btag.first->Write();
    1327   histos_btag.second->Write();
    1328   histos_tautag.first->Write();
    1329   histos_tautag.second->Write();
    1330   C_btag1->Write();
    1331   C_tautag1->Write();
    1332   C_jet->Write();
    1333   C_met->Write();
     2675
     2676    for (k = 0; k < etaVals.size()-1; k++)
     2677    {
     2678       plots_trkpi_res_pt[k].at(bin).resolHist->Write();
     2679       plots_trkele_res_pt[k].at(bin).resolHist->Write();
     2680           plots_trkmu_res_pt[k].at(bin).resolHist->Write();
     2681       plots_ecal_res_e[k].at(bin).resolHist->Write();
     2682       plots_hcal_res_e[k].at(bin).resolHist->Write();
     2683       plots_pfele_res_e[k].at(bin).resolHist->Write();
     2684       plots_pfpi_res_e[k].at(bin).resolHist->Write();
     2685       plots_pfjet_res_e[k].at(bin).resolHist->Write();
     2686       plots_cajet_res_e[k].at(bin).resolHist->Write();
     2687
     2688    }
     2689
     2690    for (k = 0; k < ptVals.size(); k++)
     2691    {
     2692       plots_trkpi_res_eta[k].at(bin).resolHist->Write();
     2693       plots_trkele_res_eta[k].at(bin).resolHist->Write();
     2694       plots_trkmu_res_eta[k].at(bin).resolHist->Write();
     2695       plots_ecal_res_eta[k].at(bin).resolHist->Write();
     2696       plots_hcal_res_eta[k].at(bin).resolHist->Write();
     2697       plots_pfele_res_eta[k].at(bin).resolHist->Write();
     2698       plots_pfpi_res_eta[k].at(bin).resolHist->Write();
     2699       plots_pfjet_res_eta[k].at(bin).resolHist->Write();
     2700       plots_cajet_res_eta[k].at(bin).resolHist->Write();
     2701
     2702    }
     2703
     2704    plots_pfmet.at(bin).resolHist->Write();
     2705    plots_camet.at(bin).resolHist->Write();
     2706
     2707  }
    13342708
    13352709  fout->Write();
    … …  
    13572731  char *appName = "Validation";
    13582732
    1359   if(argc != 8)
     2733  if(argc != 11)
    13602734  {
    13612735    cout << " Usage: " << appName << " input_file_electron input_file_muon input_file_photon input_file_jet input_file_bjet input_file_taujet output_file" << endl;
     2736    cout << " input_file_pion  - input file in ROOT format ('Delphes' tree)," << endl;
    13622737    cout << " input_file_electron  - input file in ROOT format ('Delphes' tree)," << endl;
    13632738    cout << " input_file_muon - input file in ROOT format ('Delphes' tree)," << endl;
    13642739    cout << " input_file_photon - input file in ROOT format ('Delphes' tree)," << endl;
     2740    cout << " input_file_neutralhadron - input file in ROOT format ('Delphes' tree)," << endl;
    13652741    cout << " input_file_jet - input file in ROOT format ('Delphes' tree)," << endl;
    13662742    cout << " input_file_bjet - input file in ROOT format ('Delphes' tree)," << endl;
     2743    cout << " input_file_cjet - input file in ROOT format ('Delphes' tree)," << endl;
    13672744    cout << " input_file_taujet - input file in ROOT format ('Delphes' tree)," << endl;
    13682745    cout << " output_file - output file in ROOT format" << endl;
    … …  
    13762753  TApplication app(appName, &appargc, appargv);
    13772754
    1378   Validation(argv[1], argv[2], argv[3], argv[4], argv[5], argv[6], argv[7]);
     2755  Validation(argv[1], argv[2], argv[3], argv[4], argv[5], argv[6], argv[7], argv[8], argv[9], argv[10]);
    13792756}
    13802757
    13812758
     2759
Note: See TracChangeset for help on using the changeset viewer.