Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-11-15 23:40:49

0001 // --------------------------------------------------------------------------------------------------------
0002 // ROOT script for making tracking performance plots using the ntuples produced by L1TrackNtupleMaker.cc
0003 //
0004 // e.g. in ROOT do: .L L1TrackNtuplePlot.C++, L1TrackNtuplePlot("L1TrkNtuple")
0005 //
0006 // By Louise Skinnari, June 2013
0007 // --------------------------------------------------------------------------------------------------------
0008 
0009 #include "TROOT.h"
0010 #include "TStyle.h"
0011 #include "TLatex.h"
0012 #include "TFile.h"
0013 #include "TTree.h"
0014 #include "TChain.h"
0015 #include "TBranch.h"
0016 #include "TLeaf.h"
0017 #include "TCanvas.h"
0018 #include "TLegend.h"
0019 #include "TH1.h"
0020 #include "TH2.h"
0021 #include "TF1.h"
0022 #include "TProfile.h"
0023 #include "TProfile2D.h"
0024 #include "TMath.h"
0025 #include <TError.h>
0026 #include "TSystem.h"
0027 
0028 #include <iostream>
0029 #include <string>
0030 #include <vector>
0031 #include <algorithm>
0032 
0033 using namespace std;
0034 
0035 void SetPlotStyle();
0036 void mySmallText(Double_t x, Double_t y, Color_t color, char* text);
0037 
0038 double getIntervalContainingFractionOfEntries(TH1* histogram, double interval, int minEntries = 5);
0039 void makeResidualIntervalPlot(
0040     TString type, TString dir, TString variable, TH1F* h_68, TH1F* h_90, TH1F* h_99, double minY, double maxY);
0041 
0042 // ----------------------------------------------------------------------------------------------------------------
0043 // Main script
0044 // ----------------------------------------------------------------------------------------------------------------
0045 
0046 void L1TrackNtuplePlot(TString type,
0047                        TString type_dir = "",
0048                        TString treeName = "",
0049                        int TP_select_injet = 0,
0050                        int TP_select_pdgid = 0,
0051                        int TP_select_eventid = 0,
0052                        bool useTightCuts = false,
0053                        bool useDeadRegion = false,
0054                        float TP_minPt = 2.0,
0055                        float TP_maxPt = 100.0,
0056                        float TP_maxEta = 2.4,
0057                        float TP_maxDxy = 1.0,
0058                        float TP_maxD0 = 1.0,
0059                        bool doDetailedPlots = false) {
0060   // type:              this is the name of the input file you want to process (minus ".root" extension)
0061   // type_dir:          this is the directory containing the input file you want to process. Note that this must end with a "/", as in "EventSets/"
0062   // TP_select_pdgid:   if non-zero, only select TPs with a given PDG ID
0063   // TP_select_eventid: if zero, only look at TPs from primary interaction, else, include TPs from pileup
0064   // TP_minPt:          only look at TPs with pt > X GeV
0065   // TP_maxPt:          only look at TPs with pt < X GeV
0066   // TP_maxEta:         only look at TPs with |eta| < X
0067   // doDetailedPlots:   includes extra plots, such as  performance vs d0.
0068 
0069   // TP_select_injet: only look at TPs that are within a jet with pt > 30 GeV (==1) or within a jet with pt > 100 GeV (==2), >200 GeV (==3) or all TPs (==0)
0070 
0071   //--  N.B. For standard displaced tracking plots, set TP_minPt=3.0, TP_maxEta=2.0, TP_maxDxy=10.0,
0072   //--  TO_maxD0=10.0, doDetailedPlots=true. (Efficiency plots vs eta also usually made for d0 < 5).
0073 
0074   gROOT->SetBatch();
0075   gErrorIgnoreLevel = kWarning;
0076 
0077   SetPlotStyle();
0078 
0079   cout.setf(ios::fixed);
0080   cout.precision(2);
0081 
0082   // ----------------------------------------------------------------------------------------------------------------
0083   // define input options
0084 
0085   // these are the LOOSE cuts, baseline scenario for efficiency and rate plots ==> configure as appropriate
0086   int L1Tk_minNstub = 4;
0087   float L1Tk_maxChi2 = 999999;
0088   float L1Tk_maxChi2dof = 999999.;
0089 
0090   // TIGHT cuts (separate plots / rates) ==> configure as appropriate
0091   // this is currently set up as an either or for performance plots, to not duplicate a ton of code.
0092   int L1Tk_TIGHT_minNstub = 4;
0093   float L1Tk_TIGHT_maxChi2 = 999999;
0094   float L1Tk_TIGHT_maxChi2dof = 999999.;
0095   if (useTightCuts) {
0096     L1Tk_minNstub = L1Tk_TIGHT_minNstub;
0097     L1Tk_maxChi2 = L1Tk_TIGHT_maxChi2;
0098     L1Tk_maxChi2dof = L1Tk_TIGHT_maxChi2dof;
0099   }
0100 
0101   bool doGausFit = false;     //do gaussian fit for resolution vs eta/pt plots
0102   bool doLooseMatch = false;  //looser MC truth matching
0103 
0104   // tracklet variables
0105   int L1Tk_seed = 0;
0106 
0107   //some counters for integrated efficiencies
0108   int n_all_eta2p5 = 0;
0109   int n_all_eta1p75 = 0;
0110   int n_all_eta1p0 = 0;
0111   int n_match_eta2p5 = 0;
0112   int n_match_eta1p75 = 0;
0113   int n_match_eta1p0 = 0;
0114   int n_all_ptg2 = 0;
0115   int n_all_ptg8 = 0;
0116   int n_all_pt2to8 = 0;
0117   int n_all_ptg40 = 0;
0118   int n_match_ptg2 = 0;
0119   int n_match_ptg8 = 0;
0120   int n_match_pt2to8 = 0;
0121   int n_match_ptg40 = 0;
0122 
0123   // counters for total track rates
0124   int ntrk = 0;
0125   int ntrk_pt2 = 0;
0126   int ntrk_pt3 = 0;
0127   int ntrk_pt10 = 0;
0128   int ntp_pt2 = 0;
0129   int ntp_pt3 = 0;
0130   int ntp_pt10 = 0;
0131   int ntrk_genuine_pt2 = 0;
0132   int ntp_nmatch = 0;
0133   // ----------------------------------------------------------------------------------------------------------------
0134   // read ntuples
0135   TChain* tree = new TChain("L1TrackNtuple" + treeName + "/eventTree");
0136   tree->Add(type_dir + type + ".root");
0137 
0138   if (tree->GetEntries() == 0) {
0139     cout << "File doesn't exist or is empty, returning..."
0140          << endl;  //cout's kept in this file as it is an example standalone plotting script, not running in central CMSSW
0141     return;
0142   }
0143 
0144   // ----------------------------------------------------------------------------------------------------------------
0145   // define leafs & branches
0146 
0147   // tracking particles
0148   vector<float>* tp_pt;
0149   vector<float>* tp_eta;
0150   vector<float>* tp_phi;
0151   vector<float>* tp_dxy;
0152   vector<float>* tp_z0;
0153   vector<float>* tp_d0;
0154   vector<int>* tp_pdgid;
0155   vector<int>* tp_nmatch;
0156   vector<int>* tp_nstub;
0157   vector<int>* tp_eventid;
0158   vector<int>* tp_injet;
0159   vector<int>* tp_injet_highpt;
0160   vector<int>* tp_injet_vhighpt;
0161 
0162   // *L1 track* properties, for tracking particles matched to a L1 track
0163   vector<float>* matchtrk_pt;
0164   vector<float>* matchtrk_eta;
0165   vector<float>* matchtrk_phi;
0166   vector<float>* matchtrk_d0;
0167   vector<float>* matchtrk_z0;
0168   vector<float>* matchtrk_chi2;
0169   vector<float>* matchtrk_chi2_dof;
0170   vector<float>* matchtrk_chi2rphi;
0171   vector<float>* matchtrk_chi2rphi_dof;
0172   vector<float>* matchtrk_chi2rz;
0173   vector<float>* matchtrk_chi2rz_dof;
0174   vector<int>* matchtrk_nstub;
0175   vector<int>* matchtrk_lhits;
0176   vector<int>* matchtrk_dhits;
0177   vector<int>* matchtrk_seed;
0178   vector<int>* matchtrk_hitpattern;
0179   vector<int>* matchtrk_injet;
0180   vector<int>* matchtrk_injet_highpt;
0181   vector<int>* matchtrk_injet_vhighpt;
0182 
0183   // all L1 tracks
0184   vector<float>* trk_pt;
0185   vector<float>* trk_eta;
0186   vector<float>* trk_phi;
0187   vector<float>* trk_chi2;
0188   vector<float>* trk_chi2_dof;
0189   vector<float>* trk_chi2rphi;
0190   vector<float>* trk_chi2rphi_dof;
0191   vector<float>* trk_chi2rz;
0192   vector<float>* trk_chi2rz_dof;
0193   vector<int>* trk_nstub;
0194   vector<int>* trk_lhits;
0195   vector<int>* trk_dhits;
0196   vector<int>* trk_seed;
0197   vector<int>* trk_hitpattern;
0198   vector<unsigned int>* trk_phiSector;
0199   vector<int>* trk_injet;
0200   vector<int>* trk_injet_highpt;
0201   vector<int>* trk_injet_vhighpt;
0202   vector<int>* trk_fake;
0203   vector<int>* trk_genuine;
0204   vector<int>* trk_loose;
0205 
0206   TBranch* b_tp_pt;
0207   TBranch* b_tp_eta;
0208   TBranch* b_tp_phi;
0209   TBranch* b_tp_dxy;
0210   TBranch* b_tp_z0;
0211   TBranch* b_tp_d0;
0212   TBranch* b_tp_pdgid;
0213   TBranch* b_tp_nmatch;
0214   TBranch* b_tp_nstub;
0215   TBranch* b_tp_eventid;
0216   TBranch* b_tp_injet;
0217   TBranch* b_tp_injet_highpt;
0218   TBranch* b_tp_injet_vhighpt;
0219 
0220   TBranch* b_matchtrk_pt;
0221   TBranch* b_matchtrk_eta;
0222   TBranch* b_matchtrk_phi;
0223   TBranch* b_matchtrk_d0;
0224   TBranch* b_matchtrk_z0;
0225   TBranch* b_matchtrk_chi2;
0226   TBranch* b_matchtrk_chi2_dof;
0227   TBranch* b_matchtrk_chi2rphi;
0228   TBranch* b_matchtrk_chi2rphi_dof;
0229   TBranch* b_matchtrk_chi2rz;
0230   TBranch* b_matchtrk_chi2rz_dof;
0231   TBranch* b_matchtrk_nstub;
0232   TBranch* b_matchtrk_lhits;
0233   TBranch* b_matchtrk_dhits;
0234   TBranch* b_matchtrk_seed;
0235   TBranch* b_matchtrk_hitpattern;
0236   TBranch* b_matchtrk_injet;
0237   TBranch* b_matchtrk_injet_highpt;
0238   TBranch* b_matchtrk_injet_vhighpt;
0239 
0240   TBranch* b_trk_pt;
0241   TBranch* b_trk_eta;
0242   TBranch* b_trk_phi;
0243   TBranch* b_trk_chi2;
0244   TBranch* b_trk_chi2_dof;
0245   TBranch* b_trk_chi2rphi;
0246   TBranch* b_trk_chi2rphi_dof;
0247   TBranch* b_trk_chi2rz;
0248   TBranch* b_trk_chi2rz_dof;
0249   TBranch* b_trk_nstub;
0250   TBranch* b_trk_lhits;
0251   TBranch* b_trk_dhits;
0252   TBranch* b_trk_phiSector;
0253   TBranch* b_trk_seed;
0254   TBranch* b_trk_hitpattern;
0255   TBranch* b_trk_injet;
0256   TBranch* b_trk_injet_highpt;
0257   TBranch* b_trk_injet_vhighpt;
0258   TBranch* b_trk_fake;
0259   TBranch* b_trk_genuine;
0260   TBranch* b_trk_loose;
0261 
0262   tp_pt = 0;
0263   tp_eta = 0;
0264   tp_phi = 0;
0265   tp_dxy = 0;
0266   tp_z0 = 0;
0267   tp_d0 = 0;
0268   tp_pdgid = 0;
0269   tp_nmatch = 0;
0270   tp_nstub = 0;
0271   tp_eventid = 0;
0272   tp_injet = 0;
0273   tp_injet_highpt = 0;
0274   tp_injet_vhighpt = 0;
0275 
0276   matchtrk_pt = 0;
0277   matchtrk_eta = 0;
0278   matchtrk_phi = 0;
0279   matchtrk_d0 = 0;
0280   matchtrk_z0 = 0;
0281   matchtrk_chi2 = 0;
0282   matchtrk_chi2_dof = 0;
0283   matchtrk_chi2rphi = 0;
0284   matchtrk_chi2rphi_dof = 0;
0285   matchtrk_chi2rz = 0;
0286   matchtrk_chi2rz_dof = 0;
0287   matchtrk_nstub = 0;
0288   matchtrk_lhits = 0;
0289   matchtrk_dhits = 0;
0290   matchtrk_seed = 0;
0291   matchtrk_hitpattern = 0;
0292   matchtrk_injet = 0;
0293   matchtrk_injet_highpt = 0;
0294   matchtrk_injet_vhighpt = 0;
0295 
0296   trk_pt = 0;
0297   trk_eta = 0;
0298   trk_phi = 0;
0299   trk_chi2 = 0;
0300   trk_chi2_dof = 0;
0301   trk_chi2rphi = 0;
0302   trk_chi2rphi_dof = 0;
0303   trk_chi2rz = 0;
0304   trk_chi2rz_dof = 0;
0305   trk_nstub = 0;
0306   trk_lhits = 0;
0307   trk_dhits = 0;
0308   trk_phiSector = 0;
0309   trk_seed = 0;
0310   trk_hitpattern = 0;
0311   trk_injet = 0;
0312   trk_injet_highpt = 0;
0313   trk_injet_vhighpt = 0;
0314   trk_fake = 0;
0315   trk_genuine = 0;
0316   trk_loose = 0;
0317 
0318   tree->SetBranchAddress("tp_pt", &tp_pt, &b_tp_pt);
0319   tree->SetBranchAddress("tp_eta", &tp_eta, &b_tp_eta);
0320   tree->SetBranchAddress("tp_phi", &tp_phi, &b_tp_phi);
0321   tree->SetBranchAddress("tp_dxy", &tp_dxy, &b_tp_dxy);
0322   tree->SetBranchAddress("tp_z0", &tp_z0, &b_tp_z0);
0323   tree->SetBranchAddress("tp_d0", &tp_d0, &b_tp_d0);
0324   tree->SetBranchAddress("tp_pdgid", &tp_pdgid, &b_tp_pdgid);
0325   if (doLooseMatch)
0326     tree->SetBranchAddress("tp_nloosematch", &tp_nmatch, &b_tp_nmatch);
0327   else
0328     tree->SetBranchAddress("tp_nmatch", &tp_nmatch, &b_tp_nmatch);
0329   tree->SetBranchAddress("tp_nstub", &tp_nstub, &b_tp_nstub);
0330   tree->SetBranchAddress("tp_eventid", &tp_eventid, &b_tp_eventid);
0331   if (TP_select_injet > 0) {
0332     tree->SetBranchAddress("tp_injet", &tp_injet, &b_tp_injet);
0333     tree->SetBranchAddress("tp_injet_highpt", &tp_injet_highpt, &b_tp_injet_highpt);
0334     tree->SetBranchAddress("tp_injet_vhighpt", &tp_injet_vhighpt, &b_tp_injet_vhighpt);
0335   }
0336 
0337   if (doLooseMatch) {
0338     tree->SetBranchAddress("loosematchtrk_pt", &matchtrk_pt, &b_matchtrk_pt);
0339     tree->SetBranchAddress("loosematchtrk_eta", &matchtrk_eta, &b_matchtrk_eta);
0340     tree->SetBranchAddress("loosematchtrk_phi", &matchtrk_phi, &b_matchtrk_phi);
0341     tree->SetBranchAddress("loosematchtrk_d0", &matchtrk_d0, &b_matchtrk_d0);
0342     tree->SetBranchAddress("loosematchtrk_z0", &matchtrk_z0, &b_matchtrk_z0);
0343     tree->SetBranchAddress("loosematchtrk_chi2", &matchtrk_chi2, &b_matchtrk_chi2);
0344     tree->SetBranchAddress("loosematchtrk_chi2_dof", &matchtrk_chi2_dof, &b_matchtrk_chi2_dof);
0345     tree->SetBranchAddress("loosematchtrk_chi2rphi", &matchtrk_chi2rphi, &b_matchtrk_chi2rphi);
0346     tree->SetBranchAddress("loosematchtrk_chi2rphi_dof", &matchtrk_chi2rphi_dof, &b_matchtrk_chi2rphi_dof);
0347     tree->SetBranchAddress("loosematchtrk_chi2rz", &matchtrk_chi2rz, &b_matchtrk_chi2rz);
0348     tree->SetBranchAddress("loosematchtrk_chi2rz_dof", &matchtrk_chi2rz_dof, &b_matchtrk_chi2rz_dof);
0349     tree->SetBranchAddress("loosematchtrk_nstub", &matchtrk_nstub, &b_matchtrk_nstub);
0350     tree->SetBranchAddress("loosematchtrk_seed", &matchtrk_seed, &b_matchtrk_seed);
0351     tree->SetBranchAddress("loosematchtrk_hitpattern", &matchtrk_hitpattern, &b_matchtrk_hitpattern);
0352     if (TP_select_injet > 0) {
0353       tree->SetBranchAddress("loosematchtrk_injet", &matchtrk_injet, &b_matchtrk_injet);
0354       tree->SetBranchAddress("loosematchtrk_injet_highpt", &matchtrk_injet_highpt, &b_matchtrk_injet_highpt);
0355       tree->SetBranchAddress("loosematchtrk_injet_vhighpt", &matchtrk_injet_vhighpt, &b_matchtrk_injet_vhighpt);
0356     }
0357   } else {
0358     tree->SetBranchAddress("matchtrk_pt", &matchtrk_pt, &b_matchtrk_pt);
0359     tree->SetBranchAddress("matchtrk_eta", &matchtrk_eta, &b_matchtrk_eta);
0360     tree->SetBranchAddress("matchtrk_phi", &matchtrk_phi, &b_matchtrk_phi);
0361     tree->SetBranchAddress("matchtrk_d0", &matchtrk_d0, &b_matchtrk_d0);
0362     tree->SetBranchAddress("matchtrk_z0", &matchtrk_z0, &b_matchtrk_z0);
0363     tree->SetBranchAddress("matchtrk_chi2", &matchtrk_chi2, &b_matchtrk_chi2);
0364     tree->SetBranchAddress("matchtrk_chi2_dof", &matchtrk_chi2_dof, &b_matchtrk_chi2_dof);
0365     tree->SetBranchAddress("matchtrk_chi2rphi", &matchtrk_chi2rphi, &b_matchtrk_chi2rphi);
0366     tree->SetBranchAddress("matchtrk_chi2rphi_dof", &matchtrk_chi2rphi_dof, &b_matchtrk_chi2rphi_dof);
0367     tree->SetBranchAddress("matchtrk_chi2rz", &matchtrk_chi2rz, &b_matchtrk_chi2rz);
0368     tree->SetBranchAddress("matchtrk_chi2rz_dof", &matchtrk_chi2rz_dof, &b_matchtrk_chi2rz_dof);
0369     tree->SetBranchAddress("matchtrk_nstub", &matchtrk_nstub, &b_matchtrk_nstub);
0370     tree->SetBranchAddress("matchtrk_lhits", &matchtrk_lhits, &b_matchtrk_lhits);
0371     tree->SetBranchAddress("matchtrk_dhits", &matchtrk_dhits, &b_matchtrk_dhits);
0372     tree->SetBranchAddress("matchtrk_seed", &matchtrk_seed, &b_matchtrk_seed);
0373     tree->SetBranchAddress("matchtrk_hitpattern", &matchtrk_hitpattern, &b_matchtrk_hitpattern);
0374     if (TP_select_injet > 0) {
0375       tree->SetBranchAddress("matchtrk_injet", &matchtrk_injet, &b_matchtrk_injet);
0376       tree->SetBranchAddress("matchtrk_injet_highpt", &matchtrk_injet_highpt, &b_matchtrk_injet_highpt);
0377       tree->SetBranchAddress("matchtrk_injet_vhighpt", &matchtrk_injet_vhighpt, &b_matchtrk_injet_vhighpt);
0378     }
0379   }
0380 
0381   tree->SetBranchAddress("trk_pt", &trk_pt, &b_trk_pt);
0382   tree->SetBranchAddress("trk_eta", &trk_eta, &b_trk_eta);
0383   tree->SetBranchAddress("trk_phi", &trk_phi, &b_trk_phi);
0384   tree->SetBranchAddress("trk_chi2", &trk_chi2, &b_trk_chi2);
0385   tree->SetBranchAddress("trk_chi2_dof", &trk_chi2_dof, &b_trk_chi2_dof);
0386   tree->SetBranchAddress("trk_chi2rphi", &trk_chi2rphi, &b_trk_chi2rphi);
0387   tree->SetBranchAddress("trk_chi2rphi_dof", &trk_chi2rphi_dof, &b_trk_chi2rphi_dof);
0388   tree->SetBranchAddress("trk_chi2rz", &trk_chi2rz, &b_trk_chi2rz);
0389   tree->SetBranchAddress("trk_chi2rz_dof", &trk_chi2rz_dof, &b_trk_chi2rz_dof);
0390   tree->SetBranchAddress("trk_nstub", &trk_nstub, &b_trk_nstub);
0391   tree->SetBranchAddress("trk_lhits", &trk_lhits, &b_trk_lhits);
0392   tree->SetBranchAddress("trk_dhits", &trk_dhits, &b_trk_dhits);
0393   tree->SetBranchAddress("trk_phiSector", &trk_phiSector, &b_trk_phiSector);
0394   tree->SetBranchAddress("trk_seed", &trk_seed, &b_trk_seed);
0395   tree->SetBranchAddress("trk_hitpattern", &trk_hitpattern, &b_trk_hitpattern);
0396   tree->SetBranchAddress("trk_fake", &trk_fake, &b_trk_fake);
0397   tree->SetBranchAddress("trk_genuine", &trk_genuine, &b_trk_genuine);
0398   tree->SetBranchAddress("trk_loose", &trk_loose, &b_trk_loose);
0399   if (TP_select_injet > 0) {
0400     tree->SetBranchAddress("trk_injet", &trk_injet, &b_trk_injet);
0401     tree->SetBranchAddress("trk_injet_highpt", &trk_injet_highpt, &b_trk_injet_highpt);
0402     tree->SetBranchAddress("trk_injet_vhighpt", &trk_injet_vhighpt, &b_trk_injet_vhighpt);
0403   }
0404 
0405   // ----------------------------------------------------------------------------------------------------------------
0406   // histograms
0407   // ----------------------------------------------------------------------------------------------------------------
0408 
0409   /////////////////////////////////////////////////
0410   // NOTATION:                                   //
0411   // 'C' - Central eta range, |eta|<0.8          //
0412   // 'I' - Intermediate eta range, 0.8<|eta|<1.6 //
0413   // 'F' - Forward eta range, |eta|>1.6          //
0414   //                                             //
0415   // 'L' - Low pt range,  pt<8 GeV               //
0416   // 'H' - High pt range, pt>8 GeV               //
0417   /////////////////////////////////////////////////
0418 
0419   // ----------------------------------------------------------------------------------------------------------------
0420   // for efficiencies
0421 
0422   TH1F* h_tp_pt = new TH1F("tp_pt", ";Tracking particle p_{T} [GeV]; Tracking particles / 1.0 GeV", 100, 0, 100.0);
0423   TH1F* h_tp_pt_L = new TH1F("tp_pt_L", ";Tracking particle p_{T} [GeV]; Tracking particles / 0.1 GeV", 80, 0, 8.0);
0424   TH1F* h_tp_pt_LC = new TH1F("tp_pt_LC", ";Tracking particle p_{T} [GeV]; Tracking particles / 0.1 GeV", 80, 0, 8.0);
0425   TH1F* h_tp_pt_H = new TH1F("tp_pt_H", ";Tracking particle p_{T} [GeV]; Tracking particles / 1.0 GeV", 92, 8.0, 100.0);
0426   TH1F* h_tp_eta = new TH1F("tp_eta", ";Tracking particle #eta; Tracking particles / 0.1", 50, -2.5, 2.5);
0427   TH1F* h_tp_eta_L = new TH1F("tp_eta_L", ";Tracking particle #eta; Tracking particles / 0.1", 50, -2.5, 2.5);
0428   TH1F* h_tp_eta_H = new TH1F("tp_eta_H", ";Tracking particle #eta; Tracking particles / 0.1", 50, -2.5, 2.5);
0429   TH1F* h_tp_eta_23 = new TH1F("tp_eta_23", ";Tracking particle #eta; Tracking particles / 0.1", 50, -2.5, 2.5);
0430   TH1F* h_tp_eta_35 = new TH1F("tp_eta_35", ";Tracking particle #eta; Tracking particles / 0.1", 50, -2.5, 2.5);
0431   TH1F* h_tp_eta_5 = new TH1F("tp_eta_5", ";Tracking particle #eta; Tracking particles / 0.1", 50, -2.5, 2.5);
0432 
0433   TH1F* h_match_tp_pt =
0434       new TH1F("match_tp_pt", ";Tracking particle p_{T} [GeV]; Tracking particles / 1.0 GeV", 100, 0, 100.0);
0435   TH1F* h_match_tp_pt_L =
0436       new TH1F("match_tp_pt_L", ";Tracking particle p_{T} [GeV]; Tracking particles / 0.1 GeV", 80, 0, 8.0);
0437   TH1F* h_match_tp_pt_LC =
0438       new TH1F("match_tp_pt_LC", ";Tracking particle p_{T} [GeV]; Tracking particles / 0.1 GeV", 80, 0, 8.0);
0439   TH1F* h_match_tp_pt_H =
0440       new TH1F("match_tp_pt_H", ";Tracking particle p_{T} [GeV]; Tracking particles / 0.1 GeV", 92, 8.0, 100.0);
0441   TH1F* h_match_tp_eta = new TH1F("match_tp_eta", ";Tracking particle #eta; Tracking particles / 0.1", 50, -2.5, 2.5);
0442   TH1F* h_match_tp_eta_L =
0443       new TH1F("match_tp_eta_L", ";Tracking particle #eta; Tracking particles / 0.1", 50, -2.5, 2.5);
0444   TH1F* h_match_tp_eta_H =
0445       new TH1F("match_tp_eta_H", ";Tracking particle #eta; Tracking particles / 0.1", 50, -2.5, 2.5);
0446   TH1F* h_match_tp_eta_23 =
0447       new TH1F("match_tp_eta_23", ";Tracking particle #eta; Tracking particles / 0.1", 50, -2.5, 2.5);
0448   TH1F* h_match_tp_eta_35 =
0449       new TH1F("match_tp_eta_35", ";Tracking particle #eta; Tracking particles / 0.1", 50, -2.5, 2.5);
0450   TH1F* h_match_tp_eta_5 =
0451       new TH1F("match_tp_eta_5", ";Tracking particle #eta; Tracking particles / 0.1", 50, -2.5, 2.5);
0452 
0453   // ----------------------------------------------------------------------------------------------------------------
0454   // Tracklet propogation efficiencies vs. eta for seeding layers
0455 
0456   int trackletEffEtaBins = 24;
0457   double trackletEffMaxEta = 2.4;
0458   int numLayers = 11;
0459   TH2F* h_trk_tracklet_hits = new TH2F("trk_tracklet_hits",
0460                                        ";Track |#eta|; Layer index (0-5 = L1-6, 6-10 = D1-5)",
0461                                        trackletEffEtaBins,
0462                                        0,
0463                                        trackletEffMaxEta,
0464                                        11,
0465                                        0,
0466                                        11);  //used to create below hist
0467   TH2F* h_trk_tracklet_eff = new TH2F("trk_tracklet_eff",
0468                                       ";Track |#eta|; Layer index (0-5 = L1-6, 6-10 = D1-5)",
0469                                       trackletEffEtaBins,
0470                                       0,
0471                                       trackletEffMaxEta,
0472                                       11,
0473                                       0,
0474                                       11);
0475 
0476   // ----------------------------------------------------------------------------------------------------------------
0477   // resolution vs. pt histograms
0478 
0479   // ----------------------------------------------
0480   // for new versions of resolution vs pt/eta plots
0481   unsigned int nBinsPtRes = 500;
0482   double maxPtRes = 30.;
0483 
0484   unsigned int nBinsPtRelRes = 1000;
0485   double maxPtRelRes = 10.;
0486 
0487   unsigned int nBinsEtaRes = 500;
0488   double maxEtaRes = 0.1;
0489 
0490   unsigned int nBinsPhiRes = 500;
0491   double maxPhiRes = 0.2;
0492 
0493   unsigned int nBinsZ0Res = 100;
0494   double maxZ0Res = 4.0;
0495   // ----------------------------------------------
0496 
0497   const int nRANGE = 20;
0498   TString ptrange[nRANGE] = {"0-5",   "5-10",  "10-15", "15-20", "20-25", "25-30", "30-35", "35-40", "40-45", "45-50",
0499                              "50-55", "55-60", "60-65", "65-70", "70-75", "75-80", "80-85", "85-90", "90-95", "95-100"};
0500 
0501   const float pt_resmin = 1.5;
0502   const int nRANGE_L = 13;
0503   TString ptrange_L[nRANGE_L] = {"1.5-2",
0504                                  "2-2.5",
0505                                  "2.5-3",
0506                                  "3-3.5",
0507                                  "3.5-4",
0508                                  "4-4.5",
0509                                  "4.5-5",
0510                                  "5-5.5",
0511                                  "5.5-6",
0512                                  "6-6.5",
0513                                  "6.5-7",
0514                                  "7-7.5",
0515                                  "7.5-8"};
0516 
0517   TH1F* h_absResVsPt_pt[nRANGE];
0518   TH1F* h_absResVsPt_ptRel[nRANGE];
0519   TH1F* h_absResVsPt_z0[nRANGE];
0520   TH1F* h_absResVsPt_phi[nRANGE];
0521   TH1F* h_absResVsPt_eta[nRANGE];
0522   TH1F* h_absResVsPt_d0[nRANGE];
0523 
0524   TH1F* h_absResVsPt_pt_L[nRANGE_L];
0525   TH1F* h_absResVsPt_ptRel_L[nRANGE_L];
0526   TH1F* h_absResVsPt_z0_L[nRANGE_L];
0527   TH1F* h_absResVsPt_phi_L[nRANGE_L];
0528   TH1F* h_absResVsPt_eta_L[nRANGE_L];
0529   TH1F* h_absResVsPt_d0_L[nRANGE_L];
0530 
0531   for (int i = 0; i < nRANGE; i++) {
0532     h_absResVsPt_pt[i] = new TH1F(
0533         "absResVsPt_pt_" + ptrange[i], ";p_{T} residual (L1 - sim) [GeV]; L1 tracks / 0.1", nBinsPtRes, 0, maxPtRes);
0534     h_absResVsPt_ptRel[i] = new TH1F("absResVsPt_ptRel_" + ptrange[i],
0535                                      ";p_{T} residual (L1 - sim) / p_{T}; L1 tracks / 0.02",
0536                                      nBinsPtRelRes,
0537                                      0,
0538                                      maxPtRelRes);
0539     h_absResVsPt_z0[i] = new TH1F(
0540         "absResVsPt_z0_" + ptrange[i], ";z_{0} residual (L1 - sim) [GeV]; L1 tracks / 0.1", nBinsZ0Res, 0, maxZ0Res);
0541     h_absResVsPt_phi[i] = new TH1F(
0542         "absResVsPt_phi_" + ptrange[i], ";#phi residual (L1 - sim) [GeV]; L1 tracks / 0.1", nBinsPhiRes, 0, maxPhiRes);
0543     h_absResVsPt_eta[i] = new TH1F(
0544         "absResVsPt_eta_" + ptrange[i], ";#eta residual (L1 - sim) [GeV]; L1 tracks / 0.1", nBinsEtaRes, 0, maxEtaRes);
0545     h_absResVsPt_d0[i] =
0546         new TH1F("absResVsPt_d0_" + ptrange[i], ";d_{0}residual (L1 - sim) [GeV]; L1 tracks / 0.1", 100, 0, 0.02);
0547   }
0548 
0549   for (int i = 0; i < nRANGE_L; i++) {
0550     h_absResVsPt_pt_L[i] = new TH1F(
0551         "absResVsPt_L_pt_" + ptrange_L[i], ";p_{T} residual (L1 - sim) [GeV]; L1 tracks / 0.1", nBinsPtRes, 0, maxPtRes);
0552     h_absResVsPt_ptRel_L[i] = new TH1F("absResVsPt_L_ptRel_" + ptrange_L[i],
0553                                        ";p_{T} residual (L1 - sim) / p_{T}; L1 tracks / 0.02",
0554                                        nBinsPtRelRes,
0555                                        0,
0556                                        maxPtRelRes);
0557     h_absResVsPt_z0_L[i] = new TH1F(
0558         "absResVsPt_L_z0_" + ptrange_L[i], ";z_{0} residual (L1 - sim) [GeV]; L1 tracks / 0.1", nBinsZ0Res, 0, maxZ0Res);
0559     h_absResVsPt_phi_L[i] = new TH1F("absResVsPt_L_phi_" + ptrange_L[i],
0560                                      ";#phi residual (L1 - sim) [GeV]; L1 tracks / 0.1",
0561                                      nBinsPhiRes,
0562                                      0,
0563                                      maxPhiRes);
0564     h_absResVsPt_eta_L[i] = new TH1F("absResVsPt_L_eta_" + ptrange_L[i],
0565                                      ";#eta residual (L1 - sim) [GeV]; L1 tracks / 0.1",
0566                                      nBinsEtaRes,
0567                                      0,
0568                                      maxEtaRes);
0569     h_absResVsPt_d0_L[i] =
0570         new TH1F("absResVsPt_L_d0_" + ptrange_L[i], ";d_{0}residual (L1 - sim) [GeV]; L1 tracks / 0.1", 100, 0, 0.02);
0571   }
0572 
0573   // resolution vs. eta histograms
0574 
0575   const float eta_resmax = 2.5;
0576   const int nETARANGE = 25;
0577   TString etarange[nETARANGE] = {"0.1", "0.2", "0.3", "0.4", "0.5", "0.6", "0.7", "0.8", "0.9",
0578                                  "1.0", "1.1", "1.2", "1.3", "1.4", "1.5", "1.6", "1.7", "1.8",
0579                                  "1.9", "2.0", "2.1", "2.2", "2.3", "2.4", "2.5"};
0580   /*
0581   const int nETARANGE = 12;
0582   TString etarange[nETARANGE] = {"0.2","0.4","0.6","0.8","1.0",
0583                  "1.2","1.4","1.6","1.8","2.0",
0584                  "2.2","2.4"};
0585   */
0586 
0587   TH1F* h_absResVsEta_eta[nETARANGE];
0588   TH1F* h_absResVsEta_z0[nETARANGE];
0589   TH1F* h_absResVsEta_phi[nETARANGE];
0590   TH1F* h_absResVsEta_ptRel[nETARANGE];
0591   TH1F* h_absResVsEta_d0[nETARANGE];
0592 
0593   TH1F* h_absResVsEta_eta_L[nETARANGE];
0594   TH1F* h_absResVsEta_z0_L[nETARANGE];
0595   TH1F* h_absResVsEta_phi_L[nETARANGE];
0596   TH1F* h_absResVsEta_ptRel_L[nETARANGE];
0597   TH1F* h_absResVsEta_d0_L[nETARANGE];
0598 
0599   TH1F* h_absResVsEta_eta_H[nETARANGE];
0600   TH1F* h_absResVsEta_z0_H[nETARANGE];
0601   TH1F* h_absResVsEta_phi_H[nETARANGE];
0602   TH1F* h_absResVsEta_ptRel_H[nETARANGE];
0603   TH1F* h_absResVsEta_d0_H[nETARANGE];
0604 
0605   for (int i = 0; i < nETARANGE; i++) {
0606     h_absResVsEta_eta[i] = new TH1F(
0607         "absResVsEta_eta_" + etarange[i], ";#eta residual (L1 - sim) [GeV]; L1 tracks / 0.1", nBinsEtaRes, 0, maxEtaRes);
0608     h_absResVsEta_z0[i] = new TH1F(
0609         "absResVsEta_z0_" + etarange[i], ";|z_{0} residual (L1 - sim)| [cm]; L1 tracks / 0.01", nBinsZ0Res, 0, maxZ0Res);
0610     h_absResVsEta_phi[i] = new TH1F(
0611         "absResVsEta_phi_" + etarange[i], ";#phi residual (L1 - sim) [GeV]; L1 tracks / 0.1", nBinsPhiRes, 0, maxPhiRes);
0612     h_absResVsEta_ptRel[i] = new TH1F("absResVsEta_ptRel_" + etarange[i],
0613                                       ";p_{T} residual (L1 - sim) / p_{T}; L1 tracks / 0.02",
0614                                       nBinsPtRelRes,
0615                                       0,
0616                                       maxPtRelRes);
0617     h_absResVsEta_d0[i] =
0618         new TH1F("absResVsEta_d0_" + etarange[i], ";d_{0}residual (L1 - sim) [GeV]; L1 tracks / 0.1", 100, 0, 0.02);
0619 
0620     h_absResVsEta_eta_L[i] = new TH1F("absResVsEta_eta_L_" + etarange[i],
0621                                       ";#eta residual (L1 - sim) [GeV]; L1 tracks / 0.1",
0622                                       nBinsEtaRes,
0623                                       0,
0624                                       maxEtaRes);
0625     h_absResVsEta_z0_L[i] = new TH1F("absResVsEta_z0_L_" + etarange[i],
0626                                      ";|z_{0} residual (L1 - sim)| [cm]; L1 tracks / 0.01",
0627                                      nBinsZ0Res,
0628                                      0,
0629                                      maxZ0Res);
0630     h_absResVsEta_phi_L[i] = new TH1F("absResVsEta_phi_L_" + etarange[i],
0631                                       ";#phi residual (L1 - sim) [GeV]; L1 tracks / 0.1",
0632                                       nBinsPhiRes,
0633                                       0,
0634                                       maxPhiRes);
0635     h_absResVsEta_ptRel_L[i] = new TH1F("absResVsEta_ptRel_L_" + etarange[i],
0636                                         ";p_{T} residual (L1 - sim) / p_{T}; L1 tracks / 0.02",
0637                                         nBinsPtRelRes,
0638                                         0,
0639                                         maxPtRelRes);
0640     h_absResVsEta_d0_L[i] =
0641         new TH1F("absResVsEta_d0_L_" + etarange[i], ";d_{0}residual (L1 - sim) [GeV]; L1 tracks / 0.1", 100, 0, 0.02);
0642 
0643     h_absResVsEta_eta_H[i] = new TH1F("absResVsEta_eta_H_" + etarange[i],
0644                                       ";#eta residual (L1 - sim) [GeV]; L1 tracks / 0.1",
0645                                       nBinsEtaRes,
0646                                       0,
0647                                       maxEtaRes);
0648     h_absResVsEta_z0_H[i] = new TH1F("absResVsEta_z0_H_" + etarange[i],
0649                                      ";|z_{0} residual (L1 - sim)| [cm]; L1 tracks / 0.01",
0650                                      nBinsZ0Res,
0651                                      0,
0652                                      maxZ0Res);
0653     h_absResVsEta_phi_H[i] = new TH1F("absResVsEta_phi_H_" + etarange[i],
0654                                       ";#phi residual (L1 - sim) [GeV]; L1 tracks / 0.1",
0655                                       nBinsPhiRes,
0656                                       0,
0657                                       maxPhiRes);
0658     h_absResVsEta_ptRel_H[i] = new TH1F("absResVsEta_ptRel_H_" + etarange[i],
0659                                         ";p_{T} residual (L1 - sim) / p_{T}; L1 tracks / 0.02",
0660                                         nBinsPtRelRes,
0661                                         0,
0662                                         maxPtRelRes);
0663     h_absResVsEta_d0_H[i] =
0664         new TH1F("absResVsEta_d0_H_" + etarange[i], ";d_{0}residual (L1 - sim) [GeV]; L1 tracks / 0.1", 100, 0, 0.02);
0665   }
0666 
0667   // resolution vs phi
0668 
0669   const int nPHIRANGE = 32;
0670   TString phirange[nPHIRANGE] = {"-3.0", "-2.8", "-2.6", "-2.4", "-2.2", "-2.0", "-1.8", "-1.6", "-1.4", "-1.2", "-1.0",
0671                                  "-0.8", "-0.6", "-0.4", "-0.2", "0.0",  "0.2",  "0.4",  "0.6",  "0.8",  "1.0",  "1.2",
0672                                  "1.4",  "1.6",  "1.8",  "2.0",  "2.2",  "2.4",  "2.6",  "2.8",  "3.0",  "3.2"};
0673 
0674   TH1F* h_absResVsPhi_pt[nPHIRANGE];
0675   TH1F* h_absResVsPhi_ptRel[nPHIRANGE];
0676 
0677   for (int i = 0; i < nPHIRANGE; i++) {
0678     h_absResVsPhi_pt[i] = new TH1F(
0679         "absResVsPt_pt_" + phirange[i], ";p_{T} residual (L1 - sim) [GeV]; L1 tracks / 0.1", nBinsPtRes, 0, maxPtRes);
0680     h_absResVsPhi_ptRel[i] = new TH1F("absResVsPt_ptRel_" + phirange[i],
0681                                       ";p_{T} residual (L1 - sim) / p_{T}; L1 tracks / 0.02",
0682                                       nBinsPtRelRes,
0683                                       0,
0684                                       maxPtRelRes);
0685   }
0686 
0687   // ----------------------------------------------------------------------------------------------------------------
0688   // chi2 histograms (last bin is an overflow bin)
0689 
0690   TH1F* h_trk_chi2 = new TH1F("trk_chi2", ";#chi^{2}; L1 tracks / 1.0", 100, 0, 100);
0691   TH1F* h_trk_chi2_dof = new TH1F("trk_chi2_dof", ";#chi^{2} / D.O.F.; L1 tracks / 0.2", 100, 0, 20);
0692   TH1F* h_match_trk_chi2 = new TH1F("match_trk_chi2", ";#chi^{2}; L1 tracks / 1.0", 100, 0, 100);
0693   TH1F* h_match_trk_chi2_dof = new TH1F("match_trk_chi2_dof", ";#chi^{2} / D.O.F.; L1 tracks / 0.2", 100, 0, 20);
0694 
0695   TH1F* h_trk_chi2rphi = new TH1F("trk_chi2rphi", ";#chi^{2}_{r-#phi}; L1 tracks / 1.0", 100, 0, 100);
0696   TH1F* h_trk_chi2rphi_dof = new TH1F("trk_chi2rphi_dof", ";#chi^{2}_{r-#phi} / D.O.F.; L1 tracks / 0.2", 100, 0, 20);
0697   TH1F* h_match_trk_chi2rphi = new TH1F("match_trk_chi2rphi", ";#chi^{2}_{r-#phi}; L1 tracks / 1.0", 100, 0, 100);
0698   TH1F* h_match_trk_chi2rphi_dof =
0699       new TH1F("match_trk_chi2rphi_dof", ";#chi^{2}_{r-#phi} / D.O.F.; L1 tracks / 0.2", 100, 0, 20);
0700 
0701   TH1F* h_trk_chi2rz = new TH1F("trk_chi2rz", ";#chi^{2}_{r-z}; L1 tracks / 1.0", 100, 0, 100);
0702   TH1F* h_trk_chi2rz_dof = new TH1F("trk_chi2rz_dof", ";#chi^{2}_{r-z} / D.O.F.; L1 tracks / 0.2", 100, 0, 20);
0703   TH1F* h_match_trk_chi2rz = new TH1F("match_trk_chi2rz", ";#chi^{2}_{r-z}; L1 tracks / 1.0", 100, 0, 100);
0704   TH1F* h_match_trk_chi2rz_dof =
0705       new TH1F("match_trk_chi2rz_dof", ";#chi^{2}_{r-z} / D.O.F.; L1 tracks / 0.2", 100, 0, 20);
0706 
0707   // ----------------------------------------------------------------------------------------------------------------
0708   // total track rates
0709 
0710   TH1F* h_trk_all_vspt = new TH1F("trk_all_vspt", ";Track p_{T} [GeV]; ", 50, 0, 25);
0711   TH1F* h_trk_loose_vspt = new TH1F("trk_loose_vspt", ";Track p_{T} [GeV]; ", 50, 0, 25);
0712   TH1F* h_trk_genuine_vspt = new TH1F("trk_genuine_vspt", ";Track p_{T} [GeV]; ", 50, 0, 25);
0713   TH1F* h_trk_notloose_vspt = new TH1F(
0714       "trk_notloose_vspt", ";Track p_{T} [GeV]; ", 50, 0, 25);  //(same as "fake" according to the trk_fake labeling)
0715   TH1F* h_trk_notgenuine_vspt = new TH1F("trk_notgenuine_vspt", ";Track p_{T} [GeV]; ", 50, 0, 25);
0716   TH1F* h_trk_duplicate_vspt = new TH1F("trk_duplicate_vspt",
0717                                         ";Track p_{T} [GeV]; ",
0718                                         50,
0719                                         0,
0720                                         25);  //where a TP is genuinely matched to more than one L1 track
0721   TH1F* h_tp_vspt = new TH1F("tp_vspt", ";TP p_{T} [GeV]; ", 50, 0, 25);
0722 
0723   // ----------------------------------------------------------------------------------------------------------------
0724 
0725   TH1F* h_tp_z0 = new TH1F("tp_z0", ";Tracking particle z_{0} [cm]; Tracking particles / 1.0 cm", 50, -25.0, 25.0);
0726   TH1F* h_tp_z0_L = new TH1F("tp_z0_L", ";Tracking particle z_{0} [cm]; Tracking particles / 1.0 cm", 50, -25.0, 25.0);
0727   TH1F* h_tp_z0_H = new TH1F("tp_z0_H", ";Tracking particle z_{0} [cm]; Tracking particles / 1.0 cm", 50, -25.0, 25.0);
0728 
0729   TH1F* h_match_tp_z0 =
0730       new TH1F("match_tp_z0", ";Tracking particle z_{0} [cm]; Tracking particles / 1.0 cm", 50, -25.0, 25.0);
0731   TH1F* h_match_tp_z0_L =
0732       new TH1F("match_tp_z0_L", ";Tracking particle z_{0} [cm]; Tracking particles / 1.0 cm", 50, -25.0, 25.0);
0733   TH1F* h_match_tp_z0_H =
0734       new TH1F("match_tp_z0_H", ";Tracking particle z_{0} [cm]; Tracking particles / 1.0 cm", 50, -25.0, 25.0);
0735 
0736   // ----------------------------------------------------------------------------------------------------------------
0737   //
0738   //       ******************   additional histograms drawn if using 'detailed' option   ******************
0739   //
0740   // ----------------------------------------------------------------------------------------------------------------
0741 
0742   const float maxD0plot = TP_maxD0;
0743 
0744   TH1F* h_tp_phi = new TH1F("tp_phi", ";Tracking particle #phi [rad]; Tracking particles / 0.1", 64, -3.2, 3.2);
0745   TH1F* h_tp_d0 =
0746       new TH1F("tp_d0", ";Tracking particle d_{0} [cm]; Tracking particles / 0.01 cm", 50, -maxD0plot, maxD0plot);
0747   TH1F* h_tp_absd0 =
0748       new TH1F("tp_absd0", ";Tracking particle |d_{0}| [cm]; Tracking particles / 0.04 cm", 50, 0, maxD0plot);
0749   TH1F* h_tp_absd0_eta2 =
0750       new TH1F("tp_absd0_eta2", ";Tracking particle |d_{0}| [cm]; Tracking particles / 0.04 cm", 50, 0, maxD0plot);
0751   TH1F* h_tp_absd0_eta2_pt3 =
0752       new TH1F("tp_absd0_eta2_pt3", ";Tracking particle |d_{0}| [cm]; Tracking particles / 0.04 cm", 50, 0, maxD0plot);
0753 
0754   TH1F* h_match_tp_phi =
0755       new TH1F("match_tp_phi", ";Tracking particle #phi [rad]; Tracking particles / 0.1", 64, -3.2, 3.2);
0756   TH1F* h_match_tp_d0 =
0757       new TH1F("match_tp_d0", ";Tracking particle d_{0} [cm]; Tracking particles / 0.01 cm", 50, -maxD0plot, maxD0plot);
0758   TH1F* h_match_tp_absd0 =
0759       new TH1F("match_tp_absd0", ";Tracking particle d_{0} [cm]; Tracking particles / 0.04 cm", 50, 0, maxD0plot);
0760   TH1F* h_match_tp_absd0_eta2 =
0761       new TH1F("match_tp_absd0_eta2", ";Tracking particle d_{0} [cm]; Tracking particles / 0.04 cm", 50, 0, maxD0plot);
0762   TH1F* h_match_tp_absd0_eta2_pt3 = new TH1F(
0763       "match_tp_absd0_eta2_pt3", ";Tracking particle d_{0} [cm]; Tracking particles / 0.04 cm", 50, 0, maxD0plot);
0764 
0765   TH1F* h_match_trk_nstub = new TH1F("match_trk_nstub", ";Number of stubs; L1 tracks / 1.0", 15, 0, 15);
0766   TH1F* h_match_trk_nstub_C = new TH1F("match_trk_nstub_C", ";Number of stubs; L1 tracks / 1.0", 15, 0, 15);
0767   TH1F* h_match_trk_nstub_I = new TH1F("match_trk_nstub_I", ";Number of stubs; L1 tracks / 1.0", 15, 0, 15);
0768   TH1F* h_match_trk_nstub_F = new TH1F("match_trk_nstub_F", ";Number of stubs; L1 tracks / 1.0", 15, 0, 15);
0769 
0770   // note that we are only making the below chi2 histograms using the joint chi2, not the separate chi2rphi and chi2rz
0771 
0772   // chi2 histograms
0773   // note: last bin is an overflow bin
0774   TH1F* h_match_trk_chi2_C_L = new TH1F("match_trk_chi2_C_L", ";#chi^{2}; L1 tracks / 1.0", 100, 0, 100);
0775   TH1F* h_match_trk_chi2_I_L = new TH1F("match_trk_chi2_I_L", ";#chi^{2}; L1 tracks / 1.0", 100, 0, 100);
0776   TH1F* h_match_trk_chi2_F_L = new TH1F("match_trk_chi2_F_L", ";#chi^{2}; L1 tracks / 1.0", 100, 0, 100);
0777   TH1F* h_match_trk_chi2_C_H = new TH1F("match_trk_chi2_C_H", ";#chi^{2}; L1 tracks / 1.0", 100, 0, 100);
0778   TH1F* h_match_trk_chi2_I_H = new TH1F("match_trk_chi2_I_H", ";#chi^{2}; L1 tracks / 1.0", 100, 0, 100);
0779   TH1F* h_match_trk_chi2_F_H = new TH1F("match_trk_chi2_F_H", ";#chi^{2}; L1 tracks / 1.0", 100, 0, 100);
0780 
0781   // chi2/dof histograms
0782   // note: lastbin is an overflow bin
0783   TH1F* h_match_trk_chi2_dof_C_L =
0784       new TH1F("match_trk_chi2_dof_C_L", ";#chi^{2} / D.O.F.; L1 tracks / 0.2", 100, 0, 20);
0785   TH1F* h_match_trk_chi2_dof_I_L =
0786       new TH1F("match_trk_chi2_dof_I_L", ";#chi^{2} / D.O.F.; L1 tracks / 0.2", 100, 0, 20);
0787   TH1F* h_match_trk_chi2_dof_F_L =
0788       new TH1F("match_trk_chi2_dof_F_L", ";#chi^{2} / D.O.F.; L1 tracks / 0.2", 100, 0, 20);
0789   TH1F* h_match_trk_chi2_dof_C_H =
0790       new TH1F("match_trk_chi2_dof_C_H", ";#chi^{2} / D.O.F.; L1 tracks / 0.2", 100, 0, 20);
0791   TH1F* h_match_trk_chi2_dof_I_H =
0792       new TH1F("match_trk_chi2_dof_I_H", ";#chi^{2} / D.O.F.; L1 tracks / 0.2", 100, 0, 20);
0793   TH1F* h_match_trk_chi2_dof_F_H =
0794       new TH1F("match_trk_chi2_dof_F_H", ";#chi^{2} / D.O.F.; L1 tracks / 0.2", 100, 0, 20);
0795 
0796   // ----------------------------------------------------------------------------------------------------------------
0797   // resolution histograms
0798   TH1F* h_res_pt = new TH1F("res_pt", ";p_{T} residual (L1 - sim) [GeV]; L1 tracks / 0.05", 200, -5.0, 5.0);
0799   TH1F* h_res_ptRel = new TH1F("res_ptRel", ";p_{T} residual (L1 - sim) / p_{T}; L1 tracks / 0.01", 200, -1.0, 1.0);
0800   TH1F* h_res_eta = new TH1F("res_eta", ";#eta residual (L1 - sim); L1 tracks / 0.0002", 100, -0.01, 0.01);
0801   TH1F* h_res_phi = new TH1F("res_phi", ";#phi residual (L1 - sim) [rad]; L1 tracks / 0.0001", 100, -0.005, 0.005);
0802 
0803   TH1F* h_res_z0 = new TH1F("res_z0", ";z_{0} residual (L1 - sim) [cm]; L1 tracks / 0.02", 100, -1.0, 1.0);
0804   TH1F* h_res_z0_C = new TH1F("res_z0_C", ";z_{0} residual (L1 - sim) [cm]; L1 tracks / 0.02", 100, -1.0, 1.0);
0805   TH1F* h_res_z0_I = new TH1F("res_z0_I", ";z_{0} residual (L1 - sim) [cm]; L1 tracks / 0.02", 100, -1.0, 1.0);
0806   TH1F* h_res_z0_F = new TH1F("res_z0_F", ";z_{0} residual (L1 - sim) [cm]; L1 tracks / 0.02", 100, -1.0, 1.0);
0807   TH1F* h_res_z0_L = new TH1F("res_z0_L", ";z_{0} residual (L1 - sim) [cm]; L1 tracks / 0.02", 100, -1.0, 1.0);
0808   TH1F* h_res_z0_H = new TH1F("res_z0_H", ";z_{0} residual (L1 - sim) [cm]; L1 tracks / 0.02", 100, -1.0, 1.0);
0809 
0810   TH1F* h_res_z0_C_L =
0811       new TH1F("res_z0_C_L", ";z_{0} residual (L1 - sim) [cm]; L1 tracks / 0.02", 100, (-1) * 1.0, 1.0);
0812   TH1F* h_res_z0_I_L =
0813       new TH1F("res_z0_I_L", ";z_{0} residual (L1 - sim) [cm]; L1 tracks / 0.02", 100, (-1) * 1.0, 1.0);
0814   TH1F* h_res_z0_F_L =
0815       new TH1F("res_z0_F_L", ";z_{0} residual (L1 - sim) [cm]; L1 tracks / 0.02", 100, (-1) * 1.0, 1.0);
0816   TH1F* h_res_z0_C_H =
0817       new TH1F("res_z0_C_H", ";z_{0} residual (L1 - sim) [cm]; L1 tracks / 0.02", 100, (-1) * 1.0, 1.0);
0818   TH1F* h_res_z0_I_H =
0819       new TH1F("res_z0_I_H", ";z_{0} residual (L1 - sim) [cm]; L1 tracks / 0.02", 100, (-1) * 1.0, 1.0);
0820   TH1F* h_res_z0_F_H =
0821       new TH1F("res_z0_F_H", ";z_{0} residual (L1 - sim) [cm]; L1 tracks / 0.02", 100, (-1) * 1.0, 1.0);
0822 
0823   TH1F* h_res_d0 = new TH1F("res_d0", ";d_{0} residual (L1 - sim) [cm]; L1 tracks / 0.0002 cm", 200, -0.02, 0.02);
0824   TH1F* h_res_d0_C = new TH1F("res_d0_C", ";d_{0} residual (L1 - sim) [cm]; L1 tracks / 0.0001 cm", 200, -0.05, 0.05);
0825   TH1F* h_res_d0_I = new TH1F("res_d0_I", ";d_{0} residual (L1 - sim) [cm]; L1 tracks / 0.0001 cm", 200, -0.05, 0.05);
0826   TH1F* h_res_d0_F = new TH1F("res_d0_F", ";d_{0} residual (L1 - sim) [cm]; L1 tracks / 0.0001 cm", 200, -0.05, 0.05);
0827   TH1F* h_res_d0_L = new TH1F("res_d0_L", ";d_{0} residual (L1 - sim) [cm]; L1 tracks / 0.0001 cm", 200, -0.05, 0.05);
0828   TH1F* h_res_d0_H = new TH1F("res_d0_H", ";d_{0} residual (L1 - sim) [cm]; L1 tracks / 0.0001 cm", 200, -0.05, 0.05);
0829 
0830   TH1F* h_res_d0_C_L =
0831       new TH1F("res_d0_C_L", ";d_{0} residual (L1 - sim) [cm]; L1 tracks / 0.0001 cm", 200, -0.05, 0.05);
0832   TH1F* h_res_d0_I_L =
0833       new TH1F("res_d0_I_L", ";d_{0} residual (L1 - sim) [cm]; L1 tracks / 0.0001 cm", 200, -0.05, 0.05);
0834   TH1F* h_res_d0_F_L =
0835       new TH1F("res_d0_F_L", ";d_{0} residual (L1 - sim) [cm]; L1 tracks / 0.0001 cm", 200, -0.05, 0.05);
0836   TH1F* h_res_d0_C_H =
0837       new TH1F("res_d0_C_H", ";d_{0} residual (L1 - sim) [cm]; L1 tracks / 0.0001 cm", 200, -0.05, 0.05);
0838   TH1F* h_res_d0_I_H =
0839       new TH1F("res_d0_I_H", ";d_{0} residual (L1 - sim) [cm]; L1 tracks / 0.0001 cm", 200, -0.05, 0.05);
0840   TH1F* h_res_d0_F_H =
0841       new TH1F("res_d0_F_H", ";d_{0} residual (L1 - sim) [cm]; L1 tracks / 0.0001 cm", 200, -0.05, 0.05);
0842 
0843   // ----------------------------------------------------------------------------------------------------------------
0844   // more resolution vs pt
0845 
0846   TH1F* h_resVsPt_pt[nRANGE];
0847   TH1F* h_resVsPt_pt_C[nRANGE];
0848   TH1F* h_resVsPt_pt_I[nRANGE];
0849   TH1F* h_resVsPt_pt_F[nRANGE];
0850 
0851   TH1F* h_resVsPt_ptRel[nRANGE];
0852   TH1F* h_resVsPt_ptRel_C[nRANGE];
0853   TH1F* h_resVsPt_ptRel_I[nRANGE];
0854   TH1F* h_resVsPt_ptRel_F[nRANGE];
0855 
0856   TH1F* h_resVsPt_z0[nRANGE];
0857   TH1F* h_resVsPt_z0_C[nRANGE];
0858   TH1F* h_resVsPt_z0_I[nRANGE];
0859   TH1F* h_resVsPt_z0_F[nRANGE];
0860 
0861   TH1F* h_resVsPt_phi[nRANGE];
0862   TH1F* h_resVsPt_phi_C[nRANGE];
0863   TH1F* h_resVsPt_phi_I[nRANGE];
0864   TH1F* h_resVsPt_phi_F[nRANGE];
0865 
0866   TH1F* h_resVsPt_eta[nRANGE];
0867   TH1F* h_resVsPt_d0[nRANGE];
0868 
0869   for (int i = 0; i < nRANGE; i++) {
0870     h_resVsPt_pt[i] =
0871         new TH1F("resVsPt_pt_" + ptrange[i], ";p_{T} residual (L1 - sim) [GeV]; L1 tracks / 0.1", 100, -5.0, 5.0);
0872     h_resVsPt_pt_C[i] =
0873         new TH1F("resVsPt_pt_C_" + ptrange[i], ";p_{T} residual (L1 - sim) [GeV]; L1 tracks / 0.1", 100, -5.0, 5.0);
0874     h_resVsPt_pt_I[i] =
0875         new TH1F("resVsPt_pt_I_" + ptrange[i], ";p_{T} residual (L1 - sim) [GeV]; L1 tracks / 0.1", 100, -5.0, 5.0);
0876     h_resVsPt_pt_F[i] =
0877         new TH1F("resVsPt_pt_F_" + ptrange[i], ";p_{T} residual (L1 - sim) [GeV]; L1 tracks / 0.1", 100, -5.0, 5.0);
0878 
0879     h_resVsPt_ptRel[i] = new TH1F(
0880         "resVsPt_ptRel_" + ptrange[i], ";p_{T} residual (L1 - sim) / p_{T}; L1 tracks / 0.02", 300, -0.15, 0.15);
0881     h_resVsPt_ptRel_C[i] = new TH1F(
0882         "resVsPt_ptRel_c_" + ptrange[i], ";p_{T} residual (L1 - sim) / p_{T}; L1 tracks / 0.02", 300, -0.15, 0.15);
0883     h_resVsPt_ptRel_I[i] = new TH1F(
0884         "resVsPt_ptRel_I_" + ptrange[i], ";p_{T} residual (L1 - sim) / p_{T}; L1 tracks / 0.02", 300, -0.15, 0.15);
0885     h_resVsPt_ptRel_F[i] = new TH1F(
0886         "resVsPt_ptRel_F_" + ptrange[i], ";p_{T} residual (L1 - sim) / p_{T}; L1 tracks / 0.02", 300, -0.15, 0.15);
0887 
0888     h_resVsPt_z0[i] =
0889         new TH1F("resVsPt_z0_" + ptrange[i], ";z_{0} residual (L1 - sim) [cm]; L1 tracks / 0.02", 100, -1, 1);
0890     h_resVsPt_z0_C[i] =
0891         new TH1F("resVsPt_z0_C_" + ptrange[i], ";z_{0} residual (L1 - sim) [cm]; L1 tracks / 0.02", 100, -1, 1);
0892     h_resVsPt_z0_I[i] =
0893         new TH1F("resVsPt_z0_I_" + ptrange[i], ";z_{0} residual (L1 - sim) [cm]; L1 tracks / 0.02", 100, -1, 1);
0894     h_resVsPt_z0_F[i] =
0895         new TH1F("resVsPt_z0_F_" + ptrange[i], ";z_{0} residual (L1 - sim) [cm]; L1 tracks / 0.02", 100, -1, 1);
0896 
0897     h_resVsPt_phi[i] = new TH1F(
0898         "resVsPt_phi_" + ptrange[i], ";#phi residual (L1 - sim) [rad]; L1 tracks / 0.0001", 100, -0.005, 0.005);
0899     h_resVsPt_phi_C[i] = new TH1F(
0900         "resVsPt_phi_C_" + ptrange[i], ";#phi residual (L1 - sim) [rad]; L1 tracks / 0.0001", 100, -0.005, 0.005);
0901     h_resVsPt_phi_I[i] = new TH1F(
0902         "resVsPt_phi_I_" + ptrange[i], ";#phi residual (L1 - sim) [rad]; L1 tracks / 0.0001", 100, -0.005, 0.005);
0903     h_resVsPt_phi_F[i] = new TH1F(
0904         "resVsPt_phi_F_" + ptrange[i], ";#phi residual (L1 - sim) [rad]; L1 tracks / 0.0001", 100, -0.005, 0.005);
0905 
0906     h_resVsPt_eta[i] =
0907         new TH1F("resVsPt_eta_" + ptrange[i], ";#eta residual (L1 - sim); L1 tracks / 0.0002", 100, -0.01, 0.01);
0908 
0909     h_resVsPt_d0[i] =
0910         new TH1F("resVsPt_d0_" + ptrange[i], ";d_{0} residual (L1 - sim) [cm]; L1 tracks / 0.0004", 100, -0.02, 0.02);
0911   }
0912 
0913   // ----------------------------------------------------------------------------------------------------------------
0914   // more resolution vs eta
0915 
0916   TH1F* h_resVsEta_eta[nETARANGE];
0917   TH1F* h_resVsEta_eta_L[nETARANGE];
0918   TH1F* h_resVsEta_eta_H[nETARANGE];
0919 
0920   TH1F* h_resVsEta_z0[nETARANGE];
0921   TH1F* h_resVsEta_z0_L[nETARANGE];
0922   TH1F* h_resVsEta_z0_H[nETARANGE];
0923 
0924   TH1F* h_resVsEta_phi[nETARANGE];
0925   TH1F* h_resVsEta_phi_L[nETARANGE];
0926   TH1F* h_resVsEta_phi_H[nETARANGE];
0927 
0928   TH1F* h_resVsEta_ptRel[nETARANGE];
0929   TH1F* h_resVsEta_ptRel_L[nETARANGE];
0930   TH1F* h_resVsEta_ptRel_H[nETARANGE];
0931 
0932   TH1F* h_resVsEta_d0[nETARANGE];
0933   TH1F* h_resVsEta_d0_L[nETARANGE];
0934   TH1F* h_resVsEta_d0_H[nETARANGE];
0935 
0936   for (int i = 0; i < nETARANGE; i++) {
0937     h_resVsEta_eta[i] =
0938         new TH1F("resVsEta2_eta_" + etarange[i], ";#eta residual (L1 - sim); L1 tracks / 0.0002", 100, -0.01, 0.01);
0939     h_resVsEta_eta_L[i] =
0940         new TH1F("resVsEta2_eta_L_" + etarange[i], ";#eta residual (L1 - sim); L1 tracks / 0.0002", 100, -0.01, 0.01);
0941     h_resVsEta_eta_H[i] =
0942         new TH1F("resVsEta2_eta_H_" + etarange[i], ";#eta residual (L1 - sim); L1 tracks / 0.0002", 100, -0.01, 0.01);
0943 
0944     h_resVsEta_z0[i] =
0945         new TH1F("resVsEta2_z0_" + etarange[i], ";z_{0} residual (L1 - sim) [cm]; L1 tracks / 0.01", 100, -1, 1);
0946     h_resVsEta_z0_L[i] =
0947         new TH1F("resVsEta2_z0_L_" + etarange[i], ";z_{0} residual (L1 - sim) [cm]; L1 tracks / 0.01", 100, -1, 1);
0948     h_resVsEta_z0_H[i] =
0949         new TH1F("resVsEta2_z0_H_" + etarange[i], ";z_{0} residual (L1 - sim) [cm]; L1 tracks / 0.01", 100, -1, 1);
0950 
0951     h_resVsEta_phi[i] = new TH1F(
0952         "resVsEta2_phi_" + etarange[i], ";#phi residual (L1 - sim) [rad]; L1 tracks / 0.0001", 100, -0.005, 0.005);
0953     h_resVsEta_phi_L[i] = new TH1F(
0954         "resVsEta2_phi_L_" + etarange[i], ";#phi residual (L1 - sim) [rad]; L1 tracks / 0.0001", 100, -0.005, 0.005);
0955     h_resVsEta_phi_H[i] = new TH1F(
0956         "resVsEta2_phi_H_" + etarange[i], ";#phi residual (L1 - sim) [rad]; L1 tracks / 0.0001", 100, -0.005, 0.005);
0957 
0958     h_resVsEta_ptRel[i] = new TH1F(
0959         "resVsEta2_ptRel_" + etarange[i], ";p_{T} residual (L1 - sim) / p_{T}; L1 tracks / 0.01", 100, -0.5, 0.5);
0960     h_resVsEta_ptRel_L[i] = new TH1F(
0961         "resVsEta2_ptRel_L_" + etarange[i], ";p_{T} residual (L1 - sim) / p_{T}; L1 tracks / 0.02", 100, -0.1, 0.1);
0962     h_resVsEta_ptRel_H[i] = new TH1F(
0963         "resVsEta2_ptRel_H_" + etarange[i], ";p_{T} residual (L1 - sim) / p_{T}; L1 tracks / 0.02", 100, -0.25, 0.25);
0964 
0965     h_resVsEta_d0[i] =
0966         new TH1F("resVsEta2_d0_" + etarange[i], ";d_{0} residual (L1 - sim) [cm]; L1 tracks / 0.004", 100, -0.02, 0.02);
0967     h_resVsEta_d0_L[i] = new TH1F(
0968         "resVsEta2_d0_L_" + etarange[i], ";d_{0} residual (L1 - sim) [cm]; L1 tracks / 0.004", 100, -0.02, 0.02);
0969     h_resVsEta_d0_H[i] = new TH1F(
0970         "resVsEta2_d0_H_" + etarange[i], ";d_{0} residual (L1 - sim) [cm]; L1 tracks / 0.004", 100, -0.02, 0.02);
0971   }
0972   // ----------------------------------------------------------------------------------------------------------------
0973 
0974   // ----------------------------------------------------------------------------------------------------------------
0975   // additional ones for sum pt in jets
0976   /*
0977   TH1F* h_jet_tp_sumpt_vspt  = new TH1F("jet_tp_sumpt_vspt",  ";sum(TP p_{T}) [GeV]; Gen jets / 5.0 GeV", 20, 0, 200.0);
0978   TH1F* h_jet_trk_sumpt_vspt = new TH1F("jet_trk_sumpt_vspt", ";sum(TP p_{T}) [GeV]; Gen jets / 5.0 GeV", 20, 0, 200.0);
0979   TH1F* h_jet_matchtrk_sumpt_vspt = new TH1F("jet_matchtrk_sumpt_vspt", ";sum(TP p_{T}) [GeV]; Gen jets / 5.0 GeV", 20, 0, 200.0);
0980 
0981   TH1F* h_jet_tp_sumpt_vseta  = new TH1F("jet_tp_sumpt_vseta",  ";Gen jet #eta; Gen jets / 0.2", 24, -2.4, 2.4);
0982   TH1F* h_jet_trk_sumpt_vseta = new TH1F("jet_trk_sumpt_vseta", ";Gen jet #eta; Gen jets / 0.2", 24, -2.4, 2.4);
0983   TH1F* h_jet_matchtrk_sumpt_vseta = new TH1F("jet_matchtrk_sumpt_vseta", ";Gen jet #eta; Gen jets / 0.2", 24, -2.4, 2.4);
0984 
0985   h_jet_tp_sumpt_vseta->Sumw2();
0986   h_jet_tp_sumpt_vspt->Sumw2();
0987   h_jet_trk_sumpt_vseta->Sumw2();
0988   h_jet_trk_sumpt_vspt->Sumw2();
0989   h_jet_matchtrk_sumpt_vseta->Sumw2();
0990   h_jet_matchtrk_sumpt_vspt->Sumw2();
0991   */
0992 
0993   // ----------------------------------------------------------------------------------------------------------------
0994   // number of tracks per event
0995 
0996   // all tracks
0997   TH1F* h_ntrk_pt2 = new TH1F("ntrk_pt2", ";# tracks (p_{T} > 2 GeV) / event; Events", 400, 0, 400.0);
0998   TH1F* h_ntrk_pt3 = new TH1F("ntrk_pt3", ";# tracks (p_{T} > 3 GeV) / event; Events", 300, 0, 300.0);
0999   TH1F* h_ntrk_pt10 = new TH1F("ntrk_pt10", ";# tracks (p_{T} > 10 GeV) / event; Events", 100, 0, 100.0);
1000 
1001   // tracks flagged as genuine (this would include duplicates (?))
1002   TH1F* h_ntrk_genuine_pt2 =
1003       new TH1F("ntrk_genuine_pt2", ";# genuine tracks (p_{T} > 2 GeV) / event; Events", 400, 0, 400.0);
1004   TH1F* h_ntrk_genuine_pt3 =
1005       new TH1F("ntrk_genuine_pt3", ";# genuine tracks (p_{T} > 3 GeV) / event; Events", 300, 0, 300.0);
1006   TH1F* h_ntrk_genuine_pt10 =
1007       new TH1F("ntrk_genuine_pt10", ";# genuine tracks (p_{T} > 10 GeV) / event; Events", 100, 0, 100.0);
1008 
1009   // Max N tracks from a sector per event
1010   TH1F* h_ntrkPerSector_all =
1011       new TH1F("ntrkPerSector_all", ";Max. # tracks from a sector / event; Events", 50, 0, 100.0);
1012   TH1F* h_ntrkPerSector_pt2 =
1013       new TH1F("ntrkPerSector_pt2", ";Max. # tracks from a sector (p_{T} > 2 GeV) / event; Events", 50, 0, 100.0);
1014   TH1F* h_ntrkPerSector_pt3 =
1015       new TH1F("ntrkPerSector_pt3", ";Max. # tracks from a sector (p_{T} > 3 GeV) / event; Events", 50, 0, 100.0);
1016   TH1F* h_ntrkPerSector_pt4 =
1017       new TH1F("ntrkPerSector_pt4", ";Max. # tracks from a sector (p_{T} > 10 GeV) / event; Events", 50, 0, 100.0);
1018 
1019   // number of tracks vs. efficiency (eta, pT)
1020   TH1F* h_trk_pt = new TH1F("trk_pt", Form(";Track p_{T} (GeV);Tracks / 0.5 GeV"), 200, 0., 100.);
1021   TH1F* h_trk_eta = new TH1F("trk_eta", Form(";Track #eta;Tracks / 0.026"), 200, -2.6, 2.6);
1022 
1023   // ----------------------------------------------------------------------------------------------------------------
1024   //        * * * * *     S T A R T   O F   A C T U A L   R U N N I N G   O N   E V E N T S     * * * * *
1025   // ----------------------------------------------------------------------------------------------------------------
1026 
1027   int nevt = tree->GetEntries();
1028 
1029   // ----------------------------------------------------------------------------------------------------------------
1030   // event loop
1031   for (int i = 0; i < nevt; i++) {
1032     tree->GetEntry(i, 0);
1033 
1034     /*
1035     // ----------------------------------------------------------------------------------------------------------------
1036     // sumpt in jets
1037     if (TP_select_injet > 0) {
1038       for (int ij=0; ij<(int)jet_tp_sumpt->size(); ij++) {
1039 
1040     float fraction = 0;
1041     float fractionMatch = 0;
1042     if (jet_tp_sumpt->at(ij) > 0) {
1043       fraction = jet_trk_sumpt->at(ij)/jet_tp_sumpt->at(ij);
1044       fractionMatch = jet_matchtrk_sumpt->at(ij)/jet_tp_sumpt->at(ij);
1045     }
1046 
1047     h_jet_tp_sumpt_vspt->Fill(jet_tp_sumpt->at(ij),1.0);
1048     h_jet_trk_sumpt_vspt->Fill(jet_tp_sumpt->at(ij),fraction);
1049     h_jet_matchtrk_sumpt_vspt->Fill(jet_tp_sumpt->at(ij),fractionMatch);
1050 
1051     h_jet_tp_sumpt_vseta->Fill(jet_eta->at(ij),1.0);
1052     h_jet_trk_sumpt_vseta->Fill(jet_eta->at(ij),fraction);
1053     h_jet_matchtrk_sumpt_vseta->Fill(jet_eta->at(ij),fractionMatch);
1054       }
1055     }
1056     */
1057 
1058     // ----------------------------------------------------------------------------------------------------------------
1059     // track loop for total rates & fake rates.
1060 
1061     int ntrkevt_pt2 = 0;
1062     int ntrkevt_pt3 = 0;
1063     int ntrkevt_pt10 = 0;
1064 
1065     int ntrkevt_genuine_pt2 = 0;
1066     int ntrkevt_genuine_pt3 = 0;
1067     int ntrkevt_genuine_pt10 = 0;
1068 
1069     vector<unsigned int> nTrksPerSector_all(9, 0);
1070     vector<unsigned int> nTrksPerSector_pt2(9, 0);
1071     vector<unsigned int> nTrksPerSector_pt3(9, 0);
1072     vector<unsigned int> nTrksPerSector_pt4(9, 0);
1073 
1074     for (int it = 0; it < (int)trk_pt->size(); it++) {  // Loop reco tracks
1075       // ----------------------------------------------------------------------------------------------------------------
1076       // track properties
1077 
1078       // ----------------------------------------------------------------------------------------------------------------
1079       // Fill number of tracks vs track param
1080       h_trk_pt->Fill(trk_pt->at(it));
1081       h_trk_eta->Fill(trk_eta->at(it));
1082 
1083       // fill all trk chi2 & chi2/dof histograms, including for chi2 r-phi and chi2 r-z
1084       float chi2 = trk_chi2->at(it);
1085       float chi2dof = trk_chi2_dof->at(it);
1086       float chi2rphi = trk_chi2rphi->at(it);
1087       float chi2rphidof = trk_chi2rphi_dof->at(it);
1088       float chi2rz = trk_chi2rz->at(it);
1089       float chi2rzdof = trk_chi2rz_dof->at(it);
1090 
1091       // create overflow bins by restricting range of chi2
1092       int chi2Overflow = 100;
1093       int chi2DOFOverflow = 20;  //apprx chi2Overflow / avg. nstubs
1094       double buffer = 0.1;
1095 
1096       if (chi2 > chi2Overflow)
1097         chi2 = chi2Overflow - buffer;
1098       if (chi2dof > chi2DOFOverflow)
1099         chi2dof = chi2DOFOverflow - buffer;
1100       if (chi2rphi > chi2Overflow)
1101         chi2rphi = chi2Overflow - buffer;
1102       if (chi2rphidof > chi2DOFOverflow)
1103         chi2rphidof = chi2DOFOverflow - buffer;
1104       if (chi2rz > chi2Overflow)
1105         chi2rz = chi2Overflow - buffer;
1106       if (chi2rzdof > chi2DOFOverflow)
1107         chi2rzdof = chi2DOFOverflow - buffer;
1108 
1109       if (trk_pt->at(it) > TP_minPt) {  //TRK pt > TP_minPt
1110 
1111         h_trk_chi2->Fill(chi2);
1112         h_trk_chi2_dof->Fill(chi2dof);
1113 
1114         h_trk_chi2rphi->Fill(chi2rphi);
1115         h_trk_chi2rphi_dof->Fill(chi2rphidof);
1116 
1117         h_trk_chi2rz->Fill(chi2rz);
1118         h_trk_chi2rz_dof->Fill(chi2rzdof);
1119 
1120       }  //end TRK pt > TP_minPt
1121 
1122       // ----------------------------------------------------------------------------------------------------------------
1123       // look at track rate and fake rate, etc.
1124 
1125       // only look at tracks in (ttbar) jets ?
1126       if (TP_select_injet > 0) {
1127         if (TP_select_injet == 1 && trk_injet->at(it) == 0)
1128           continue;
1129         if (TP_select_injet == 2 && trk_injet_highpt->at(it) == 0)
1130           continue;
1131         if (TP_select_injet == 3 && trk_injet_vhighpt->at(it) == 0)
1132           continue;
1133       }
1134       ntrk++;
1135       if (trk_pt->at(it) >= 0.0)
1136         ++nTrksPerSector_all.at(trk_phiSector->at(it) % 9);
1137       if (std::abs(trk_eta->at(it)) > TP_maxEta)
1138         continue;
1139       if (trk_pt->at(it) < TP_minPt)
1140         continue;
1141 
1142       // Uncomment these cuts to see effect on rate & fake rate.
1143       //if (trk_chi2->at(it) > L1Tk_maxChi2) continue;
1144       //if (trk_chi2_dof->at(it) > L1Tk_maxChi2dof) continue;
1145       //if (trk_nstub->at(it) < L1Tk_minNstub) continue;
1146 
1147       // Tracklet & Hybrid have 9 sectors, but TMTT has 18 (with sectors 0 & 1 in nonant 0 etc).
1148       // As don't know here with algo used, "% 9" added to prevent crash, but not correct for TMTT.
1149       if (trk_pt->at(it) > 2.0)
1150         ++nTrksPerSector_pt2.at(trk_phiSector->at(it) % 9);
1151       if (trk_pt->at(it) > 3.0)
1152         ++nTrksPerSector_pt3.at(trk_phiSector->at(it) % 9);
1153       if (trk_pt->at(it) > 4.0)
1154         ++nTrksPerSector_pt4.at(trk_phiSector->at(it) % 9);
1155 
1156       if (trk_pt->at(it) > 2.0) {
1157         ntrk_pt2++;
1158         ntrkevt_pt2++;
1159         h_trk_all_vspt->Fill(trk_pt->at(it));
1160         if (trk_genuine->at(it) == 1) {
1161           ntrk_genuine_pt2++;
1162           ntrkevt_genuine_pt2++;
1163           h_trk_genuine_vspt->Fill(trk_pt->at(it));
1164         } else
1165           h_trk_notgenuine_vspt->Fill(trk_pt->at(it));
1166         if (trk_loose->at(it) == 1)
1167           h_trk_loose_vspt->Fill(trk_pt->at(it));
1168         else
1169           h_trk_notloose_vspt->Fill(trk_pt->at(it));
1170       }
1171       if (trk_pt->at(it) > 3.0) {
1172         ntrk_pt3++;
1173         ntrkevt_pt3++;
1174         if (trk_genuine->at(it) == 1)
1175           ntrkevt_genuine_pt3++;
1176       }
1177       if (trk_pt->at(it) > 10.0) {
1178         ntrk_pt10++;
1179         ntrkevt_pt10++;
1180         if (trk_genuine->at(it) == 1)
1181           ntrkevt_genuine_pt10++;
1182       }
1183 
1184       // ----------------------------------------------------------------------------------------------------------------
1185       // Fill tracklet propogation efficiency histo
1186 
1187       // create an 11-bit long iterable from lhits and dhits
1188       int num_layers = 6;
1189       int num_discs = 5;
1190       int lhits = trk_lhits->at(it);
1191       int dhits = trk_dhits->at(it);
1192       std::vector<int> layers = {};
1193       for (int layer_index = 0; layer_index < num_layers + num_discs; layer_index++) {
1194         if (layer_index < num_layers) {
1195           layers.push_back(lhits % 10);
1196           lhits /= 10;
1197         } else {
1198           layers.push_back(dhits % 10);
1199           dhits /= 10;
1200         }
1201       }
1202 
1203       for (unsigned int layer = 0; layer < layers.size(); layer++) {
1204         if (layers.at(layer)) {                                         // if there was a hit at this layer...
1205           h_trk_tracklet_hits->Fill(std::abs(trk_eta->at(it)), layer);  // ...fill this bin with the layer of the track.
1206         }
1207       }
1208     }  // End loop reco tracks
1209 
1210     h_ntrk_pt2->Fill(ntrkevt_pt2);
1211     h_ntrk_pt3->Fill(ntrkevt_pt3);
1212     h_ntrk_pt10->Fill(ntrkevt_pt10);
1213 
1214     h_ntrk_genuine_pt2->Fill(ntrkevt_genuine_pt2);
1215     h_ntrk_genuine_pt3->Fill(ntrkevt_genuine_pt3);
1216     h_ntrk_genuine_pt10->Fill(ntrkevt_genuine_pt10);
1217 
1218     h_ntrkPerSector_all->Fill(*std::max_element(nTrksPerSector_all.begin(), nTrksPerSector_all.end()));
1219     h_ntrkPerSector_pt2->Fill(*std::max_element(nTrksPerSector_pt2.begin(), nTrksPerSector_pt2.end()));
1220     h_ntrkPerSector_pt3->Fill(*std::max_element(nTrksPerSector_pt3.begin(), nTrksPerSector_pt3.end()));
1221     h_ntrkPerSector_pt4->Fill(*std::max_element(nTrksPerSector_pt4.begin(), nTrksPerSector_pt4.end()));
1222 
1223     // ----------------------------------------------------------------------------------------------------------------
1224     // Loop tracking particles
1225     for (int it = 0; it < (int)tp_pt->size(); it++) {
1226       // only look at TPs in (ttbar) jets ?
1227       if (TP_select_injet > 0) {
1228         if (TP_select_injet == 1 && tp_injet->at(it) == 0)
1229           continue;
1230         if (TP_select_injet == 2 && tp_injet_highpt->at(it) == 0)
1231           continue;
1232         if (TP_select_injet == 3 && tp_injet_vhighpt->at(it) == 0)
1233           continue;
1234       }
1235 
1236       // cut on PDG ID at plot stage?
1237       if (TP_select_pdgid != 0) {
1238         if (abs(tp_pdgid->at(it)) != abs(TP_select_pdgid))
1239           continue;
1240       }
1241 
1242       // kinematic cuts
1243       if (std::abs(tp_dxy->at(it)) > TP_maxDxy)
1244         continue;
1245       if (std::abs(tp_d0->at(it)) > TP_maxD0)
1246         continue;
1247       if (tp_pt->at(it) < 0.2)
1248         continue;
1249       if (tp_pt->at(it) > TP_maxPt)
1250         continue;
1251       if (std::abs(tp_eta->at(it)) > TP_maxEta)
1252         continue;
1253 
1254       // total track rates
1255       if (tp_pt->at(it) > TP_minPt) {
1256         if (tp_pt->at(it) > 2.0) {
1257           ntp_pt2++;
1258           h_tp_vspt->Fill(tp_pt->at(it));
1259           // duplicate rate
1260           if (tp_nmatch->at(it) > 1) {
1261             for (int inm = 1; inm < tp_nmatch->at(it); inm++) {  // N.B. Loop doesn't start at zero.
1262               ntp_nmatch++;
1263               h_trk_duplicate_vspt->Fill(matchtrk_pt->at(it));
1264             }
1265           }
1266         }
1267         if (tp_pt->at(it) > 3.0)
1268           ntp_pt3++;
1269         if (tp_pt->at(it) > 10.0)
1270           ntp_pt10++;
1271       }
1272 
1273       // cut on event ID (eventid=0 means the TP is from the primary interaction, so *not* selecting only eventid=0 means including stuff from pileup)
1274       if (TP_select_eventid == 0 && tp_eventid->at(it) != 0)
1275         continue;
1276 
1277       // look at failure scenarios?
1278       if (useDeadRegion) {
1279         if (tp_phi->at(it) < 0 || tp_phi->at(it) > 1)
1280           continue;
1281       }
1282 
1283       h_tp_pt->Fill(tp_pt->at(it));
1284       if (tp_pt->at(it) < 8.0)
1285         h_tp_pt_L->Fill(tp_pt->at(it));
1286       else
1287         h_tp_pt_H->Fill(tp_pt->at(it));
1288       if (tp_pt->at(it) < 8.0 && std::abs(tp_eta->at(it)) < 1.0)
1289         h_tp_pt_LC->Fill(tp_pt->at(it));
1290 
1291       if (tp_pt->at(it) > TP_minPt) {
1292         if (std::abs(tp_eta->at(it)) < 1.0)
1293           n_all_eta1p0++;
1294         else if (std::abs(tp_eta->at(it)) < 1.75)
1295           n_all_eta1p75++;
1296         else
1297           n_all_eta2p5++;
1298 
1299         if (std::abs(tp_pt->at(it)) > 2.0)
1300           n_all_ptg2++;
1301         if (std::abs(tp_pt->at(it)) > 2.0 && std::abs(tp_pt->at(it)) < 8.0)
1302           n_all_pt2to8++;
1303         if (std::abs(tp_pt->at(it)) > 8.0)
1304           n_all_ptg8++;
1305         if (std::abs(tp_pt->at(it)) > 40.0)
1306           n_all_ptg40++;
1307 
1308         h_tp_eta->Fill(tp_eta->at(it));
1309         h_tp_phi->Fill(tp_phi->at(it));
1310         h_tp_z0->Fill(tp_z0->at(it));
1311         h_tp_d0->Fill(tp_d0->at(it));
1312         h_tp_absd0->Fill(std::abs(tp_d0->at(it)));
1313         if (std::abs(tp_eta->at(it)) < 2.0)
1314           h_tp_absd0_eta2->Fill(std::abs(tp_d0->at(it)));
1315         if (std::abs(tp_eta->at(it)) < 2.0 && tp_pt->at(it) > 3.0)
1316           h_tp_absd0_eta2_pt3->Fill(std::abs(tp_d0->at(it)));
1317 
1318         if (tp_pt->at(it) < 3.0)
1319           h_tp_eta_23->Fill(tp_eta->at(it));
1320         else if (tp_pt->at(it) < 5.0)
1321           h_tp_eta_35->Fill(tp_eta->at(it));
1322         else
1323           h_tp_eta_5->Fill(tp_eta->at(it));
1324 
1325         if (tp_pt->at(it) < 8.0) {
1326           h_tp_eta_L->Fill(tp_eta->at(it));
1327           h_tp_z0_L->Fill(tp_z0->at(it));
1328         } else {
1329           h_tp_eta_H->Fill(tp_eta->at(it));
1330           h_tp_z0_H->Fill(tp_z0->at(it));
1331         }
1332       }
1333 
1334       // ----------------------------------------------------------------------------------------------------------------
1335       // was the tracking particle matched to a L1 track?
1336       if (tp_nmatch->at(it) < 1)
1337         continue;
1338 
1339       // ----------------------------------------------------------------------------------------------------------------
1340       // use only tracks with min X stubs
1341       if (matchtrk_nstub->at(it) < L1Tk_minNstub)
1342         continue;
1343 
1344       int thisseed = matchtrk_seed->at(it);
1345       if (thisseed > 25)
1346         thisseed = thisseed - 20;
1347       if ((L1Tk_seed != 0) && (thisseed != L1Tk_seed))
1348         continue;
1349 
1350       // ----------------------------------------------------------------------------------------------------------------
1351       // fill matchtrk chi2 & chi2/dof histograms before making chi2 cut
1352 
1353       float chi2 = matchtrk_chi2->at(it);
1354       float chi2dof = matchtrk_chi2_dof->at(it);
1355       float chi2rphi = matchtrk_chi2rphi->at(it);
1356       float chi2rphidof = matchtrk_chi2rphi_dof->at(it);
1357       float chi2rz = matchtrk_chi2rz->at(it);
1358       float chi2rzdof = matchtrk_chi2rz_dof->at(it);
1359 
1360       // create overflow bins by restricting range of chi2
1361       int chi2Overflow = 100;
1362       int chi2DOFOverflow = 20;  //apprx chi2Overflow / avg. nstubs
1363       double buffer = 0.1;
1364 
1365       if (chi2 > chi2Overflow)
1366         chi2 = chi2Overflow - buffer;
1367       if (chi2dof > chi2DOFOverflow)
1368         chi2dof = chi2DOFOverflow - buffer;
1369       if (chi2rphi > chi2Overflow)
1370         chi2rphi = chi2Overflow - buffer;
1371       if (chi2rphidof > chi2DOFOverflow)
1372         chi2rphidof = chi2DOFOverflow - buffer;
1373       if (chi2rz > chi2Overflow)
1374         chi2rz = chi2Overflow - buffer;
1375       if (chi2rzdof > chi2DOFOverflow)
1376         chi2rzdof = chi2DOFOverflow - buffer;
1377 
1378       if (tp_pt->at(it) > TP_minPt) {  //TP pt > TP_minPt
1379 
1380         h_match_trk_chi2->Fill(chi2);
1381         h_match_trk_chi2_dof->Fill(chi2dof);
1382 
1383         h_match_trk_chi2rphi->Fill(chi2rphi);
1384         h_match_trk_chi2rphi_dof->Fill(chi2rphidof);
1385 
1386         h_match_trk_chi2rz->Fill(chi2rz);
1387         h_match_trk_chi2rz_dof->Fill(chi2rzdof);
1388 
1389         // central eta
1390         if (std::abs(tp_eta->at(it)) < 0.8) {
1391           if (tp_pt->at(it) < 8.0) {
1392             h_match_trk_chi2_C_L->Fill(chi2);
1393             h_match_trk_chi2_dof_C_L->Fill(chi2dof);
1394           } else {
1395             h_match_trk_chi2_C_H->Fill(chi2);
1396             h_match_trk_chi2_dof_C_H->Fill(chi2dof);
1397           }
1398         }
1399         // intermediate eta
1400         else if (std::abs(tp_eta->at(it)) < 1.6 && std::abs(tp_eta->at(it)) >= 0.8) {
1401           if (tp_pt->at(it) < 8.0) {
1402             h_match_trk_chi2_I_L->Fill(chi2);
1403             h_match_trk_chi2_dof_I_L->Fill(chi2dof);
1404           } else {
1405             h_match_trk_chi2_I_H->Fill(chi2);
1406             h_match_trk_chi2_dof_I_H->Fill(chi2dof);
1407           }
1408         }
1409         // forward eta
1410         else if (std::abs(tp_eta->at(it)) >= 1.6) {
1411           if (tp_pt->at(it) < 8.0) {
1412             h_match_trk_chi2_F_L->Fill(chi2);
1413             h_match_trk_chi2_dof_F_L->Fill(chi2dof);
1414           } else {
1415             h_match_trk_chi2_F_H->Fill(chi2);
1416             h_match_trk_chi2_dof_F_H->Fill(chi2dof);
1417           }
1418         }
1419       }  //end TP pt > TP_minPt
1420 
1421       // ----------------------------------------------------------------------------------------------------------------
1422       // cut on chi2?
1423       if (matchtrk_chi2->at(it) > L1Tk_maxChi2)
1424         continue;
1425       if (matchtrk_chi2_dof->at(it) > L1Tk_maxChi2dof)
1426         continue;
1427 
1428       // use tight quality cut selection?
1429       /*
1430       if (matchtrk_nstub->at(it)==4) {
1431     if (std::abs(matchtrk_eta->at(it))<2.2 && matchtrk_consistency->at(it)>10) continue;
1432     else if (std::abs(matchtrk_eta->at(it))>2.2 && chi2dof>5.0) continue;
1433       }
1434       if (matchtrk_pt->at(it)>10.0 && chi2dof>5.0) continue;
1435       */
1436 
1437       // ----------------------------------------------------------------------------------------------------------------
1438       // more plots
1439 
1440       // fill matched track histograms
1441       h_match_tp_pt->Fill(tp_pt->at(it));
1442       if (tp_pt->at(it) < 8.0)
1443         h_match_tp_pt_L->Fill(tp_pt->at(it));
1444       else
1445         h_match_tp_pt_H->Fill(tp_pt->at(it));
1446       if (tp_pt->at(it) < 8.0 && std::abs(tp_eta->at(it)) < 1.0)
1447         h_match_tp_pt_LC->Fill(tp_pt->at(it));
1448 
1449       if (tp_pt->at(it) > TP_minPt) {
1450         h_match_tp_eta->Fill(tp_eta->at(it));
1451         h_match_tp_phi->Fill(tp_phi->at(it));
1452         h_match_tp_z0->Fill(tp_z0->at(it));
1453         h_match_tp_d0->Fill(tp_d0->at(it));
1454         h_match_tp_absd0->Fill(std::abs(tp_d0->at(it)));
1455         if (std::abs(tp_eta->at(it)) < 2.0)
1456           h_match_tp_absd0_eta2->Fill(std::abs(tp_d0->at(it)));
1457         if (std::abs(tp_eta->at(it)) < 2.0 && tp_pt->at(it) > 3.0)
1458           h_match_tp_absd0_eta2_pt3->Fill(std::abs(tp_d0->at(it)));
1459 
1460         if (std::abs(tp_eta->at(it)) < 1.0)
1461           n_match_eta1p0++;
1462         else if (std::abs(tp_eta->at(it)) < 1.75)
1463           n_match_eta1p75++;
1464         else
1465           n_match_eta2p5++;
1466 
1467         if (std::abs(tp_pt->at(it)) > 2.0)
1468           n_match_ptg2++;
1469         if (std::abs(tp_pt->at(it)) > 2.0 && std::abs(tp_pt->at(it)) < 8.0)
1470           n_match_pt2to8++;
1471         if (std::abs(tp_pt->at(it)) > 8.0)
1472           n_match_ptg8++;
1473         if (std::abs(tp_pt->at(it)) > 40.0)
1474           n_match_ptg40++;
1475 
1476         if (tp_pt->at(it) < 3.0)
1477           h_match_tp_eta_23->Fill(tp_eta->at(it));
1478         else if (tp_pt->at(it) < 5.0)
1479           h_match_tp_eta_35->Fill(tp_eta->at(it));
1480         else
1481           h_match_tp_eta_5->Fill(tp_eta->at(it));
1482 
1483         if (tp_pt->at(it) < 8.0) {
1484           h_match_tp_eta_L->Fill(tp_eta->at(it));
1485           h_match_tp_z0_L->Fill(tp_z0->at(it));
1486         } else {
1487           h_match_tp_z0_H->Fill(tp_z0->at(it));
1488           h_match_tp_eta_H->Fill(tp_eta->at(it));
1489         }
1490       }
1491 
1492       // for the following, only consider TPs with pt > TP_minPt
1493       if (tp_pt->at(it) < TP_minPt)
1494         continue;
1495 
1496       // fill nstub histograms
1497       h_match_trk_nstub->Fill(matchtrk_nstub->at(it));
1498       if (std::abs(tp_eta->at(it)) < 0.8)
1499         h_match_trk_nstub_C->Fill(matchtrk_nstub->at(it));
1500       else if (std::abs(tp_eta->at(it)) < 1.6 && std::abs(tp_eta->at(it)) >= 0.8)
1501         h_match_trk_nstub_I->Fill(matchtrk_nstub->at(it));
1502       else if (std::abs(tp_eta->at(it)) >= 1.6)
1503         h_match_trk_nstub_F->Fill(matchtrk_nstub->at(it));
1504 
1505       // ----------------------------------------------------------------------------------------------------------------
1506       // fill resolution histograms
1507 
1508       h_res_pt->Fill(matchtrk_pt->at(it) - tp_pt->at(it));
1509       h_res_ptRel->Fill((matchtrk_pt->at(it) - tp_pt->at(it)) / tp_pt->at(it));
1510       h_res_eta->Fill(matchtrk_eta->at(it) - tp_eta->at(it));
1511       h_res_phi->Fill(matchtrk_phi->at(it) - tp_phi->at(it));
1512       h_res_z0->Fill(matchtrk_z0->at(it) - tp_z0->at(it));
1513       if (matchtrk_d0->at(it) < 999.)
1514         h_res_d0->Fill(matchtrk_d0->at(it) - tp_d0->at(it));
1515 
1516       if (std::abs(tp_eta->at(it)) < 0.8)
1517         h_res_z0_C->Fill(matchtrk_z0->at(it) - tp_z0->at(it));
1518       else if (std::abs(tp_eta->at(it)) < 1.6 && std::abs(tp_eta->at(it)) >= 0.8)
1519         h_res_z0_I->Fill(matchtrk_z0->at(it) - tp_z0->at(it));
1520       else if (std::abs(tp_eta->at(it)) >= 1.6)
1521         h_res_z0_F->Fill(matchtrk_z0->at(it) - tp_z0->at(it));
1522 
1523       if (tp_pt->at(it) < 8.0) {
1524         h_res_z0_L->Fill(matchtrk_z0->at(it) - tp_z0->at(it));
1525         if (std::abs(tp_eta->at(it)) < 1.0)
1526           h_res_z0_C_L->Fill(matchtrk_z0->at(it) - tp_z0->at(it));
1527         else
1528           h_res_z0_F_L->Fill(matchtrk_z0->at(it) - tp_z0->at(it));
1529       } else {
1530         h_res_z0_H->Fill(matchtrk_z0->at(it) - tp_z0->at(it));
1531         if (std::abs(tp_eta->at(it)) < 1.0)
1532           h_res_z0_C_H->Fill(matchtrk_z0->at(it) - tp_z0->at(it));
1533         else
1534           h_res_z0_F_H->Fill(matchtrk_z0->at(it) - tp_z0->at(it));
1535       }
1536 
1537       if (matchtrk_d0->at(it) < 999.) {
1538         if (std::abs(tp_eta->at(it)) < 0.8)
1539           h_res_d0_C->Fill(matchtrk_d0->at(it) - tp_d0->at(it));
1540         else if (std::abs(tp_eta->at(it)) < 1.6 && std::abs(tp_eta->at(it)) >= 0.8)
1541           h_res_d0_I->Fill(matchtrk_d0->at(it) - tp_d0->at(it));
1542         else if (std::abs(tp_eta->at(it)) >= 1.6)
1543           h_res_d0_F->Fill(matchtrk_d0->at(it) - tp_d0->at(it));
1544 
1545         if (tp_pt->at(it) < 8.0) {
1546           h_res_d0_L->Fill(matchtrk_d0->at(it) - tp_d0->at(it));
1547           if (std::abs(tp_eta->at(it)) < 1.0)
1548             h_res_d0_C_L->Fill(matchtrk_d0->at(it) - tp_d0->at(it));
1549           else
1550             h_res_d0_F_L->Fill(matchtrk_d0->at(it) - tp_d0->at(it));
1551         } else {
1552           h_res_d0_H->Fill(matchtrk_d0->at(it) - tp_d0->at(it));
1553           if (std::abs(tp_eta->at(it)) < 1.0)
1554             h_res_d0_C_H->Fill(matchtrk_d0->at(it) - tp_d0->at(it));
1555           else
1556             h_res_d0_F_H->Fill(matchtrk_d0->at(it) - tp_d0->at(it));
1557         }
1558       }
1559 
1560       // ----------------------------------------------------------------------------------------------------------------
1561       // fill resolution vs. pt histograms
1562       for (int im = 0; im < nRANGE; im++) {
1563         if ((tp_pt->at(it) > (float)im * 5.0) && (tp_pt->at(it) < (float)(im + 1) * 5.0)) {
1564           h_resVsPt_pt[im]->Fill(matchtrk_pt->at(it) - tp_pt->at(it));
1565           h_resVsPt_ptRel[im]->Fill((matchtrk_pt->at(it) - tp_pt->at(it)) / tp_pt->at(it));
1566           h_resVsPt_eta[im]->Fill(matchtrk_eta->at(it) - tp_eta->at(it));
1567           h_resVsPt_phi[im]->Fill(matchtrk_phi->at(it) - tp_phi->at(it));
1568           h_resVsPt_z0[im]->Fill(matchtrk_z0->at(it) - tp_z0->at(it));
1569 
1570           h_absResVsPt_pt[im]->Fill(std::abs(matchtrk_pt->at(it) - tp_pt->at(it)));
1571           h_absResVsPt_ptRel[im]->Fill(std::abs((matchtrk_pt->at(it) - tp_pt->at(it))) / tp_pt->at(it));
1572           h_absResVsPt_z0[im]->Fill(std::abs(matchtrk_z0->at(it) - tp_z0->at(it)));
1573           h_absResVsPt_phi[im]->Fill(std::abs(matchtrk_phi->at(it) - tp_phi->at(it)));
1574           h_absResVsPt_eta[im]->Fill(std::abs(matchtrk_eta->at(it) - tp_eta->at(it)));
1575 
1576           if (std::abs(tp_eta->at(it)) < 0.8) {
1577             h_resVsPt_pt_C[im]->Fill(matchtrk_pt->at(it) - tp_pt->at(it));
1578             h_resVsPt_ptRel_C[im]->Fill((matchtrk_pt->at(it) - tp_pt->at(it)) / tp_pt->at(it));
1579             h_resVsPt_z0_C[im]->Fill(matchtrk_z0->at(it) - tp_z0->at(it));
1580             h_resVsPt_phi_C[im]->Fill(matchtrk_phi->at(it) - tp_phi->at(it));
1581           } else if (std::abs(tp_eta->at(it)) < 1.6 && std::abs(tp_eta->at(it)) >= 0.8) {
1582             h_resVsPt_pt_I[im]->Fill(matchtrk_pt->at(it) - tp_pt->at(it));
1583             h_resVsPt_ptRel_I[im]->Fill((matchtrk_pt->at(it) - tp_pt->at(it)) / tp_pt->at(it));
1584             h_resVsPt_z0_I[im]->Fill(matchtrk_z0->at(it) - tp_z0->at(it));
1585             h_resVsPt_phi_I[im]->Fill(matchtrk_phi->at(it) - tp_phi->at(it));
1586           } else if (std::abs(tp_eta->at(it)) >= 1.6) {
1587             h_resVsPt_pt_F[im]->Fill(matchtrk_pt->at(it) - tp_pt->at(it));
1588             h_resVsPt_ptRel_F[im]->Fill((matchtrk_pt->at(it) - tp_pt->at(it)) / tp_pt->at(it));
1589             h_resVsPt_z0_F[im]->Fill(matchtrk_z0->at(it) - tp_z0->at(it));
1590             h_resVsPt_phi_F[im]->Fill(matchtrk_phi->at(it) - tp_phi->at(it));
1591           }
1592 
1593           if (matchtrk_d0->at(it) < 999) {
1594             h_resVsPt_d0[im]->Fill(matchtrk_d0->at(it) - tp_d0->at(it));
1595             h_absResVsPt_d0[im]->Fill(std::abs(matchtrk_d0->at(it) - tp_d0->at(it)));
1596           }
1597         }
1598       }
1599 
1600       for (int im = 3; im < nRANGE_L + 3; im++) {
1601         if ((tp_pt->at(it) > (float)im * 0.5) && (tp_pt->at(it) <= (float)(im + 1) * 0.5)) {
1602           h_absResVsPt_pt_L[im - 3]->Fill(std::abs(matchtrk_pt->at(it) - tp_pt->at(it)));
1603           h_absResVsPt_ptRel_L[im - 3]->Fill(std::abs((matchtrk_pt->at(it) - tp_pt->at(it))) / tp_pt->at(it));
1604           h_absResVsPt_z0_L[im - 3]->Fill(std::abs(matchtrk_z0->at(it) - tp_z0->at(it)));
1605           h_absResVsPt_phi_L[im - 3]->Fill(std::abs(matchtrk_phi->at(it) - tp_phi->at(it)));
1606           h_absResVsPt_eta_L[im - 3]->Fill(std::abs(matchtrk_eta->at(it) - tp_eta->at(it)));
1607           h_absResVsPt_d0_L[im - 3]->Fill(std::abs(matchtrk_d0->at(it) - tp_d0->at(it)));
1608         }
1609       }
1610 
1611       // ----------------------------------------------------------------------------------------------------------------
1612       // fill resolution vs. eta histograms
1613       for (int im = 0; im < nETARANGE; im++) {
1614         if ((std::abs(tp_eta->at(it)) > (float)im * 0.1) && (std::abs(tp_eta->at(it)) < (float)(im + 1) * 0.1)) {
1615           h_resVsEta_ptRel[im]->Fill((matchtrk_pt->at(it) - tp_pt->at(it)) / tp_pt->at(it));
1616           h_resVsEta_eta[im]->Fill(matchtrk_eta->at(it) - tp_eta->at(it));
1617           h_resVsEta_phi[im]->Fill(matchtrk_phi->at(it) - tp_phi->at(it));
1618           h_resVsEta_z0[im]->Fill(matchtrk_z0->at(it) - tp_z0->at(it));
1619 
1620           h_absResVsEta_ptRel[im]->Fill(std::abs((matchtrk_pt->at(it) - tp_pt->at(it))) / tp_pt->at(it));
1621           h_absResVsEta_eta[im]->Fill(std::abs(matchtrk_eta->at(it) - tp_eta->at(it)));
1622           h_absResVsEta_phi[im]->Fill(std::abs(matchtrk_phi->at(it) - tp_phi->at(it)));
1623           h_absResVsEta_z0[im]->Fill(std::abs(matchtrk_z0->at(it) - tp_z0->at(it)));
1624 
1625           if (tp_pt->at(it) < 8.0) {
1626             h_resVsEta_ptRel_L[im]->Fill((matchtrk_pt->at(it) - tp_pt->at(it)) / tp_pt->at(it));
1627             h_resVsEta_eta_L[im]->Fill(matchtrk_eta->at(it) - tp_eta->at(it));
1628             h_resVsEta_z0_L[im]->Fill(matchtrk_z0->at(it) - tp_z0->at(it));
1629             h_resVsEta_phi_L[im]->Fill(matchtrk_phi->at(it) - tp_phi->at(it));
1630 
1631             h_absResVsEta_ptRel_L[im]->Fill(std::abs((matchtrk_pt->at(it) - tp_pt->at(it))) / tp_pt->at(it));
1632             h_absResVsEta_eta_L[im]->Fill(std::abs(matchtrk_eta->at(it) - tp_eta->at(it)));
1633             h_absResVsEta_phi_L[im]->Fill(std::abs(matchtrk_phi->at(it) - tp_phi->at(it)));
1634             h_absResVsEta_z0_L[im]->Fill(std::abs(matchtrk_z0->at(it) - tp_z0->at(it)));
1635 
1636             if (matchtrk_d0->at(it) < 999) {
1637               h_resVsEta_d0_L[im]->Fill(matchtrk_d0->at(it) - tp_d0->at(it));
1638               h_absResVsEta_d0_L[im]->Fill(std::abs(matchtrk_d0->at(it) - tp_d0->at(it)));
1639             }
1640           } else {
1641             h_resVsEta_ptRel_H[im]->Fill((matchtrk_pt->at(it) - tp_pt->at(it)) / tp_pt->at(it));
1642             h_resVsEta_eta_H[im]->Fill(matchtrk_eta->at(it) - tp_eta->at(it));
1643             h_resVsEta_z0_H[im]->Fill(matchtrk_z0->at(it) - tp_z0->at(it));
1644             h_resVsEta_phi_H[im]->Fill(matchtrk_phi->at(it) - tp_phi->at(it));
1645 
1646             h_absResVsEta_ptRel_H[im]->Fill(std::abs((matchtrk_pt->at(it) - tp_pt->at(it))) / tp_pt->at(it));
1647             h_absResVsEta_eta_H[im]->Fill(std::abs(matchtrk_eta->at(it) - tp_eta->at(it)));
1648             h_absResVsEta_phi_H[im]->Fill(std::abs(matchtrk_phi->at(it) - tp_phi->at(it)));
1649             h_absResVsEta_z0_H[im]->Fill(std::abs(matchtrk_z0->at(it) - tp_z0->at(it)));
1650 
1651             if (matchtrk_d0->at(it) < 999) {
1652               h_resVsEta_d0_H[im]->Fill(matchtrk_d0->at(it) - tp_d0->at(it));
1653               h_absResVsEta_d0_H[im]->Fill(std::abs(matchtrk_d0->at(it) - tp_d0->at(it)));
1654             }
1655           }
1656 
1657           if (matchtrk_d0->at(it) < 999) {
1658             h_resVsEta_d0[im]->Fill(matchtrk_d0->at(it) - tp_d0->at(it));
1659             h_absResVsEta_d0[im]->Fill(std::abs(matchtrk_d0->at(it) - tp_d0->at(it)));
1660           }
1661         }
1662       }
1663 
1664       // ----------------------------------------------------------------------------------------------------------------
1665       // fill resolution vs. phi histograms
1666       for (int im = 0; im < nPHIRANGE; im++) {
1667         if ((tp_phi->at(it) > (float)im * 0.2 - 3.2) && (tp_phi->at(it) < (float)(im + 1) * 0.2 - 3.2)) {
1668           h_absResVsPhi_pt[im]->Fill(std::abs(matchtrk_pt->at(it) - tp_pt->at(it)));
1669           h_absResVsPhi_ptRel[im]->Fill(std::abs((matchtrk_pt->at(it) - tp_pt->at(it))) / tp_pt->at(it));
1670         }
1671       }
1672 
1673     }  // end of matched track loop
1674 
1675   }  // end of event loop
1676 
1677   // ----------------------------------------------------------------------------------------------------------------
1678 
1679   // ----------------------------------------------------------------------------------------------------------------
1680   // Post-event-loop histogram normalization(s)
1681 
1682   // Normalize tracklet efficiency by each eta slice
1683   double maxBinContents;
1684   double etaRes = trackletEffMaxEta / trackletEffEtaBins;
1685   for (double etaBin = 0; etaBin < trackletEffEtaBins;
1686        etaBin++) {  //loop through eta bin values (constants defined with relevant hist defs)
1687     maxBinContents = 0;
1688     std::vector<double> binContents = {};
1689     for (int layer = 0; layer < numLayers; layer++) {
1690       binContents.push_back(h_trk_tracklet_hits->GetBinContent(etaBin + 1, layer + 1));
1691       maxBinContents = std::max(maxBinContents, binContents.back());
1692     }
1693     float binWeight;
1694     for (int layer = 0; layer < numLayers; layer++) {
1695       binWeight = (maxBinContents == 0) ? 0 : binContents.at(layer) / maxBinContents;
1696       h_trk_tracklet_eff->Fill((etaBin + 0.5) * etaRes, layer, binWeight);
1697     }
1698   }
1699 
1700   // Adjust tracklet hits and efficiency histograms for aesthetics
1701 
1702   // ----------------------------------------------------------------------------------------------------------------
1703   // 2D plots
1704   // ----------------------------------------------------------------------------------------------------------------
1705 
1706   TH1F* h2_resVsPt_pt = new TH1F("resVsPt2_pt", ";Tracking particle p_{T} [GeV]; p_{T} resolution [GeV]", 20, 0, 100);
1707   TH1F* h2_resVsPt_pt_C =
1708       new TH1F("resVsPt2_pt_C", ";Tracking particle p_{T} [GeV]; p_{T} resolution [GeV]", 20, 0, 100);
1709   TH1F* h2_resVsPt_pt_I =
1710       new TH1F("resVsPt2_pt_I", ";Tracking particle p_{T} [GeV]; p_{T} resolution [GeV]", 20, 0, 100);
1711   TH1F* h2_resVsPt_pt_F =
1712       new TH1F("resVsPt2_pt_F", ";Tracking particle p_{T} [GeV]; p_{T} resolution [GeV]", 20, 0, 100);
1713 
1714   TH1F* h2_resVsPt_ptRel =
1715       new TH1F("resVsPt2_ptRel", ";Tracking particle p_{T} [GeV]; p_{T} resolution / p_{T}", 20, 0, 100);
1716   TH1F* h2_resVsPt_ptRel_C =
1717       new TH1F("resVsPt2_ptRel_C", ";Tracking particle p_{T} [GeV]; p_{T} resolution / p_{T}", 20, 0, 100);
1718   TH1F* h2_resVsPt_ptRel_I =
1719       new TH1F("resVsPt2_ptRel_I", ";Tracking particle p_{T} [GeV]; p_{T} resolution / p_{T}", 20, 0, 100);
1720   TH1F* h2_resVsPt_ptRel_F =
1721       new TH1F("resVsPt2_ptRel_F", ";Tracking particle p_{T} [GeV]; p_{T} resolution / p_{T}", 20, 0, 100);
1722 
1723   TH1F* h2_mresVsPt_pt =
1724       new TH1F("mresVsPt2_pt", ";Tracking particle p_{T} [GeV]; Mean(p_{T} residual) [GeV]", 20, 0, 100);
1725   TH1F* h2_mresVsPt_pt_C =
1726       new TH1F("mresVsPt2_pt_C", ";Tracking particle p_{T} [GeV]; Mean(p_{T} residual) [GeV]", 20, 0, 100);
1727   TH1F* h2_mresVsPt_pt_I =
1728       new TH1F("mresVsPt2_pt_I", ";Tracking particle p_{T} [GeV]; Mean(p_{T} residual) [GeV]", 20, 0, 100);
1729   TH1F* h2_mresVsPt_pt_F =
1730       new TH1F("mresVsPt2_pt_F", ";Tracking particle p_{T} [GeV]; Mean(p_{T} residual) [GeV]", 20, 0, 100);
1731 
1732   TH1F* h2_resVsPt_z0 = new TH1F("resVsPt2_z0", ";Tracking particle p_{T} [GeV]; z_{0} resolution [cm]", 20, 0, 100);
1733   TH1F* h2_resVsPt_z0_C =
1734       new TH1F("resVsPt2_z0_C", ";Tracking particle p_{T} [GeV]; z_{0} resolution [cm]", 20, 0, 100);
1735   TH1F* h2_resVsPt_z0_I =
1736       new TH1F("resVsPt2_z0_I", ";Tracking particle p_{T} [GeV]; z_{0} resolution [cm]", 20, 0, 100);
1737   TH1F* h2_resVsPt_z0_F =
1738       new TH1F("resVsPt2_z0_F", ";Tracking particle p_{T} [GeV]; z_{0} resolution [cm]", 20, 0, 100);
1739 
1740   TH1F* h2_resVsPt_phi = new TH1F("resVsPt2_phi", ";Tracking particle p_{T} [GeV]; #phi resolution [rad]", 20, 0, 100);
1741   TH1F* h2_resVsPt_phi_C =
1742       new TH1F("resVsPt2_phi_C", ";Tracking particle p_{T} [GeV]; #phi resolution [rad]", 20, 0, 100);
1743   TH1F* h2_resVsPt_phi_I =
1744       new TH1F("resVsPt2_phi_I", ";Tracking particle p_{T} [GeV]; #phi resolution [rad]", 20, 0, 100);
1745   TH1F* h2_resVsPt_phi_F =
1746       new TH1F("resVsPt2_phi_F", ";Tracking particle p_{T} [GeV]; #phi resolution [rad]", 20, 0, 100);
1747 
1748   TH1F* h2_resVsPt_eta = new TH1F("resVsPt2_eta", ";Tracking particle p_{T} [GeV]; #eta resolution", 20, 0, 100);
1749 
1750   TH1F* h2_resVsPt_d0 = new TH1F("resVsPt2_d0", ";Tracking particle p_{T} [GeV]; d_{0} resolution [cm]", 20, 0, 100);
1751 
1752   TH1F* h2_resVsPt_pt_68 =
1753       new TH1F("resVsPt2_pt_68", ";Tracking particle p_{T} [GeV]; p_{T} resolution [GeV]", 20, 0, 100);
1754   TH1F* h2_resVsPt_ptRel_68 =
1755       new TH1F("resVsPt2_ptRel_68", ";Tracking particle p_{T} [GeV]; p_{T} resolution / p_{T}", 20, 0, 100);
1756   TH1F* h2_resVsPt_z0_68 =
1757       new TH1F("resVsPt2_z0_68", ";Tracking particle p_{T} [GeV]; z_{0} resolution [cm]", 20, 0, 100);
1758   TH1F* h2_resVsPt_phi_68 =
1759       new TH1F("resVsPt2_phi_68", ";Tracking particle p_{T} [GeV]; #phi resolution [rad]", 20, 0, 100);
1760   TH1F* h2_resVsPt_eta_68 = new TH1F("resVsPt2_eta_68", ";Tracking particle p_{T} [GeV]; #eta resolution", 20, 0, 100);
1761   TH1F* h2_resVsPt_d0_68 =
1762       new TH1F("resVsPt2_d0_68", ";Tracking particle p_{T} [GeV]; d_{0} resolution [cm]", 20, 0, 100);
1763 
1764   TH1F* h2_resVsPt_pt_90 =
1765       new TH1F("resVsPt2_pt_90", ";Tracking particle p_{T} [GeV]; p_{T} resolution [GeV]", 20, 0, 100);
1766   TH1F* h2_resVsPt_ptRel_90 =
1767       new TH1F("resVsPt2_ptRel_90", ";Tracking particle p_{T} [GeV]; p_{T} resolution / p_{T}", 20, 0, 100);
1768   TH1F* h2_resVsPt_z0_90 =
1769       new TH1F("resVsPt2_z0_90", ";Tracking particle p_{T} [GeV]; z_{0} resolution [cm]", 20, 0, 100);
1770   TH1F* h2_resVsPt_phi_90 =
1771       new TH1F("resVsPt2_phi_90", ";Tracking particle p_{T} [GeV]; #phi resolution [rad]", 20, 0, 100);
1772   TH1F* h2_resVsPt_eta_90 = new TH1F("resVsPt2_eta_90", ";Tracking particle p_{T} [GeV]; #eta resolution", 20, 0, 100);
1773   TH1F* h2_resVsPt_d0_90 =
1774       new TH1F("resVsPt2_d0_90", ";Tracking particle p_{T} [GeV]; d_{0} resolution [cm]", 20, 0, 100);
1775 
1776   TH1F* h2_resVsPt_pt_99 =
1777       new TH1F("resVsPt2_pt_99", ";Tracking particle p_{T} [GeV]; p_{T} resolution [GeV]", 20, 0, 100);
1778   TH1F* h2_resVsPt_ptRel_99 =
1779       new TH1F("resVsPt2_ptRel_99", ";Tracking particle p_{T} [GeV]; p_{T} resolution / p_{T}", 20, 0, 100);
1780   TH1F* h2_resVsPt_z0_99 =
1781       new TH1F("resVsPt2_z0_99", ";Tracking particle p_{T} [GeV]; z_{0} resolution [cm]", 20, 0, 100);
1782   TH1F* h2_resVsPt_phi_99 =
1783       new TH1F("resVsPt2_phi_99", ";Tracking particle p_{T} [GeV]; #phi resolution [rad]", 20, 0, 100);
1784   TH1F* h2_resVsPt_eta_99 = new TH1F("resVsPt2_eta_99", ";Tracking particle p_{T} [GeV]; #eta resolution", 20, 0, 100);
1785   TH1F* h2_resVsPt_d0_99 =
1786       new TH1F("resVsPt2_d0_99", ";Tracking particle p_{T} [GeV]; d_{0} resolution [cm]", 20, 0, 100);
1787 
1788   TH1F* h2_resVsPt_pt_L_68 =
1789       new TH1F("resVsPt2_pt_L_68", ";Tracking particle p_{T} [GeV]; p_{T} resolution [GeV]", nRANGE_L, pt_resmin, 8);
1790   TH1F* h2_resVsPt_ptRel_L_68 = new TH1F(
1791       "resVsPt2_ptRel_L_68", ";Tracking particle p_{T} [GeV]; p_{T} resolution / p_{T}", nRANGE_L, pt_resmin, 8);
1792   TH1F* h2_resVsPt_z0_L_68 =
1793       new TH1F("resVsPt2_z0_L_68", ";Tracking particle p_{T} [GeV]; z_{0} resolution [cm]", nRANGE_L, pt_resmin, 8);
1794   TH1F* h2_resVsPt_phi_L_68 =
1795       new TH1F("resVsPt2_phi_L_68", ";Tracking particle p_{T} [GeV]; #phi resolution [rad]", nRANGE_L, pt_resmin, 8);
1796   TH1F* h2_resVsPt_eta_L_68 =
1797       new TH1F("resVsPt2_eta_L_68", ";Tracking particle p_{T} [GeV]; #eta resolution", nRANGE_L, pt_resmin, 8);
1798   TH1F* h2_resVsPt_d0_L_68 =
1799       new TH1F("resVsPt2_d0_L_68", ";Tracking particle p_{T} [GeV]; d_{0} resolution [cm]", nRANGE_L, pt_resmin, 8);
1800 
1801   TH1F* h2_resVsPt_pt_L_90 =
1802       new TH1F("resVsPt2_pt_L_90", ";Tracking particle p_{T} [GeV]; p_{T} resolution [GeV]", nRANGE_L, pt_resmin, 8);
1803   TH1F* h2_resVsPt_ptRel_L_90 = new TH1F(
1804       "resVsPt2_ptRel_L_90", ";Tracking particle p_{T} [GeV]; p_{T} resolution / p_{T}", nRANGE_L, pt_resmin, 8);
1805   TH1F* h2_resVsPt_z0_L_90 =
1806       new TH1F("resVsPt2_z0_L_90", ";Tracking particle p_{T} [GeV]; z_{0} resolution [cm]", nRANGE_L, pt_resmin, 8);
1807   TH1F* h2_resVsPt_phi_L_90 =
1808       new TH1F("resVsPt2_phi_L_90", ";Tracking particle p_{T} [GeV]; #phi resolution [rad]", nRANGE_L, pt_resmin, 8);
1809   TH1F* h2_resVsPt_eta_L_90 =
1810       new TH1F("resVsPt2_eta_L_90", ";Tracking particle p_{T} [GeV]; #eta resolution", nRANGE_L, pt_resmin, 8);
1811   TH1F* h2_resVsPt_d0_L_90 =
1812       new TH1F("resVsPt2_d0_L_90", ";Tracking particle p_{T} [GeV]; d_{0} resolution [cm]", nRANGE_L, pt_resmin, 8);
1813 
1814   TH1F* h2_resVsPt_pt_L_99 =
1815       new TH1F("resVsPt2_pt_L_99", ";Tracking particle p_{T} [GeV]; p_{T} resolution [cm]", nRANGE_L, pt_resmin, 8);
1816   TH1F* h2_resVsPt_ptRel_L_99 = new TH1F(
1817       "resVsPt2_ptRel_L_99", ";Tracking particle p_{T} [GeV]; p_{T} resolution / p_{T}", nRANGE_L, pt_resmin, 8);
1818   TH1F* h2_resVsPt_z0_L_99 =
1819       new TH1F("resVsPt2_z0_L_99", ";Tracking particle p_{T} [GeV]; z_{0} resolution [cm]", nRANGE_L, pt_resmin, 8);
1820   TH1F* h2_resVsPt_phi_L_99 =
1821       new TH1F("resVsPt2_phi_L_99", ";Tracking particle p_{T} [GeV]; #phi resolution [rad]", nRANGE_L, pt_resmin, 8);
1822   TH1F* h2_resVsPt_eta_L_99 =
1823       new TH1F("resVsPt2_eta_L_99", ";Tracking particle p_{T} [GeV]; #eta resolution", nRANGE_L, pt_resmin, 8);
1824   TH1F* h2_resVsPt_d0_L_99 =
1825       new TH1F("resVsPt2_d0_L_99", ";Tracking particle p_{T} [GeV]; d_{0} resolution [cm]", nRANGE_L, pt_resmin, 8);
1826 
1827   for (int i = 0; i < nRANGE; i++) {
1828     // set bin content and error
1829     h2_resVsPt_pt->SetBinContent(i + 1, h_resVsPt_pt[i]->GetRMS());
1830     h2_resVsPt_pt->SetBinError(i + 1, h_resVsPt_pt[i]->GetRMSError());
1831     h2_resVsPt_pt_C->SetBinContent(i + 1, h_resVsPt_pt_C[i]->GetRMS());
1832     h2_resVsPt_pt_C->SetBinError(i + 1, h_resVsPt_pt_C[i]->GetRMSError());
1833     h2_resVsPt_pt_I->SetBinContent(i + 1, h_resVsPt_pt_I[i]->GetRMS());
1834     h2_resVsPt_pt_I->SetBinError(i + 1, h_resVsPt_pt_I[i]->GetRMSError());
1835     h2_resVsPt_pt_F->SetBinContent(i + 1, h_resVsPt_pt_F[i]->GetRMS());
1836     h2_resVsPt_pt_F->SetBinError(i + 1, h_resVsPt_pt_F[i]->GetRMSError());
1837 
1838     h2_resVsPt_ptRel->SetBinContent(i + 1, h_resVsPt_ptRel[i]->GetRMS());
1839     h2_resVsPt_ptRel->SetBinError(i + 1, h_resVsPt_ptRel[i]->GetRMSError());
1840     h2_resVsPt_ptRel_C->SetBinContent(i + 1, h_resVsPt_ptRel_C[i]->GetRMS());
1841     h2_resVsPt_ptRel_C->SetBinError(i + 1, h_resVsPt_ptRel_C[i]->GetRMSError());
1842     h2_resVsPt_ptRel_I->SetBinContent(i + 1, h_resVsPt_ptRel_I[i]->GetRMS());
1843     h2_resVsPt_ptRel_I->SetBinError(i + 1, h_resVsPt_ptRel_I[i]->GetRMSError());
1844     h2_resVsPt_ptRel_F->SetBinContent(i + 1, h_resVsPt_ptRel_F[i]->GetRMS());
1845     h2_resVsPt_ptRel_F->SetBinError(i + 1, h_resVsPt_ptRel_F[i]->GetRMSError());
1846 
1847     h2_mresVsPt_pt->SetBinContent(i + 1, h_resVsPt_pt[i]->GetMean());
1848     h2_mresVsPt_pt->SetBinError(i + 1, h_resVsPt_pt[i]->GetMeanError());
1849     h2_mresVsPt_pt_C->SetBinContent(i + 1, h_resVsPt_pt_C[i]->GetMean());
1850     h2_mresVsPt_pt_C->SetBinError(i + 1, h_resVsPt_pt_C[i]->GetMeanError());
1851     h2_mresVsPt_pt_I->SetBinContent(i + 1, h_resVsPt_pt_I[i]->GetMean());
1852     h2_mresVsPt_pt_I->SetBinError(i + 1, h_resVsPt_pt_I[i]->GetMeanError());
1853     h2_mresVsPt_pt_F->SetBinContent(i + 1, h_resVsPt_pt_F[i]->GetMean());
1854     h2_mresVsPt_pt_F->SetBinError(i + 1, h_resVsPt_pt_F[i]->GetMeanError());
1855 
1856     h2_resVsPt_z0->SetBinContent(i + 1, h_resVsPt_z0[i]->GetRMS());
1857     h2_resVsPt_z0->SetBinError(i + 1, h_resVsPt_z0[i]->GetRMSError());
1858     h2_resVsPt_z0_C->SetBinContent(i + 1, h_resVsPt_z0_C[i]->GetRMS());
1859     h2_resVsPt_z0_C->SetBinError(i + 1, h_resVsPt_z0_C[i]->GetRMSError());
1860     h2_resVsPt_z0_I->SetBinContent(i + 1, h_resVsPt_z0_I[i]->GetRMS());
1861     h2_resVsPt_z0_I->SetBinError(i + 1, h_resVsPt_z0_I[i]->GetRMSError());
1862     h2_resVsPt_z0_F->SetBinContent(i + 1, h_resVsPt_z0_F[i]->GetRMS());
1863     h2_resVsPt_z0_F->SetBinError(i + 1, h_resVsPt_z0_F[i]->GetRMSError());
1864 
1865     h2_resVsPt_phi->SetBinContent(i + 1, h_resVsPt_phi[i]->GetRMS());
1866     h2_resVsPt_phi->SetBinError(i + 1, h_resVsPt_phi[i]->GetRMSError());
1867     h2_resVsPt_phi_C->SetBinContent(i + 1, h_resVsPt_phi_C[i]->GetRMS());
1868     h2_resVsPt_phi_C->SetBinError(i + 1, h_resVsPt_phi_C[i]->GetRMSError());
1869     h2_resVsPt_phi_I->SetBinContent(i + 1, h_resVsPt_phi_I[i]->GetRMS());
1870     h2_resVsPt_phi_I->SetBinError(i + 1, h_resVsPt_phi_I[i]->GetRMSError());
1871     h2_resVsPt_phi_F->SetBinContent(i + 1, h_resVsPt_phi_F[i]->GetRMS());
1872     h2_resVsPt_phi_F->SetBinError(i + 1, h_resVsPt_phi_F[i]->GetRMSError());
1873 
1874     h2_resVsPt_eta->SetBinContent(i + 1, h_resVsPt_eta[i]->GetRMS());
1875     h2_resVsPt_eta->SetBinError(i + 1, h_resVsPt_eta[i]->GetRMSError());
1876 
1877     h2_resVsPt_d0->SetBinContent(i + 1, h_resVsPt_d0[i]->GetRMS());
1878     h2_resVsPt_d0->SetBinError(i + 1, h_resVsPt_d0[i]->GetRMSError());
1879 
1880     h2_resVsPt_pt_68->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_pt[i], 0.68));
1881     h2_resVsPt_pt_90->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_pt[i], 0.90));
1882     h2_resVsPt_pt_99->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_pt[i], 0.99));
1883 
1884     h2_resVsPt_ptRel_68->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_ptRel[i], 0.68));
1885     h2_resVsPt_ptRel_90->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_ptRel[i], 0.90));
1886     h2_resVsPt_ptRel_99->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_ptRel[i], 0.99));
1887 
1888     h2_resVsPt_eta_68->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_eta[i], 0.68));
1889     h2_resVsPt_eta_90->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_eta[i], 0.90));
1890     h2_resVsPt_eta_99->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_eta[i], 0.99));
1891 
1892     h2_resVsPt_z0_68->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_z0[i], 0.68));
1893     h2_resVsPt_z0_90->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_z0[i], 0.90));
1894     h2_resVsPt_z0_99->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_z0[i], 0.99));
1895 
1896     h2_resVsPt_phi_68->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_phi[i], 0.68));
1897     h2_resVsPt_phi_90->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_phi[i], 0.90));
1898     h2_resVsPt_phi_99->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_phi[i], 0.99));
1899 
1900     h2_resVsPt_d0_68->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_d0[i], 0.68));
1901     h2_resVsPt_d0_90->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_d0[i], 0.90));
1902     h2_resVsPt_d0_99->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_d0[i], 0.99));
1903   }
1904 
1905   for (int i = 0; i < nRANGE_L; i++) {
1906     h2_resVsPt_pt_L_68->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_pt_L[i], 0.68));
1907     h2_resVsPt_pt_L_90->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_pt_L[i], 0.90));
1908     h2_resVsPt_pt_L_99->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_pt_L[i], 0.99));
1909 
1910     h2_resVsPt_ptRel_L_68->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_ptRel_L[i], 0.68));
1911     h2_resVsPt_ptRel_L_90->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_ptRel_L[i], 0.90));
1912     h2_resVsPt_ptRel_L_99->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_ptRel_L[i], 0.99));
1913 
1914     h2_resVsPt_eta_L_68->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_eta_L[i], 0.68));
1915     h2_resVsPt_eta_L_90->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_eta_L[i], 0.90));
1916     h2_resVsPt_eta_L_99->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_eta_L[i], 0.99));
1917 
1918     h2_resVsPt_z0_L_68->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_z0_L[i], 0.68));
1919     h2_resVsPt_z0_L_90->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_z0_L[i], 0.90));
1920     h2_resVsPt_z0_L_99->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_z0_L[i], 0.99));
1921 
1922     h2_resVsPt_phi_L_68->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_phi_L[i], 0.68));
1923     h2_resVsPt_phi_L_90->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_phi_L[i], 0.90));
1924     h2_resVsPt_phi_L_99->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_phi_L[i], 0.99));
1925 
1926     h2_resVsPt_d0_L_68->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_d0_L[i], 0.68));
1927     h2_resVsPt_d0_L_90->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_d0_L[i], 0.90));
1928     h2_resVsPt_d0_L_99->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPt_d0_L[i], 0.99));
1929   }
1930 
1931   // resolution vs. eta histograms
1932   TH1F* h2_resVsEta_eta =
1933       new TH1F("resVsEta_eta", ";Tracking particle |#eta|; #eta resolution", nETARANGE, 0, eta_resmax);
1934   TH1F* h2_resVsEta_eta_L =
1935       new TH1F("resVsEta_eta_L", ";Tracking particle |#eta|; #eta resolution", nETARANGE, 0, eta_resmax);
1936   TH1F* h2_resVsEta_eta_H =
1937       new TH1F("resVsEta_eta_H", ";Tracking particle |#eta|; #eta resolution", nETARANGE, 0, eta_resmax);
1938 
1939   TH1F* h2_resVsEta_z0 =
1940       new TH1F("resVsEta_z0", ";Tracking particle |#eta|; z_{0} resolution [cm]", nETARANGE, 0, eta_resmax);
1941   TH1F* h2_resVsEta_z0_L =
1942       new TH1F("resVsEta_z0_L", ";Tracking particle |#eta|; z_{0} resolution [cm]", nETARANGE, 0, eta_resmax);
1943   TH1F* h2_resVsEta_z0_H =
1944       new TH1F("resVsEta_z0_H", ";Tracking particle |#eta|; z_{0} resolution [cm]", nETARANGE, 0, eta_resmax);
1945 
1946   TH1F* h2_resVsEta_phi =
1947       new TH1F("resVsEta_phi", ";Tracking particle |#eta|; #phi resolution [rad]", nETARANGE, 0, eta_resmax);
1948   TH1F* h2_resVsEta_phi_L =
1949       new TH1F("resVsEta_phi_L", ";Tracking particle |#eta|; #phi resolution [rad]", nETARANGE, 0, eta_resmax);
1950   TH1F* h2_resVsEta_phi_H =
1951       new TH1F("resVsEta_phi_H", ";Tracking particle |#eta|; #phi resolution [rad]", nETARANGE, 0, eta_resmax);
1952 
1953   TH1F* h2_resVsEta_ptRel =
1954       new TH1F("resVsEta_ptRel", ";Tracking particle |#eta|; p_{T} resolution / p_{T}", nETARANGE, 0, eta_resmax);
1955   TH1F* h2_resVsEta_ptRel_L =
1956       new TH1F("resVsEta_ptRel_L", ";Tracking particle |#eta|; p_{T} resolution / p_{T}", nETARANGE, 0, eta_resmax);
1957   TH1F* h2_resVsEta_ptRel_H =
1958       new TH1F("resVsEta_ptRel_H", ";Tracking particle |#eta|; p_{T} resolution / p_{T}", nETARANGE, 0, eta_resmax);
1959 
1960   TH1F* h2_resVsEta_d0 =
1961       new TH1F("resVsEta_d0", ";Tracking particle |#eta|; d_{0} resolution [cm]", nETARANGE, 0, eta_resmax);
1962   TH1F* h2_resVsEta_d0_L =
1963       new TH1F("resVsEta_d0_L", ";Tracking particle |#eta|; d_{0} resolution [cm]", nETARANGE, 0, eta_resmax);
1964   TH1F* h2_resVsEta_d0_H =
1965       new TH1F("resVsEta_d0_H", ";Tracking particle |#eta|; d_{0} resolution [cm]", nETARANGE, 0, eta_resmax);
1966 
1967   // mean of residuals
1968   TH1F* h2_mresVsEta_eta =
1969       new TH1F("mresVsEta_eta", ";Tracking particle |#eta|; Mean(#eta residual)", nETARANGE, 0, eta_resmax);
1970   TH1F* h2_mresVsEta_z0 =
1971       new TH1F("mresVsEta_z0", ";Tracking particle |#eta|; Mean(z_{0} residual) [cm]", nETARANGE, 0, eta_resmax);
1972   TH1F* h2_mresVsEta_phi =
1973       new TH1F("mresVsEta_phi", ";Tracking particle |#eta|; Mean(phi residual) [rad]", nETARANGE, 0, eta_resmax);
1974   TH1F* h2_mresVsEta_ptRel =
1975       new TH1F("mresVsEta_ptRel", ";Tracking particle |#eta|; Mean(ptrel residual)", nETARANGE, 0, eta_resmax);
1976 
1977   // 68 / 90 / 99% residuals
1978   TH1F* h2_resVsEta_eta_68 =
1979       new TH1F("resVsEta_eta_68", ";Tracking particle |#eta|; #eta resolution", nETARANGE, 0, eta_resmax);
1980   TH1F* h2_resVsEta_eta_90 =
1981       new TH1F("resVsEta_eta_90", ";Tracking particle |#eta|; #eta resolution", nETARANGE, 0, eta_resmax);
1982   TH1F* h2_resVsEta_eta_99 =
1983       new TH1F("resVsEta_eta_99", ";Tracking particle |#eta|; #eta resolution", nETARANGE, 0, eta_resmax);
1984 
1985   TH1F* h2_resVsEta_z0_68 =
1986       new TH1F("resVsEta_z0_68", ";Tracking particle |#eta|; z_{0} resolution [cm]", nETARANGE, 0, eta_resmax);
1987   TH1F* h2_resVsEta_z0_90 =
1988       new TH1F("resVsEta_z0_90", ";Tracking particle |#eta|; z_{0} resolution [cm]", nETARANGE, 0, eta_resmax);
1989   TH1F* h2_resVsEta_z0_99 =
1990       new TH1F("resVsEta_z0_99", ";Tracking particle |#eta|; z_{0} resolution [cm]", nETARANGE, 0, eta_resmax);
1991 
1992   TH1F* h2_resVsEta_phi_68 =
1993       new TH1F("resVsEta_phi_68", ";Tracking particle |#eta|; #phi resolution [rad]", nETARANGE, 0, eta_resmax);
1994   TH1F* h2_resVsEta_phi_90 =
1995       new TH1F("resVsEta_phi_90", ";Tracking particle |#eta|; #phi resolution [rad]", nETARANGE, 0, eta_resmax);
1996   TH1F* h2_resVsEta_phi_99 =
1997       new TH1F("resVsEta_phi_99", ";Tracking particle |#eta|; #phi resolution [rad]", nETARANGE, 0, eta_resmax);
1998 
1999   TH1F* h2_resVsEta_ptRel_68 =
2000       new TH1F("resVsEta_ptRel_68", ";Tracking particle |#eta|; p_{T} resolution / p_{T}", nETARANGE, 0, eta_resmax);
2001   TH1F* h2_resVsEta_ptRel_90 =
2002       new TH1F("resVsEta_ptRel_90", ";Tracking particle |#eta|; p_{T} resolution / p_{T}", nETARANGE, 0, eta_resmax);
2003   TH1F* h2_resVsEta_ptRel_99 =
2004       new TH1F("resVsEta_ptRel_99", ";Tracking particle |#eta|; p_{T} resolution / p_{T}", nETARANGE, 0, eta_resmax);
2005 
2006   TH1F* h2_resVsEta_d0_68 =
2007       new TH1F("resVsEta_d0_68", ";Tracking particle |#eta|; d_{0} resolution [cm]", nETARANGE, 0, eta_resmax);
2008   TH1F* h2_resVsEta_d0_90 =
2009       new TH1F("resVsEta_d0_90", ";Tracking particle |#eta|; d_{0} resolution [cm]", nETARANGE, 0, eta_resmax);
2010   TH1F* h2_resVsEta_d0_99 =
2011       new TH1F("resVsEta_d0_99", ";Tracking particle |#eta|; d_{0} resolution [cm]", nETARANGE, 0, eta_resmax);
2012 
2013   TH1F* h2_resVsEta_eta_L_68 =
2014       new TH1F("resVsEta_eta_L_68", ";Tracking particle |#eta|; #eta resolution", nETARANGE, 0, eta_resmax);
2015   TH1F* h2_resVsEta_z0_L_68 =
2016       new TH1F("resVsEta_z0_L_68", ";Tracking particle |#eta|; z_{0} resolution [cm]", nETARANGE, 0, eta_resmax);
2017   TH1F* h2_resVsEta_phi_L_68 =
2018       new TH1F("resVsEta_phi_L_68", ";Tracking particle |#eta|; #phi resolution [rad]", nETARANGE, 0, eta_resmax);
2019   TH1F* h2_resVsEta_ptRel_L_68 =
2020       new TH1F("resVsEta_ptRel_L_68", ";Tracking particle |#eta|; p_{T} resolution / p_{T}", nETARANGE, 0, eta_resmax);
2021   TH1F* h2_resVsEta_d0_L_68 =
2022       new TH1F("resVsEta_d0_L_68", ";Tracking particle |#eta|; d_{0} resolution [cm]", nETARANGE, 0, eta_resmax);
2023 
2024   TH1F* h2_resVsEta_eta_L_90 =
2025       new TH1F("resVsEta_eta_L_90", ";Tracking particle |#eta|; #eta resolution", nETARANGE, 0, eta_resmax);
2026   TH1F* h2_resVsEta_z0_L_90 =
2027       new TH1F("resVsEta_z0_L_90", ";Tracking particle |#eta|; z_{0} resolution [cm]", nETARANGE, 0, eta_resmax);
2028   TH1F* h2_resVsEta_phi_L_90 =
2029       new TH1F("resVsEta_phi_L_90", ";Tracking particle |#eta|; #phi resolution [rad]", nETARANGE, 0, eta_resmax);
2030   TH1F* h2_resVsEta_ptRel_L_90 =
2031       new TH1F("resVsEta_ptRel_L_90", ";Tracking particle |#eta|; p_{T} resolution / p_{T}", nETARANGE, 0, eta_resmax);
2032   TH1F* h2_resVsEta_d0_L_90 =
2033       new TH1F("resVsEta_d0_L_90", ";Tracking particle |#eta|; d_{0} resolution [cm]", nETARANGE, 0, eta_resmax);
2034 
2035   TH1F* h2_resVsEta_eta_L_99 =
2036       new TH1F("resVsEta_eta_L_99", ";Tracking particle |#eta|; #eta resolution", nETARANGE, 0, eta_resmax);
2037   TH1F* h2_resVsEta_z0_L_99 =
2038       new TH1F("resVsEta_z0_L_99", ";Tracking particle |#eta|; z_{0} resolution [cm]", nETARANGE, 0, eta_resmax);
2039   TH1F* h2_resVsEta_phi_L_99 =
2040       new TH1F("resVsEta_phi_L_99", ";Tracking particle |#eta|; #phi resolution [rad]", nETARANGE, 0, eta_resmax);
2041   TH1F* h2_resVsEta_ptRel_L_99 =
2042       new TH1F("resVsEta_ptRel_L_99", ";Tracking particle |#eta|; p_{T} resolution / p_{T}", nETARANGE, 0, eta_resmax);
2043   TH1F* h2_resVsEta_d0_L_99 =
2044       new TH1F("resVsEta_d0_L_99", ";Tracking particle |#eta|; d_{0} resolution [cm]", nETARANGE, 0, eta_resmax);
2045 
2046   TH1F* h2_resVsEta_eta_H_68 =
2047       new TH1F("resVsEta_eta_H_68", ";Tracking particle |#eta|; #eta resolution", nETARANGE, 0, eta_resmax);
2048   TH1F* h2_resVsEta_z0_H_68 =
2049       new TH1F("resVsEta_z0_H_68", ";Tracking particle |#eta|; z_{0} resolution [cm]", nETARANGE, 0, eta_resmax);
2050   TH1F* h2_resVsEta_phi_H_68 =
2051       new TH1F("resVsEta_phi_H_68", ";Tracking particle |#eta|; #phi resolution [rad]", nETARANGE, 0, eta_resmax);
2052   TH1F* h2_resVsEta_ptRel_H_68 =
2053       new TH1F("resVsEta_ptRel_H_68", ";Tracking particle |#eta|; p_{T} resolution / p_{T}", nETARANGE, 0, eta_resmax);
2054   TH1F* h2_resVsEta_d0_H_68 =
2055       new TH1F("resVsEta_d0_H_68", ";Tracking particle |#eta|; d_{0} resolution [cm]", nETARANGE, 0, eta_resmax);
2056 
2057   TH1F* h2_resVsEta_eta_H_90 =
2058       new TH1F("resVsEta_eta_H_90", ";Tracking particle |#eta|; #eta resolution", nETARANGE, 0, eta_resmax);
2059   TH1F* h2_resVsEta_z0_H_90 =
2060       new TH1F("resVsEta_z0_H_90", ";Tracking particle |#eta|; z_{0} resolution [cm]", nETARANGE, 0, eta_resmax);
2061   TH1F* h2_resVsEta_phi_H_90 =
2062       new TH1F("resVsEta_phi_H_90", ";Tracking particle |#eta|; #phi resolution [rad]", nETARANGE, 0, eta_resmax);
2063   TH1F* h2_resVsEta_ptRel_H_90 =
2064       new TH1F("resVsEta_ptRel_H_90", ";Tracking particle |#eta|; p_{T} resolution / p_{T}", nETARANGE, 0, eta_resmax);
2065   TH1F* h2_resVsEta_d0_H_90 =
2066       new TH1F("resVsEta_d0_H_90", ";Tracking particle |#eta|; d_{0} resolution [cm]", nETARANGE, 0, eta_resmax);
2067 
2068   TH1F* h2_resVsEta_eta_H_99 =
2069       new TH1F("resVsEta_eta_H_99", ";Tracking particle |#eta|; #eta resolution", nETARANGE, 0, eta_resmax);
2070   TH1F* h2_resVsEta_z0_H_99 =
2071       new TH1F("resVsEta_z0_H_99", ";Tracking particle |#eta|; z_{0} resolution [cm]", nETARANGE, 0, eta_resmax);
2072   TH1F* h2_resVsEta_phi_H_99 =
2073       new TH1F("resVsEta_phi_H_99", ";Tracking particle |#eta|; #phi resolution [rad]", nETARANGE, 0, eta_resmax);
2074   TH1F* h2_resVsEta_ptRel_H_99 =
2075       new TH1F("resVsEta_ptRel_H_99", ";Tracking particle |#eta|; p_{T} resolution / p_{T}", nETARANGE, 0, eta_resmax);
2076   TH1F* h2_resVsEta_d0_H_99 =
2077       new TH1F("resVsEta_d0_H_99", ";Tracking particle |#eta|; d_{0} resolution [cm]", nETARANGE, 0, eta_resmax);
2078 
2079   // resolution vs. eta histograms (gaussian fit)
2080   TH1F* h3_resVsEta_eta_L = new TH1F("resVsEta_eta_L_gaus", ";|#eta|; #sigma(#eta)", nETARANGE, 0, eta_resmax);
2081   TH1F* h3_resVsEta_eta_H = new TH1F("resVsEta_eta_H_gaus", ";|#eta|; #sigma(#eta)", nETARANGE, 0, eta_resmax);
2082 
2083   TH1F* h3_resVsEta_z0_L = new TH1F("resVsEta_z0_L_gaus", ";|#eta|; #sigma(z_{0}) [cm]", nETARANGE, 0, eta_resmax);
2084   TH1F* h3_resVsEta_z0_H = new TH1F("resVsEta_z0_H_gaus", ";|#eta|; #sigma(z_{0}) [cm]", nETARANGE, 0, eta_resmax);
2085 
2086   TH1F* h3_resVsEta_phi_L = new TH1F("resVsEta_phi_L_gaus", ";|#eta|; #sigma(#phi) [rad]", nETARANGE, 0, eta_resmax);
2087   TH1F* h3_resVsEta_phi_H = new TH1F("resVsEta_phi_H_gaus", ";|#eta|; #sigma(#phi) [rad]", nETARANGE, 0, eta_resmax);
2088 
2089   TH1F* h3_resVsEta_ptRel_L =
2090       new TH1F("resVsEta_ptRel_L_gaus", ";|#eta|; #sigma(p_{T}) / p_{T}", nETARANGE, 0, eta_resmax);
2091   TH1F* h3_resVsEta_ptRel_H =
2092       new TH1F("resVsEta_ptRel_H_gaus", ";|#eta|; #sigma(p_{T}) / p_{T}", nETARANGE, 0, eta_resmax);
2093 
2094   gSystem->mkdir("FitResults");
2095   TString fitdir = "FitResults/";
2096 
2097   for (int i = 0; i < nETARANGE; i++) {
2098     // set bin content and error
2099     h2_resVsEta_eta->SetBinContent(i + 1, h_resVsEta_eta[i]->GetRMS());
2100     h2_resVsEta_eta->SetBinError(i + 1, h_resVsEta_eta[i]->GetRMSError());
2101     h2_resVsEta_eta_L->SetBinContent(i + 1, h_resVsEta_eta_L[i]->GetRMS());
2102     h2_resVsEta_eta_L->SetBinError(i + 1, h_resVsEta_eta_L[i]->GetRMSError());
2103     h2_resVsEta_eta_H->SetBinContent(i + 1, h_resVsEta_eta_H[i]->GetRMS());
2104     h2_resVsEta_eta_H->SetBinError(i + 1, h_resVsEta_eta_H[i]->GetRMSError());
2105 
2106     h2_resVsEta_z0->SetBinContent(i + 1, h_resVsEta_z0[i]->GetRMS());
2107     h2_resVsEta_z0->SetBinError(i + 1, h_resVsEta_z0[i]->GetRMSError());
2108     h2_resVsEta_z0_L->SetBinContent(i + 1, h_resVsEta_z0_L[i]->GetRMS());
2109     h2_resVsEta_z0_L->SetBinError(i + 1, h_resVsEta_z0_L[i]->GetRMSError());
2110     h2_resVsEta_z0_H->SetBinContent(i + 1, h_resVsEta_z0_H[i]->GetRMS());
2111     h2_resVsEta_z0_H->SetBinError(i + 1, h_resVsEta_z0_H[i]->GetRMSError());
2112 
2113     h2_resVsEta_phi->SetBinContent(i + 1, h_resVsEta_phi[i]->GetRMS());
2114     h2_resVsEta_phi->SetBinError(i + 1, h_resVsEta_phi[i]->GetRMSError());
2115     h2_resVsEta_phi_L->SetBinContent(i + 1, h_resVsEta_phi_L[i]->GetRMS());
2116     h2_resVsEta_phi_L->SetBinError(i + 1, h_resVsEta_phi_L[i]->GetRMSError());
2117     h2_resVsEta_phi_H->SetBinContent(i + 1, h_resVsEta_phi_H[i]->GetRMS());
2118     h2_resVsEta_phi_H->SetBinError(i + 1, h_resVsEta_phi_H[i]->GetRMSError());
2119 
2120     h2_resVsEta_ptRel->SetBinContent(i + 1, h_resVsEta_ptRel[i]->GetRMS());
2121     h2_resVsEta_ptRel->SetBinError(i + 1, h_resVsEta_ptRel[i]->GetRMSError());
2122     h2_resVsEta_ptRel_L->SetBinContent(i + 1, h_resVsEta_ptRel_L[i]->GetRMS());
2123     h2_resVsEta_ptRel_L->SetBinError(i + 1, h_resVsEta_ptRel_L[i]->GetRMSError());
2124     h2_resVsEta_ptRel_H->SetBinContent(i + 1, h_resVsEta_ptRel_H[i]->GetRMS());
2125     h2_resVsEta_ptRel_H->SetBinError(i + 1, h_resVsEta_ptRel_H[i]->GetRMSError());
2126 
2127     h2_mresVsEta_eta->SetBinContent(i + 1, h_resVsEta_eta[i]->GetMean());
2128     h2_mresVsEta_eta->SetBinError(i + 1, h_resVsEta_eta[i]->GetMeanError());
2129     h2_mresVsEta_z0->SetBinContent(i + 1, h_resVsEta_z0[i]->GetMean());
2130     h2_mresVsEta_z0->SetBinError(i + 1, h_resVsEta_z0[i]->GetMeanError());
2131     h2_mresVsEta_phi->SetBinContent(i + 1, h_resVsEta_phi[i]->GetMean());
2132     h2_mresVsEta_phi->SetBinError(i + 1, h_resVsEta_phi[i]->GetMeanError());
2133     h2_mresVsEta_ptRel->SetBinContent(i + 1, h_resVsEta_ptRel[i]->GetMean());
2134     h2_mresVsEta_ptRel->SetBinError(i + 1, h_resVsEta_ptRel[i]->GetMeanError());
2135 
2136     h2_resVsEta_d0->SetBinContent(i + 1, h_resVsEta_d0[i]->GetRMS());
2137     h2_resVsEta_d0->SetBinError(i + 1, h_resVsEta_d0[i]->GetRMSError());
2138     h2_resVsEta_d0_L->SetBinContent(i + 1, h_resVsEta_d0_L[i]->GetRMS());
2139     h2_resVsEta_d0_L->SetBinError(i + 1, h_resVsEta_d0_L[i]->GetRMSError());
2140     h2_resVsEta_d0_H->SetBinContent(i + 1, h_resVsEta_d0_H[i]->GetRMS());
2141     h2_resVsEta_d0_H->SetBinError(i + 1, h_resVsEta_d0_H[i]->GetRMSError());
2142 
2143     h2_resVsEta_eta_68->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_eta[i], 0.68));
2144     h2_resVsEta_eta_90->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_eta[i], 0.90));
2145     h2_resVsEta_eta_99->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_eta[i], 0.99));
2146 
2147     h2_resVsEta_z0_68->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_z0[i], 0.68));
2148     h2_resVsEta_z0_90->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_z0[i], 0.90));
2149     h2_resVsEta_z0_99->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_z0[i], 0.99));
2150 
2151     h2_resVsEta_phi_68->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_phi[i], 0.68));
2152     h2_resVsEta_phi_90->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_phi[i], 0.90));
2153     h2_resVsEta_phi_99->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_phi[i], 0.99));
2154 
2155     h2_resVsEta_ptRel_68->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_ptRel[i], 0.68));
2156     h2_resVsEta_ptRel_90->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_ptRel[i], 0.90));
2157     h2_resVsEta_ptRel_99->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_ptRel[i], 0.99));
2158 
2159     h2_resVsEta_d0_68->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_d0[i], 0.68));
2160     h2_resVsEta_d0_90->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_d0[i], 0.90));
2161     h2_resVsEta_d0_99->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_d0[i], 0.99));
2162 
2163     h2_resVsEta_eta_L_68->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_eta_L[i], 0.68));
2164     h2_resVsEta_z0_L_68->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_z0_L[i], 0.68));
2165     h2_resVsEta_phi_L_68->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_phi_L[i], 0.68));
2166     h2_resVsEta_ptRel_L_68->SetBinContent(i + 1,
2167                                           getIntervalContainingFractionOfEntries(h_absResVsEta_ptRel_L[i], 0.68));
2168     h2_resVsEta_d0_L_68->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_d0_L[i], 0.68));
2169 
2170     h2_resVsEta_eta_L_90->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_eta_L[i], 0.90));
2171     h2_resVsEta_z0_L_90->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_z0_L[i], 0.90));
2172     h2_resVsEta_phi_L_90->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_phi_L[i], 0.90));
2173     h2_resVsEta_ptRel_L_90->SetBinContent(i + 1,
2174                                           getIntervalContainingFractionOfEntries(h_absResVsEta_ptRel_L[i], 0.90));
2175     h2_resVsEta_d0_L_90->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_d0_L[i], 0.90));
2176 
2177     h2_resVsEta_eta_L_99->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_eta_L[i], 0.99));
2178     h2_resVsEta_z0_L_99->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_z0_L[i], 0.99));
2179     h2_resVsEta_phi_L_99->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_phi_L[i], 0.99));
2180     h2_resVsEta_ptRel_L_99->SetBinContent(i + 1,
2181                                           getIntervalContainingFractionOfEntries(h_absResVsEta_ptRel_L[i], 0.99));
2182     h2_resVsEta_d0_L_99->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_d0_L[i], 0.99));
2183 
2184     h2_resVsEta_eta_H_68->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_eta_H[i], 0.68));
2185     h2_resVsEta_z0_H_68->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_z0_H[i], 0.68));
2186     h2_resVsEta_phi_H_68->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_phi_H[i], 0.68));
2187     h2_resVsEta_ptRel_H_68->SetBinContent(i + 1,
2188                                           getIntervalContainingFractionOfEntries(h_absResVsEta_ptRel_H[i], 0.68));
2189     h2_resVsEta_d0_H_68->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_d0_H[i], 0.68));
2190 
2191     h2_resVsEta_eta_H_90->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_eta_H[i], 0.90));
2192     h2_resVsEta_z0_H_90->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_z0_H[i], 0.90));
2193     h2_resVsEta_phi_H_90->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_phi_H[i], 0.90));
2194     h2_resVsEta_ptRel_H_90->SetBinContent(i + 1,
2195                                           getIntervalContainingFractionOfEntries(h_absResVsEta_ptRel_H[i], 0.90));
2196     h2_resVsEta_d0_H_90->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_d0_H[i], 0.90));
2197 
2198     h2_resVsEta_eta_H_99->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_eta_H[i], 0.99));
2199     h2_resVsEta_z0_H_99->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_z0_H[i], 0.99));
2200     h2_resVsEta_phi_H_99->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_phi_H[i], 0.99));
2201     h2_resVsEta_ptRel_H_99->SetBinContent(i + 1,
2202                                           getIntervalContainingFractionOfEntries(h_absResVsEta_ptRel_H[i], 0.99));
2203     h2_resVsEta_d0_H_99->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsEta_d0_H[i], 0.99));
2204 
2205     // ---------------------------------------------------------------------------------------------------
2206     // gaussian fit instead
2207     // ---------------------------------------------------------------------------------------------------
2208 
2209     if (doGausFit) {
2210       TCanvas cfit;
2211       char text[500];
2212 
2213       float sigma = 0;
2214       float esigma = 0;
2215       TF1* fit;
2216 
2217       float rms = 0;
2218       float erms = 0;
2219 
2220       fit = new TF1("fit", "gaus", -0.01, 0.01);
2221       h_resVsEta_eta_L[i]->Fit("fit", "R");
2222       sigma = fit->GetParameter(2);
2223       esigma = fit->GetParError(2);
2224       rms = h_resVsEta_eta_L[i]->GetRMS();
2225       erms = h_resVsEta_eta_L[i]->GetRMSError();
2226       h3_resVsEta_eta_L->SetBinContent(i + 1, sigma);
2227       h3_resVsEta_eta_L->SetBinError(i + 1, esigma);
2228       h_resVsEta_eta_L[i]->Draw();
2229       sprintf(text, "RMS: %.4f +/- %.4f", rms, erms);
2230       mySmallText(0.2, 0.86, 1, text);
2231       sprintf(text, "Fit: %.4f +/- %.4f", sigma, esigma);
2232       mySmallText(0.2, 0.8, 2, text);
2233       sprintf(text, "p_{T} < 5 GeV");
2234       mySmallText(0.2, 0.7, 1, text);
2235       cfit.SaveAs(fitdir + "resVsEta_eta_L_" + etarange[i] + ".pdf");
2236       delete fit;
2237 
2238       fit = new TF1("fit", "gaus", -0.01, 0.01);
2239       h_resVsEta_eta_H[i]->Fit("fit", "R");
2240       sigma = fit->GetParameter(2);
2241       esigma = fit->GetParError(2);
2242       rms = h_resVsEta_eta_H[i]->GetRMS();
2243       erms = h_resVsEta_eta_H[i]->GetRMSError();
2244       h3_resVsEta_eta_H->SetBinContent(i + 1, sigma);
2245       h3_resVsEta_eta_H->SetBinError(i + 1, esigma);
2246       h_resVsEta_eta_H[i]->Draw();
2247       sprintf(text, "RMS: %.4f +/- %.4f", rms, erms);
2248       mySmallText(0.2, 0.86, 1, text);
2249       sprintf(text, "Fit: %.4f +/- %.4f", sigma, esigma);
2250       mySmallText(0.2, 0.8, 2, text);
2251       sprintf(text, "p_{T} > 15 GeV");
2252       mySmallText(0.2, 0.7, 1, text);
2253       cfit.SaveAs(fitdir + "resVsEta_eta_H_" + etarange[i] + ".pdf");
2254       delete fit;
2255 
2256       fit = new TF1("fit", "gaus", -1, 1);
2257       h_resVsEta_z0_L[i]->Fit("fit", "R");
2258       sigma = fit->GetParameter(2);
2259       esigma = fit->GetParError(2);
2260       rms = h_resVsEta_z0_L[i]->GetRMS();
2261       erms = h_resVsEta_z0_L[i]->GetRMSError();
2262       h3_resVsEta_z0_L->SetBinContent(i + 1, sigma);
2263       h3_resVsEta_z0_L->SetBinError(i + 1, esigma);
2264       h_resVsEta_z0_L[i]->Draw();
2265       sprintf(text, "RMS: %.4f +/- %.4f", rms, erms);
2266       mySmallText(0.2, 0.86, 1, text);
2267       sprintf(text, "Fit: %.4f +/- %.4f", sigma, esigma);
2268       mySmallText(0.2, 0.8, 2, text);
2269       sprintf(text, "p_{T} < 5 GeV");
2270       mySmallText(0.2, 0.7, 1, text);
2271       cfit.SaveAs(fitdir + "resVsEta_z0_L_" + etarange[i] + ".pdf");
2272       delete fit;
2273 
2274       fit = new TF1("fit", "gaus", -1, 1);
2275       h_resVsEta_z0_H[i]->Fit("fit", "R");
2276       sigma = fit->GetParameter(2);
2277       esigma = fit->GetParError(2);
2278       rms = h_resVsEta_z0_H[i]->GetRMS();
2279       erms = h_resVsEta_z0_H[i]->GetRMSError();
2280       h3_resVsEta_z0_H->SetBinContent(i + 1, sigma);
2281       h3_resVsEta_z0_H->SetBinError(i + 1, esigma);
2282       h_resVsEta_z0_H[i]->Draw();
2283       sprintf(text, "RMS: %.4f +/- %.4f", rms, erms);
2284       mySmallText(0.2, 0.86, 1, text);
2285       sprintf(text, "Fit: %.4f +/- %.4f", sigma, esigma);
2286       mySmallText(0.2, 0.8, 2, text);
2287       sprintf(text, "p_{T} > 15 GeV");
2288       mySmallText(0.2, 0.7, 1, text);
2289       cfit.SaveAs(fitdir + "resVsEta_z0_H_" + etarange[i] + ".pdf");
2290       delete fit;
2291 
2292       fit = new TF1("fit", "gaus", -0.005, 0.005);
2293       h_resVsEta_phi_L[i]->Fit("fit", "R");
2294       sigma = fit->GetParameter(2);
2295       esigma = fit->GetParError(2);
2296       rms = h_resVsEta_phi_L[i]->GetRMS();
2297       erms = h_resVsEta_phi_L[i]->GetRMSError();
2298       h3_resVsEta_phi_L->SetBinContent(i + 1, sigma);
2299       h3_resVsEta_phi_L->SetBinError(i + 1, esigma);
2300       h_resVsEta_phi_L[i]->Draw();
2301       sprintf(text, "RMS: %.4f +/- %.4f", rms, erms);
2302       mySmallText(0.2, 0.86, 1, text);
2303       sprintf(text, "Fit: %.4f +/- %.4f", sigma, esigma);
2304       mySmallText(0.2, 0.8, 2, text);
2305       sprintf(text, "p_{T} < 5 GeV");
2306       mySmallText(0.2, 0.7, 1, text);
2307       cfit.SaveAs(fitdir + "resVsEta_phi_L_" + etarange[i] + ".pdf");
2308       delete fit;
2309 
2310       fit = new TF1("fit", "gaus", -0.005, 0.005);
2311       h_resVsEta_phi_H[i]->Fit("fit", "R");
2312       sigma = fit->GetParameter(2);
2313       esigma = fit->GetParError(2);
2314       rms = h_resVsEta_phi_H[i]->GetRMS();
2315       erms = h_resVsEta_phi_H[i]->GetRMSError();
2316       h3_resVsEta_phi_H->SetBinContent(i + 1, sigma);
2317       h3_resVsEta_phi_H->SetBinError(i + 1, esigma);
2318       h_resVsEta_phi_H[i]->Draw();
2319       sprintf(text, "RMS: %.4f +/- %.4f", rms, erms);
2320       mySmallText(0.2, 0.86, 1, text);
2321       sprintf(text, "Fit: %.4f +/- %.4f", sigma, esigma);
2322       mySmallText(0.2, 0.8, 2, text);
2323       sprintf(text, "p_{T} > 15 GeV");
2324       mySmallText(0.2, 0.7, 1, text);
2325       cfit.SaveAs(fitdir + "resVsEta_phi_H_" + etarange[i] + ".pdf");
2326       delete fit;
2327 
2328       fit = new TF1("fit", "gaus", -0.5, 0.5);
2329       h_resVsEta_ptRel_L[i]->Fit("fit", "R");
2330       sigma = fit->GetParameter(2);
2331       esigma = fit->GetParError(2);
2332       rms = h_resVsEta_ptRel_L[i]->GetRMS();
2333       erms = h_resVsEta_ptRel_L[i]->GetRMSError();
2334       h3_resVsEta_ptRel_L->SetBinContent(i + 1, sigma);
2335       h3_resVsEta_ptRel_L->SetBinError(i + 1, esigma);
2336       h_resVsEta_ptRel_L[i]->Draw();
2337       sprintf(text, "RMS: %.4f +/- %.4f", rms, erms);
2338       mySmallText(0.2, 0.86, 1, text);
2339       sprintf(text, "Fit: %.4f +/- %.4f", sigma, esigma);
2340       mySmallText(0.2, 0.8, 2, text);
2341       sprintf(text, "p_{T} < 5 GeV");
2342       mySmallText(0.2, 0.7, 1, text);
2343       cfit.SaveAs(fitdir + "resVsEta_ptRel_L_" + etarange[i] + ".pdf");
2344       delete fit;
2345 
2346       fit = new TF1("fit", "gaus", -0.5, 0.5);
2347       h_resVsEta_ptRel_H[i]->Fit("fit", "R");
2348       sigma = fit->GetParameter(2);
2349       esigma = fit->GetParError(2);
2350       rms = h_resVsEta_ptRel_H[i]->GetRMS();
2351       erms = h_resVsEta_ptRel_H[i]->GetRMSError();
2352       h3_resVsEta_ptRel_H->SetBinContent(i + 1, sigma);
2353       h3_resVsEta_ptRel_H->SetBinError(i + 1, esigma);
2354       h_resVsEta_ptRel_H[i]->Draw();
2355       sprintf(text, "RMS: %.4f +/- %.4f", rms, erms);
2356       mySmallText(0.2, 0.86, 1, text);
2357       sprintf(text, "Fit: %.4f +/- %.4f", sigma, esigma);
2358       mySmallText(0.2, 0.8, 2, text);
2359       sprintf(text, "p_{T} > 15 GeV");
2360       mySmallText(0.2, 0.7, 1, text);
2361       cfit.SaveAs(fitdir + "resVsEta_ptRel_H_" + etarange[i] + ".pdf");
2362       delete fit;
2363 
2364     }  //end doGausFit
2365   }
2366 
2367   TH1F* h2_resVsPhi_pt_68 =
2368       new TH1F("resVsPhi2_pt_68", ";Tracking particle #phi; p_{T} resolution [GeV]", 32, -3.2, 3.2);
2369   TH1F* h2_resVsPhi_pt_90 =
2370       new TH1F("resVsPhi2_pt_90", ";Tracking particle #phi; p_{T} resolution [GeV]", 32, -3.2, 3.2);
2371   TH1F* h2_resVsPhi_pt_99 =
2372       new TH1F("resVsPhi2_pt_99", ";Tracking particle #phi; p_{T} resolution [GeV]", 32, -3.2, 3.2);
2373   TH1F* h2_resVsPhi_ptRel_68 =
2374       new TH1F("resVsPhi2_ptRel_68", ";Tracking particle #phi; p_{T} resolution / p_{T}", 32, -3.2, 3.2);
2375   TH1F* h2_resVsPhi_ptRel_90 =
2376       new TH1F("resVsPhi2_ptRel_90", ";Tracking particle #phi; p_{T} resolution / p_{T}", 32, -3.2, 3.2);
2377   TH1F* h2_resVsPhi_ptRel_99 =
2378       new TH1F("resVsPhi2_ptRel_99", ";Tracking particle #phi; p_{T} resolution / p_{T}", 32, -3.2, 3.2);
2379 
2380   for (int i = 0; i < nPHIRANGE; i++) {
2381     h2_resVsPhi_pt_68->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPhi_pt[i], 0.68));
2382     h2_resVsPhi_pt_90->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPhi_pt[i], 0.90));
2383     h2_resVsPhi_pt_99->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPhi_pt[i], 0.99));
2384 
2385     h2_resVsPhi_ptRel_68->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPhi_ptRel[i], 0.68));
2386     h2_resVsPhi_ptRel_90->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPhi_ptRel[i], 0.90));
2387     h2_resVsPhi_ptRel_99->SetBinContent(i + 1, getIntervalContainingFractionOfEntries(h_absResVsPhi_ptRel[i], 0.99));
2388   }
2389 
2390   // set minimum to zero
2391   h2_resVsPt_pt->SetMinimum(0);
2392   h2_resVsPt_pt_C->SetMinimum(0);
2393   h2_resVsPt_pt_I->SetMinimum(0);
2394   h2_resVsPt_pt_F->SetMinimum(0);
2395 
2396   h2_resVsPt_ptRel->SetMinimum(0);
2397   h2_resVsPt_ptRel_C->SetMinimum(0);
2398   h2_resVsPt_ptRel_I->SetMinimum(0);
2399   h2_resVsPt_ptRel_F->SetMinimum(0);
2400 
2401   h2_resVsPt_z0->SetMinimum(0);
2402   h2_resVsPt_z0_C->SetMinimum(0);
2403   h2_resVsPt_z0_I->SetMinimum(0);
2404   h2_resVsPt_z0_F->SetMinimum(0);
2405 
2406   h2_resVsPt_phi->SetMinimum(0);
2407   h2_resVsPt_phi_C->SetMinimum(0);
2408   h2_resVsPt_phi_I->SetMinimum(0);
2409   h2_resVsPt_phi_F->SetMinimum(0);
2410 
2411   h2_resVsPt_eta->SetMinimum(0);
2412 
2413   h2_resVsEta_eta->SetMinimum(0);
2414   h2_resVsEta_eta_L->SetMinimum(0);
2415   h2_resVsEta_eta_H->SetMinimum(0);
2416 
2417   h2_resVsEta_z0->SetMinimum(0);
2418   h2_resVsEta_z0_L->SetMinimum(0);
2419   h2_resVsEta_z0_H->SetMinimum(0);
2420 
2421   h2_resVsEta_phi->SetMinimum(0);
2422   h2_resVsEta_phi_L->SetMinimum(0);
2423   h2_resVsEta_phi_H->SetMinimum(0);
2424 
2425   h2_resVsEta_ptRel->SetMinimum(0);
2426   h2_resVsEta_ptRel_L->SetMinimum(0);
2427   h2_resVsEta_ptRel_H->SetMinimum(0);
2428 
2429   h2_resVsPt_d0->SetMinimum(0);
2430   h2_resVsEta_d0->SetMinimum(0);
2431   h2_resVsEta_d0_L->SetMinimum(0);
2432   h2_resVsEta_d0_H->SetMinimum(0);
2433 
2434   // -------------------------------------------------------------------------------------------
2435   // output file for histograms
2436   // -------------------------------------------------------------------------------------------
2437 
2438   if (TP_select_pdgid != 0) {
2439     char pdgidtxt[500];
2440     sprintf(pdgidtxt, "_pdgid%i", TP_select_pdgid);
2441     type = type + pdgidtxt;
2442   } else if (TP_select_injet == 1)
2443     type = type + "_injet";
2444   else if (TP_select_injet == 2)
2445     type = type + "_injet_highpt";
2446   else if (TP_select_injet == 3)
2447     type = type + "_injet_vhighpt";
2448 
2449   if (TP_select_eventid != 0)
2450     type = type + "_wpu";
2451 
2452   if (useTightCuts)
2453     type = type + "_tight";
2454   if (useDeadRegion)
2455     type = type + "_dead";
2456 
2457   if (L1Tk_seed != 0) {
2458     char seedtxt[500];
2459     sprintf(seedtxt, "_seed%i", L1Tk_seed);
2460     type = type + seedtxt;
2461   }
2462 
2463   if (TP_minPt > 2.0) {
2464     char pttxt[500];
2465     sprintf(pttxt, "_pt%.0f", TP_minPt);
2466     type = type + pttxt;
2467   }
2468 
2469   TFile* fout;
2470   if (doLooseMatch)
2471     fout = new TFile("output_looseMatch_" + type + treeName + ".root", "recreate");
2472   else
2473     fout = new TFile(type_dir + "output_" + type + treeName + ".root", "recreate");
2474 
2475   // -------------------------------------------------------------------------------------------
2476   // draw and save plots
2477   // -------------------------------------------------------------------------------------------
2478 
2479   char ctxt[500];
2480   TCanvas c;
2481 
2482   gSystem->mkdir("TrkPlots");
2483   TString DIR = "TrkPlots/";
2484 
2485   // plots overlaying 68, 90, 99% confidence levels]
2486 
2487   float max_eta_ptRel = 0.2;
2488   float max_pt_ptRel = 0.2;
2489   float max_pt_pt = 20;
2490   float max_z0 = 2.0;
2491   float max_phi = 0.01;
2492   float max_eta = 0.03;
2493 
2494   if (type.Contains("El")) {
2495     max_pt_ptRel = 1.0;
2496     max_eta_ptRel = 1.0;
2497     max_phi = 0.1;
2498   }
2499 
2500   // makeResidualIntervalPlot will save the individual plots to the root file
2501   makeResidualIntervalPlot(
2502       type, DIR, "resVsPt_ptRel", h2_resVsPt_ptRel_68, h2_resVsPt_ptRel_90, h2_resVsPt_ptRel_99, 0, max_pt_ptRel);
2503   makeResidualIntervalPlot(type, DIR, "resVsPt_pt", h2_resVsPt_pt_68, h2_resVsPt_pt_90, h2_resVsPt_pt_99, 0, max_pt_pt);
2504   makeResidualIntervalPlot(type, DIR, "resVsPt_z0", h2_resVsPt_z0_68, h2_resVsPt_z0_90, h2_resVsPt_z0_99, 0, max_z0);
2505   makeResidualIntervalPlot(
2506       type, DIR, "resVsPt_phi", h2_resVsPt_phi_68, h2_resVsPt_phi_90, h2_resVsPt_phi_99, 0, max_phi);
2507   makeResidualIntervalPlot(
2508       type, DIR, "resVsPt_eta", h2_resVsPt_eta_68, h2_resVsPt_eta_90, h2_resVsPt_eta_99, 0, max_eta);
2509   //makeResidualIntervalPlot( type, DIR, "resVsPt_d0", h2_resVsPt_d0_68, h2_resVsPt_d0_90, h2_resVsPt_d0_99, 0, 0.02 );
2510 
2511   makeResidualIntervalPlot(type,
2512                            DIR,
2513                            "resVsPt_L_ptRel",
2514                            h2_resVsPt_ptRel_L_68,
2515                            h2_resVsPt_ptRel_L_90,
2516                            h2_resVsPt_ptRel_L_99,
2517                            0,
2518                            max_pt_ptRel);
2519   makeResidualIntervalPlot(type, DIR, "resVsPt_L_pt", h2_resVsPt_pt_L_68, h2_resVsPt_pt_L_90, h2_resVsPt_pt_L_99, 0, 4);
2520   makeResidualIntervalPlot(
2521       type, DIR, "resVsPt_L_z0", h2_resVsPt_z0_L_68, h2_resVsPt_z0_L_90, h2_resVsPt_z0_L_99, 0, max_z0);
2522   makeResidualIntervalPlot(
2523       type, DIR, "resVsPt_L_phi", h2_resVsPt_phi_L_68, h2_resVsPt_phi_L_90, h2_resVsPt_phi_L_99, 0, max_phi);
2524   makeResidualIntervalPlot(
2525       type, DIR, "resVsPt_L_eta", h2_resVsPt_eta_L_68, h2_resVsPt_eta_L_90, h2_resVsPt_eta_L_99, 0, max_eta);
2526   //makeResidualIntervalPlot( type, DIR, "resVsPt_L_d0", h2_resVsPt_d0_L_68, h2_resVsPt_d0_L_90, h2_resVsPt_d0_L_99, 0, 0.02 );
2527 
2528   makeResidualIntervalPlot(
2529       type, DIR, "resVsEta_eta", h2_resVsEta_eta_68, h2_resVsEta_eta_90, h2_resVsEta_eta_99, 0, max_eta);
2530   makeResidualIntervalPlot(
2531       type, DIR, "resVsEta_z0", h2_resVsEta_z0_68, h2_resVsEta_z0_90, h2_resVsEta_z0_99, 0, max_z0);
2532   makeResidualIntervalPlot(
2533       type, DIR, "resVsEta_phi", h2_resVsEta_phi_68, h2_resVsEta_phi_90, h2_resVsEta_phi_99, 0, max_phi);
2534   makeResidualIntervalPlot(
2535       type, DIR, "resVsEta_ptRel", h2_resVsEta_ptRel_68, h2_resVsEta_ptRel_90, h2_resVsEta_ptRel_99, 0, max_eta_ptRel);
2536   //makeResidualIntervalPlot( type, DIR, "resVsEta_d0", h2_resVsEta_d0_68, h2_resVsEta_d0_90, h2_resVsEta_d0_99, 0, 0.02 );
2537 
2538   makeResidualIntervalPlot(
2539       type, DIR, "resVsEta_L_eta", h2_resVsEta_eta_L_68, h2_resVsEta_eta_L_90, h2_resVsEta_eta_L_99, 0, max_eta);
2540   makeResidualIntervalPlot(
2541       type, DIR, "resVsEta_L_z0", h2_resVsEta_z0_L_68, h2_resVsEta_z0_L_90, h2_resVsEta_z0_L_99, 0, max_z0);
2542   makeResidualIntervalPlot(
2543       type, DIR, "resVsEta_L_phi", h2_resVsEta_phi_L_68, h2_resVsEta_phi_L_90, h2_resVsEta_phi_L_99, 0, max_phi);
2544   makeResidualIntervalPlot(type,
2545                            DIR,
2546                            "resVsEta_L_ptRel",
2547                            h2_resVsEta_ptRel_L_68,
2548                            h2_resVsEta_ptRel_L_90,
2549                            h2_resVsEta_ptRel_L_99,
2550                            0,
2551                            max_eta_ptRel);
2552   //makeResidualIntervalPlot( type, DIR, "resVsEta_L_d0", h2_resVsEta_d0_L_68, h2_resVsEta_d0_L_90, h2_resVsEta_d0_L_99, 0, 0.02 );
2553 
2554   makeResidualIntervalPlot(
2555       type, DIR, "resVsEta_H_eta", h2_resVsEta_eta_H_68, h2_resVsEta_eta_H_90, h2_resVsEta_eta_H_99, 0, max_eta);
2556   makeResidualIntervalPlot(
2557       type, DIR, "resVsEta_H_z0", h2_resVsEta_z0_H_68, h2_resVsEta_z0_H_90, h2_resVsEta_z0_H_99, 0, max_z0);
2558   makeResidualIntervalPlot(
2559       type, DIR, "resVsEta_H_phi", h2_resVsEta_phi_H_68, h2_resVsEta_phi_H_90, h2_resVsEta_phi_H_99, 0, max_phi);
2560   makeResidualIntervalPlot(type,
2561                            DIR,
2562                            "resVsEta_H_ptRel",
2563                            h2_resVsEta_ptRel_H_68,
2564                            h2_resVsEta_ptRel_H_90,
2565                            h2_resVsEta_ptRel_H_99,
2566                            0,
2567                            max_eta_ptRel);
2568   //makeResidualIntervalPlot( type, DIR, "resVsEta_H_d0", h2_resVsEta_d0_H_68, h2_resVsEta_d0_H_90, h2_resVsEta_d0_H_99, 0, 0.02 );
2569 
2570   if (doDetailedPlots) {
2571     makeResidualIntervalPlot(
2572         type, DIR, "resVsPhi_ptRel", h2_resVsPhi_ptRel_68, h2_resVsPhi_ptRel_90, h2_resVsPhi_ptRel_99, 0, 0.5);
2573     makeResidualIntervalPlot(type, DIR, "resVsPhi_pt", h2_resVsPhi_pt_68, h2_resVsPhi_pt_90, h2_resVsPhi_pt_99, 0, 20);
2574   }
2575 
2576   // ----------------------------------------------------------------------------------------------------------
2577   // resoultion vs pt
2578   // ----------------------------------------------------------------------------------------------------------
2579 
2580   h2_resVsPt_pt_90->SetMinimum(0);
2581   h2_resVsPt_pt_90->SetMarkerStyle(20);
2582   h2_resVsPt_pt_90->Draw("p");
2583   c.SaveAs(DIR + type + "_resVsPt_pt_90.pdf");
2584 
2585   h2_resVsPt_ptRel_90->SetMinimum(0);
2586   h2_resVsPt_ptRel_90->SetMarkerStyle(20);
2587   h2_resVsPt_ptRel_90->Draw("p");
2588   c.SaveAs(DIR + type + "_resVsPt_ptRel_90.pdf");
2589 
2590   h2_resVsPt_z0_90->SetMinimum(0);
2591   h2_resVsPt_z0_90->SetMarkerStyle(20);
2592   h2_resVsPt_z0_90->Draw("p");
2593   c.SaveAs(DIR + type + "_resVsPt_z0_90.pdf");
2594 
2595   h2_resVsPt_phi_90->SetMinimum(0);
2596   h2_resVsPt_phi_90->SetMarkerStyle(20);
2597   h2_resVsPt_phi_90->Draw("p");
2598   c.SaveAs(DIR + type + "_resVsPt_phi_90.pdf");
2599 
2600   h2_resVsPt_eta_90->SetMinimum(0);
2601   h2_resVsPt_eta_90->SetMarkerStyle(20);
2602   h2_resVsPt_eta_90->Draw("p");
2603   c.SaveAs(DIR + type + "_resVsPt_eta_90.pdf");
2604 
2605   /*
2606   h2_resVsPt_phi_90->SetMinimum(0);
2607   h2_resVsPt_d0_90->SetMarkerStyle(20);
2608   h2_resVsPt_d0_90->Draw("p");
2609   c.SaveAs(DIR+type+"_resVsPt_d0_90.pdf");
2610   */
2611 
2612   // Limit decimal places for doubles in normalized tracklet eff. graphs
2613   gStyle->SetPaintTextFormat("4.2f");
2614 
2615   h_trk_tracklet_eff->SetMinimum(0);
2616   h_trk_tracklet_eff->Draw("col text");
2617   c.SaveAs(DIR + type + "_trk_tracklet_eff.pdf");
2618 
2619   // Remove decimal places for ints in in tracklet hits graphs
2620   gStyle->SetPaintTextFormat("4.0f");
2621 
2622   h_trk_tracklet_hits->SetMinimum(0);
2623   h_trk_tracklet_hits->Draw("col text");
2624   c.SaveAs(DIR + type + "_trk_tracklet_hits.pdf");
2625 
2626   if (doDetailedPlots) {
2627     h2_resVsPt_eta->Write();
2628 
2629     h2_resVsPt_pt->Write();
2630     h2_resVsPt_pt_C->Write();
2631     h2_resVsPt_pt_I->Write();
2632     h2_resVsPt_pt_F->Write();
2633 
2634     h2_resVsPt_ptRel->Write();
2635     h2_resVsPt_ptRel_C->Write();
2636     h2_resVsPt_ptRel_I->Write();
2637     h2_resVsPt_ptRel_F->Write();
2638 
2639     h2_mresVsPt_pt->Write();
2640     h2_mresVsPt_pt_C->Write();
2641     h2_mresVsPt_pt_I->Write();
2642     h2_mresVsPt_pt_F->Write();
2643 
2644     h2_resVsPt_d0->Write();
2645 
2646     h2_resVsPt_z0_C->Write();
2647     h2_resVsPt_z0_I->Write();
2648     h2_resVsPt_z0_F->Write();
2649 
2650     h2_resVsPt_phi->Write();
2651     h2_resVsPt_phi_C->Write();
2652     h2_resVsPt_phi_I->Write();
2653     h2_resVsPt_phi_F->Write();
2654   }
2655 
2656   // ----------------------------------------------------------------------------------------------------------
2657   // resolution vs eta
2658   // ----------------------------------------------------------------------------------------------------------
2659 
2660   h2_resVsEta_eta_90->SetMinimum(0);
2661   h2_resVsEta_eta_90->SetMarkerStyle(20);
2662   h2_resVsEta_eta_90->Draw("p");
2663   c.SaveAs(DIR + type + "_resVsEta_eta_90.pdf");
2664 
2665   h2_resVsEta_eta_68->SetMinimum(0);
2666   h2_resVsEta_eta_68->SetMarkerStyle(20);
2667   h2_resVsEta_eta_68->Draw("p");
2668   c.SaveAs(DIR + type + "_resVsEta_eta_68.pdf");
2669 
2670   if (doDetailedPlots) {
2671     h2_resVsEta_eta_L_90->Draw("p");
2672     sprintf(ctxt, "p_{T} < 8 GeV");
2673     mySmallText(0.22, 0.82, 1, ctxt);
2674     c.SaveAs(DIR + type + "_resVsEta_eta_L_90.pdf");
2675 
2676     h2_resVsEta_eta_H_90->Draw("p");
2677     sprintf(ctxt, "p_{T} > 8 GeV");
2678     mySmallText(0.22, 0.82, 1, ctxt);
2679     c.SaveAs(DIR + type + "_resVsEta_eta_H_90.pdf");
2680   }
2681 
2682   h2_resVsEta_z0_90->SetMinimum(0);
2683   h2_resVsEta_z0_90->SetMarkerStyle(20);
2684   h2_resVsEta_z0_90->Draw("p");
2685   c.SaveAs(DIR + type + "_resVsEta_z0_90.pdf");
2686 
2687   h2_resVsEta_z0_68->SetMinimum(0);
2688   h2_resVsEta_z0_68->SetMarkerStyle(20);
2689   h2_resVsEta_z0_68->Draw("p");
2690   c.SaveAs(DIR + type + "_resVsEta_z0_68.pdf");
2691 
2692   if (doDetailedPlots) {
2693     h2_resVsEta_z0_L_90->Draw();
2694     sprintf(ctxt, "p_{T} < 8 GeV");
2695     mySmallText(0.22, 0.82, 1, ctxt);
2696     c.SaveAs(DIR + type + "_resVsEta_z0_L_90.pdf");
2697 
2698     h2_resVsEta_z0_H_90->Draw();
2699     sprintf(ctxt, "p_{T} > 8 GeV");
2700     mySmallText(0.22, 0.82, 1, ctxt);
2701     c.SaveAs(DIR + type + "_resVsEta_z0_H_90.pdf");
2702   }
2703 
2704   /*
2705   h2_resVsEta_d0_90->Draw();
2706   c.SaveAs(DIR+type+"_resVsEta_d0_90.pdf");
2707 
2708   h2_resVsEta_d0_L_90->Draw();
2709   sprintf(ctxt,"p_{T} < 8 GeV");
2710   mySmallText(0.22,0.82,1,ctxt);
2711   c.SaveAs(DIR+type+"_resVsEta_d0_L_90.pdf");
2712 
2713   h2_resVsEta_d0_H_90->Draw();
2714   sprintf(ctxt,"p_{T} > 8 GeV");
2715   mySmallText(0.22,0.82,1,ctxt);
2716   c.SaveAs(DIR+type+"_resVsEta_d0_H_90.pdf");
2717   */
2718 
2719   h2_resVsEta_phi_90->SetMinimum(0);
2720   h2_resVsEta_phi_90->SetMarkerStyle(20);
2721   h2_resVsEta_phi_90->Draw("p");
2722   c.SaveAs(DIR + type + "_resVsEta_phi_90.pdf");
2723 
2724   h2_resVsEta_phi_68->SetMinimum(0);
2725   h2_resVsEta_phi_68->SetMarkerStyle(20);
2726   h2_resVsEta_phi_68->Draw("p");
2727   c.SaveAs(DIR + type + "_resVsEta_phi_68.pdf");
2728 
2729   if (doDetailedPlots) {
2730     h2_resVsEta_phi_L_90->Draw();
2731     sprintf(ctxt, "p_{T} < 8 GeV");
2732     mySmallText(0.22, 0.82, 1, ctxt);
2733     c.SaveAs(DIR + type + "_resVsEta_phi_L_90.pdf");
2734 
2735     h2_resVsEta_phi_H_90->Draw();
2736     sprintf(ctxt, "p_{T} > 8 GeV");
2737     mySmallText(0.22, 0.82, 1, ctxt);
2738     c.SaveAs(DIR + type + "_resVsEta_phi_H_90.pdf");
2739   }
2740 
2741   h2_resVsEta_ptRel_90->SetMinimum(0);
2742   h2_resVsEta_ptRel_90->SetMarkerStyle(20);
2743   h2_resVsEta_ptRel_90->Draw("p");
2744   c.SaveAs(DIR + type + "_resVsEta_ptRel_90.pdf");
2745 
2746   h2_resVsEta_ptRel_68->SetMinimum(0);
2747   h2_resVsEta_ptRel_68->SetMarkerStyle(20);
2748   h2_resVsEta_ptRel_68->Draw("p");
2749   c.SaveAs(DIR + type + "_resVsEta_ptRel_68.pdf");
2750 
2751   if (doDetailedPlots) {
2752     h2_resVsEta_ptRel_L_90->Draw();
2753     sprintf(ctxt, "p_{T} < 8 GeV");
2754     mySmallText(0.22, 0.82, 1, ctxt);
2755     c.SaveAs(DIR + type + "_resVsEta_ptRel_L_90.pdf");
2756 
2757     h2_resVsEta_ptRel_H_90->Draw();
2758     sprintf(ctxt, "p_{T} > 8 GeV");
2759     mySmallText(0.22, 0.82, 1, ctxt);
2760     c.SaveAs(DIR + type + "_resVsEta_ptRel_H_90.pdf");
2761 
2762     h2_resVsEta_eta->Write();
2763     h2_resVsEta_eta_L->Write();
2764     h2_resVsEta_eta_H->Write();
2765 
2766     h2_resVsEta_z0->Draw();
2767     c.SaveAs(DIR + type + "_resVsEta_z0_rms.pdf");
2768     h2_resVsEta_eta->Draw();
2769     c.SaveAs(DIR + type + "_resVsEta_eta_rms.pdf");
2770     h2_resVsEta_ptRel->Draw();
2771     c.SaveAs(DIR + type + "_resVsEta_ptRel_rms.pdf");
2772     h2_resVsEta_phi->Draw();
2773     c.SaveAs(DIR + type + "_resVsEta_phi_rms.pdf");
2774 
2775     // check residuals
2776     h2_mresVsEta_z0->Draw();
2777     h2_mresVsEta_z0->Write();
2778     c.SaveAs(DIR + type + "_mean-resVsEta_z0.pdf");
2779 
2780     h2_mresVsEta_eta->Draw();
2781     h2_mresVsEta_eta->Write();
2782     c.SaveAs(DIR + type + "_mean-resVsEta_eta.pdf");
2783 
2784     h2_mresVsEta_ptRel->Draw();
2785     h2_mresVsEta_ptRel->Write();
2786     c.SaveAs(DIR + type + "_mean-resVsEta_ptRel.pdf");
2787 
2788     h2_mresVsEta_phi->Draw();
2789     h2_mresVsEta_phi->Write();
2790     c.SaveAs(DIR + type + "_mean-resVsEta_phi.pdf");
2791 
2792     h2_resVsEta_z0->Write();
2793     h2_resVsEta_z0_L->Write();
2794     h2_resVsEta_z0_H->Write();
2795 
2796     h2_resVsEta_d0->Write();
2797     h2_resVsEta_d0_L->Write();
2798     h2_resVsEta_d0_H->Write();
2799 
2800     h2_resVsEta_phi->Write();
2801     h2_resVsEta_phi_L->Write();
2802     h2_resVsEta_phi_H->Write();
2803 
2804     h2_resVsEta_ptRel->Write();
2805     h2_resVsEta_ptRel_L->Write();
2806     h2_resVsEta_ptRel_H->Write();
2807   }
2808 
2809   if (doGausFit) {
2810     h3_resVsEta_eta_L->Write();
2811     h3_resVsEta_z0_L->Write();
2812     h3_resVsEta_phi_L->Write();
2813     h3_resVsEta_ptRel_L->Write();
2814 
2815     h3_resVsEta_eta_H->Write();
2816     h3_resVsEta_z0_H->Write();
2817     h3_resVsEta_phi_H->Write();
2818     h3_resVsEta_ptRel_H->Write();
2819   }
2820 
2821   // resolution vs phi
2822   if (doDetailedPlots) {
2823     h2_resVsPhi_pt_90->SetMinimum(0);
2824     h2_resVsPhi_pt_90->SetMarkerStyle(20);
2825     h2_resVsPhi_pt_90->Draw("p");
2826     c.SaveAs(DIR + type + "_resVsPhi_pt_90.pdf");
2827 
2828     h2_resVsPhi_ptRel_90->SetMinimum(0);
2829     h2_resVsPhi_ptRel_90->SetMarkerStyle(20);
2830     h2_resVsPhi_ptRel_90->Draw("p");
2831     c.SaveAs(DIR + type + "_resVsPhi_ptRel_90.pdf");
2832   }
2833 
2834   // ----------------------------------------------------------------------------------------------------------------
2835   // track quality plots
2836   // ----------------------------------------------------------------------------------------------------------------
2837 
2838   if (doDetailedPlots) {
2839     h_match_trk_nstub->Write();
2840     h_match_trk_nstub_C->Write();
2841     h_match_trk_nstub_I->Write();
2842     h_match_trk_nstub_F->Write();
2843   }
2844 
2845   h_trk_chi2->Draw();
2846   sprintf(ctxt, "|eta| < 2.4");
2847   mySmallText(0.52, 0.82, 1, ctxt);
2848   c.SaveAs(DIR + type + "_trk_chi2.pdf");
2849 
2850   h_trk_chi2_dof->Draw();
2851   sprintf(ctxt, "|eta| < 2.4");
2852   mySmallText(0.52, 0.82, 1, ctxt);
2853   c.SaveAs(DIR + type + "_trk_chi2_dof.pdf");
2854 
2855   h_trk_chi2rphi->Draw();
2856   sprintf(ctxt, "|eta| < 2.4");
2857   mySmallText(0.52, 0.82, 1, ctxt);
2858   c.SaveAs(DIR + type + "_trk_chi2rphi.pdf");
2859 
2860   h_trk_chi2rphi_dof->Draw();
2861   sprintf(ctxt, "|eta| < 2.4");
2862   mySmallText(0.52, 0.82, 1, ctxt);
2863   c.SaveAs(DIR + type + "_trk_chi2rphi_dof.pdf");
2864 
2865   h_trk_chi2rz->Draw();
2866   sprintf(ctxt, "|eta| < 2.4");
2867   mySmallText(0.52, 0.82, 1, ctxt);
2868   c.SaveAs(DIR + type + "_trk_chi2rz.pdf");
2869 
2870   h_trk_chi2rz_dof->Draw();
2871   sprintf(ctxt, "|eta| < 2.4");
2872   mySmallText(0.52, 0.82, 1, ctxt);
2873   c.SaveAs(DIR + type + "_trk_chi2rz_dof.pdf");
2874 
2875   h_trk_chi2->Write();
2876   h_trk_chi2rphi->Write();
2877   h_trk_chi2rz->Write();
2878   h_match_trk_chi2->Write();
2879   h_match_trk_chi2rphi->Write();
2880   h_match_trk_chi2rz->Write();
2881 
2882   if (doDetailedPlots) {
2883     h_match_trk_chi2_C_L->Write();
2884     h_match_trk_chi2_I_L->Write();
2885     h_match_trk_chi2_F_L->Write();
2886     h_match_trk_chi2_C_H->Write();
2887     h_match_trk_chi2_I_H->Write();
2888     h_match_trk_chi2_F_H->Write();
2889 
2890     h_match_trk_chi2_dof_C_L->Write();
2891     h_match_trk_chi2_dof_I_L->Write();
2892     h_match_trk_chi2_dof_F_L->Write();
2893     h_match_trk_chi2_dof_C_H->Write();
2894     h_match_trk_chi2_dof_I_H->Write();
2895     h_match_trk_chi2_dof_F_H->Write();
2896   }
2897 
2898   // ----------------------------------------------------------------------------------------------------------------
2899   // efficiency plots
2900   // ----------------------------------------------------------------------------------------------------------------
2901 
2902   // rebin pt/phi plots
2903   h_tp_pt->Rebin(2);
2904   h_match_tp_pt->Rebin(2);
2905   h_tp_phi->Rebin(2);
2906   h_match_tp_phi->Rebin(2);
2907 
2908   h_tp_pt_L->Rebin(2);
2909   h_match_tp_pt_L->Rebin(2);
2910   h_tp_pt_LC->Rebin(2);
2911   h_match_tp_pt_LC->Rebin(2);
2912   h_tp_pt_H->Rebin(2);
2913   h_match_tp_pt_H->Rebin(2);
2914 
2915   // h_tp_eta->Rebin(2);
2916   // h_match_tp_eta->Rebin(2);
2917   // h_tp_eta_L->Rebin(2);
2918   // h_match_tp_eta_L->Rebin(2);
2919   // h_tp_eta_H->Rebin(2);
2920   // h_match_tp_eta_H->Rebin(2);
2921 
2922   // calculate the efficiency
2923   h_match_tp_pt->Sumw2();
2924   h_tp_pt->Sumw2();
2925   TH1F* h_eff_pt = (TH1F*)h_match_tp_pt->Clone();
2926   h_eff_pt->SetName("eff_pt");
2927   h_eff_pt->GetYaxis()->SetTitle("Efficiency");
2928   h_eff_pt->Divide(h_match_tp_pt, h_tp_pt, 1.0, 1.0, "B");
2929 
2930   h_match_tp_pt_L->Sumw2();
2931   h_tp_pt_L->Sumw2();
2932   TH1F* h_eff_pt_L = (TH1F*)h_match_tp_pt_L->Clone();
2933   h_eff_pt_L->SetName("eff_pt_L");
2934   h_eff_pt_L->GetYaxis()->SetTitle("Efficiency");
2935   h_eff_pt_L->Divide(h_match_tp_pt_L, h_tp_pt_L, 1.0, 1.0, "B");
2936 
2937   h_match_tp_pt_LC->Sumw2();
2938   h_tp_pt_LC->Sumw2();
2939   TH1F* h_eff_pt_LC = (TH1F*)h_match_tp_pt_LC->Clone();
2940   h_eff_pt_LC->SetName("eff_pt_LC");
2941   h_eff_pt_LC->GetYaxis()->SetTitle("Efficiency");
2942   h_eff_pt_LC->Divide(h_match_tp_pt_LC, h_tp_pt_LC, 1.0, 1.0, "B");
2943 
2944   h_match_tp_pt_H->Sumw2();
2945   h_tp_pt_H->Sumw2();
2946   TH1F* h_eff_pt_H = (TH1F*)h_match_tp_pt_H->Clone();
2947   h_eff_pt_H->SetName("eff_pt_H");
2948   h_eff_pt_H->GetYaxis()->SetTitle("Efficiency");
2949   h_eff_pt_H->Divide(h_match_tp_pt_H, h_tp_pt_H, 1.0, 1.0, "B");
2950 
2951   h_match_tp_eta->Sumw2();
2952   h_tp_eta->Sumw2();
2953   TH1F* h_eff_eta = (TH1F*)h_match_tp_eta->Clone();
2954   h_eff_eta->SetName("eff_eta");
2955   h_eff_eta->GetYaxis()->SetTitle("Efficiency");
2956   h_eff_eta->Divide(h_match_tp_eta, h_tp_eta, 1.0, 1.0, "B");
2957 
2958   h_match_tp_eta_L->Sumw2();
2959   h_tp_eta_L->Sumw2();
2960   TH1F* h_eff_eta_L = (TH1F*)h_match_tp_eta_L->Clone();
2961   h_eff_eta_L->SetName("eff_eta_L");
2962   h_eff_eta_L->GetYaxis()->SetTitle("Efficiency");
2963   h_eff_eta_L->Divide(h_match_tp_eta_L, h_tp_eta_L, 1.0, 1.0, "B");
2964 
2965   h_match_tp_eta_H->Sumw2();
2966   h_tp_eta_H->Sumw2();
2967   TH1F* h_eff_eta_H = (TH1F*)h_match_tp_eta_H->Clone();
2968   h_eff_eta_H->SetName("eff_eta_H");
2969   h_eff_eta_H->GetYaxis()->SetTitle("Efficiency");
2970   h_eff_eta_H->Divide(h_match_tp_eta_H, h_tp_eta_H, 1.0, 1.0, "B");
2971 
2972   h_match_tp_eta_23->Sumw2();
2973   h_tp_eta_23->Sumw2();
2974   TH1F* h_eff_eta_23 = (TH1F*)h_match_tp_eta_23->Clone();
2975   h_eff_eta_23->SetName("eff_eta_23");
2976   h_eff_eta_23->GetYaxis()->SetTitle("Efficiency");
2977   h_eff_eta_23->Divide(h_match_tp_eta_23, h_tp_eta_23, 1.0, 1.0, "B");
2978 
2979   h_match_tp_eta_35->Sumw2();
2980   h_tp_eta_35->Sumw2();
2981   TH1F* h_eff_eta_35 = (TH1F*)h_match_tp_eta_35->Clone();
2982   h_eff_eta_35->SetName("eff_eta_35");
2983   h_eff_eta_35->GetYaxis()->SetTitle("Efficiency");
2984   h_eff_eta_35->Divide(h_match_tp_eta_35, h_tp_eta_35, 1.0, 1.0, "B");
2985 
2986   h_match_tp_eta_5->Sumw2();
2987   h_tp_eta_5->Sumw2();
2988   TH1F* h_eff_eta_5 = (TH1F*)h_match_tp_eta_5->Clone();
2989   h_eff_eta_5->SetName("eff_eta_5");
2990   h_eff_eta_5->GetYaxis()->SetTitle("Efficiency");
2991   h_eff_eta_5->Divide(h_match_tp_eta_5, h_tp_eta_5, 1.0, 1.0, "B");
2992 
2993   h_match_tp_phi->Sumw2();
2994   h_tp_phi->Sumw2();
2995   TH1F* h_eff_phi = (TH1F*)h_match_tp_phi->Clone();
2996   h_eff_phi->SetName("eff_phi");
2997   h_eff_phi->GetYaxis()->SetTitle("Efficiency");
2998   h_eff_phi->Divide(h_match_tp_phi, h_tp_phi, 1.0, 1.0, "B");
2999 
3000   h_match_tp_z0->Sumw2();
3001   h_tp_z0->Sumw2();
3002   TH1F* h_eff_z0 = (TH1F*)h_match_tp_z0->Clone();
3003   h_eff_z0->SetName("eff_z0");
3004   h_eff_z0->GetYaxis()->SetTitle("Efficiency");
3005   h_eff_z0->Divide(h_match_tp_z0, h_tp_z0, 1.0, 1.0, "B");
3006 
3007   h_match_tp_z0_L->Sumw2();
3008   h_tp_z0_L->Sumw2();
3009   TH1F* h_eff_z0_L = (TH1F*)h_match_tp_z0_L->Clone();
3010   h_eff_z0_L->SetName("eff_z0_L");
3011   h_eff_z0_L->GetYaxis()->SetTitle("Efficiency");
3012   h_eff_z0_L->Divide(h_match_tp_z0_L, h_tp_z0_L, 1.0, 1.0, "B");
3013 
3014   h_match_tp_z0_H->Sumw2();
3015   h_tp_z0_H->Sumw2();
3016   TH1F* h_eff_z0_H = (TH1F*)h_match_tp_z0_H->Clone();
3017   h_eff_z0_H->SetName("eff_z0_H");
3018   h_eff_z0_H->GetYaxis()->SetTitle("Efficiency");
3019   h_eff_z0_H->Divide(h_match_tp_z0_H, h_tp_z0_H, 1.0, 1.0, "B");
3020 
3021   h_match_tp_d0->Sumw2();
3022   h_tp_d0->Sumw2();
3023   TH1F* h_eff_d0 = (TH1F*)h_match_tp_d0->Clone();
3024   h_eff_d0->SetName("eff_d0");
3025   h_eff_d0->GetYaxis()->SetTitle("Efficiency");
3026   h_eff_d0->Divide(h_match_tp_d0, h_tp_d0, 1.0, 1.0, "B");
3027 
3028   h_match_tp_absd0->Sumw2();
3029   h_tp_absd0->Sumw2();
3030   TH1F* h_eff_absd0 = (TH1F*)h_match_tp_absd0->Clone();
3031   h_eff_absd0->SetName("eff_absd0");
3032   h_eff_absd0->GetYaxis()->SetTitle("Efficiency");
3033   h_eff_absd0->Divide(h_match_tp_absd0, h_tp_absd0, 1.0, 1.0, "B");
3034 
3035   h_match_tp_absd0_eta2->Sumw2();
3036   h_tp_absd0_eta2->Sumw2();
3037   TH1F* h_eff_absd0_eta2 = (TH1F*)h_match_tp_absd0_eta2->Clone();
3038   h_eff_absd0_eta2->SetName("eff_absd0_eta2");
3039   h_eff_absd0_eta2->GetYaxis()->SetTitle("Efficiency");
3040   h_eff_absd0_eta2->Divide(h_match_tp_absd0_eta2, h_tp_absd0_eta2, 1.0, 1.0, "B");
3041 
3042   h_match_tp_absd0_eta2_pt3->Sumw2();
3043   h_tp_absd0_eta2_pt3->Sumw2();
3044   TH1F* h_eff_absd0_eta2_pt3 = (TH1F*)h_match_tp_absd0_eta2_pt3->Clone();
3045   h_eff_absd0_eta2_pt3->SetName("eff_absd0_eta2_pt3");
3046   h_eff_absd0_eta2_pt3->GetYaxis()->SetTitle("Efficiency");
3047   h_eff_absd0_eta2_pt3->Divide(h_match_tp_absd0_eta2_pt3, h_tp_absd0_eta2_pt3, 1.0, 1.0, "B");
3048 
3049   // set the axis range
3050   h_eff_pt->SetAxisRange(0, 1.1, "Y");
3051   h_eff_pt_L->SetAxisRange(0, 1.1, "Y");
3052   h_eff_pt_LC->SetAxisRange(0, 1.1, "Y");
3053   h_eff_pt_H->SetAxisRange(0, 1.1, "Y");
3054   h_eff_eta->SetAxisRange(0, 1.1, "Y");
3055   h_eff_eta_L->SetAxisRange(0, 1.1, "Y");
3056   h_eff_eta_H->SetAxisRange(0, 1.1, "Y");
3057   h_eff_eta_23->SetAxisRange(0, 1.1, "Y");
3058   h_eff_eta_35->SetAxisRange(0, 1.1, "Y");
3059   h_eff_eta_5->SetAxisRange(0, 1.1, "Y");
3060   h_eff_phi->SetAxisRange(0, 1.1, "Y");
3061   h_eff_z0->SetAxisRange(0, 1.1, "Y");
3062   h_eff_z0_L->SetAxisRange(0, 1.1, "Y");
3063   h_eff_z0_H->SetAxisRange(0, 1.1, "Y");
3064   h_eff_d0->SetAxisRange(0, 1.1, "Y");
3065   h_eff_absd0->SetAxisRange(0, 1.1, "Y");
3066   h_eff_absd0_eta2->SetAxisRange(0, 1.1, "Y");
3067   h_eff_absd0_eta2_pt3->SetAxisRange(0, 1.1, "Y");
3068 
3069   gPad->SetGridx();
3070   gPad->SetGridy();
3071 
3072   // draw and save plots
3073   h_eff_pt->Draw();
3074   h_eff_pt->Write();
3075   c.SaveAs(DIR + type + "_eff_pt.pdf");
3076 
3077   if (type.Contains("Mu")) {
3078     h_eff_pt->GetYaxis()->SetRangeUser(0.8, 1.01);  // zoomed-in plot
3079     c.SaveAs(DIR + type + "_eff_pt_zoom.pdf");
3080   }
3081 
3082   h_eff_pt_L->Draw();
3083   h_eff_pt_L->Write();
3084   sprintf(ctxt, "p_{T} < 8 GeV");
3085   mySmallText(0.45, 0.5, 1, ctxt);
3086   c.SaveAs(DIR + type + "_eff_pt_L.pdf");
3087 
3088   if (doDetailedPlots) {
3089     h_eff_pt_LC->Draw();
3090     h_eff_pt_LC->Write();
3091     sprintf(ctxt, "p_{T} < 8 GeV, |#eta|<1.0");
3092     mySmallText(0.45, 0.5, 1, ctxt);
3093     c.SaveAs(DIR + type + "_eff_pt_LC.pdf");
3094   }
3095   h_eff_pt_H->Draw();
3096   h_eff_pt_H->Write();
3097   sprintf(ctxt, "p_{T} > 8 GeV");
3098   mySmallText(0.45, 0.5, 1, ctxt);
3099   c.SaveAs(DIR + type + "_eff_pt_H.pdf");
3100 
3101   h_eff_eta->Draw();
3102   h_eff_eta->Write();
3103   c.SaveAs(DIR + type + "_eff_eta.pdf");
3104 
3105   if (type.Contains("Mu")) {
3106     h_eff_eta->GetYaxis()->SetRangeUser(0.8, 1.01);  // zoomed-in plot
3107     c.SaveAs(DIR + type + "_eff_eta_zoom.pdf");
3108   }
3109 
3110   h_eff_eta_L->Draw();
3111   h_eff_eta_L->Write();
3112   sprintf(ctxt, "p_{T} < 8 GeV");
3113   mySmallText(0.45, 0.5, 1, ctxt);
3114   c.SaveAs(DIR + type + "_eff_eta_L.pdf");
3115 
3116   h_eff_eta_H->Draw();
3117   h_eff_eta_H->Write();
3118   sprintf(ctxt, "p_{T} > 8 GeV");
3119   mySmallText(0.45, 0.5, 1, ctxt);
3120   c.SaveAs(DIR + type + "_eff_eta_H.pdf");
3121 
3122   h_eff_eta_23->Write();
3123   h_eff_eta_35->Write();
3124   h_eff_eta_5->Write();
3125 
3126   if (doDetailedPlots) {
3127     h_eff_eta_23->Draw();
3128     sprintf(ctxt, "2 < p_{T} < 3 GeV");
3129     mySmallText(0.45, 0.5, 1, ctxt);
3130     c.SaveAs(DIR + type + "_eff_eta_23.pdf");
3131 
3132     h_eff_eta_35->Draw();
3133     sprintf(ctxt, "3 < p_{T} < 5 GeV");
3134     mySmallText(0.45, 0.5, 1, ctxt);
3135     c.SaveAs(DIR + type + "_eff_eta_35.pdf");
3136 
3137     h_eff_eta_5->Draw();
3138     sprintf(ctxt, "p_{T} > 5 GeV");
3139     mySmallText(0.45, 0.5, 1, ctxt);
3140     c.SaveAs(DIR + type + "_eff_eta_5.pdf");
3141 
3142     h_eff_z0->Draw();
3143     h_eff_z0->Write();
3144     c.SaveAs(DIR + type + "_eff_z0.pdf");
3145 
3146     h_eff_z0_L->Write();
3147     h_eff_z0_H->Write();
3148 
3149     h_eff_phi->Draw();
3150     h_eff_phi->Write();
3151     c.SaveAs(DIR + type + "_eff_phi.pdf");
3152 
3153     if (type.Contains("Mu")) {
3154       h_eff_phi->GetYaxis()->SetRangeUser(0.8, 1.01);  // zoomed-in plot
3155       c.SaveAs(DIR + type + "_eff_phi_zoom.pdf");
3156     }
3157   }
3158 
3159   if (doDetailedPlots || TP_maxD0 > 1.0) {
3160     h_eff_d0->Write();
3161     h_eff_absd0->Write();
3162     h_eff_absd0_eta2->Write();
3163     h_eff_absd0_eta2_pt3->Write();
3164 
3165     h_eff_d0->Draw();
3166     c.SaveAs(DIR + type + "_eff_d0.pdf");
3167 
3168     h_eff_absd0->Draw();
3169     c.SaveAs(DIR + type + "_eff_absd0.pdf");
3170 
3171     h_eff_absd0_eta2->Draw();
3172     c.SaveAs(DIR + type + "_eff_absd0_eta2.pdf");
3173 
3174     h_eff_absd0_eta2_pt3->Draw();
3175     c.SaveAs(DIR + type + "_eff_absd0_eta2_pt3.pdf");
3176   }
3177 
3178   gPad->SetGridx(0);
3179   gPad->SetGridy(0);
3180 
3181   // ----------------------------------------------------------------------------------------------------------------
3182   // more resolution plots
3183   // ----------------------------------------------------------------------------------------------------------------
3184 
3185   float rms = 0;
3186 
3187   if (doDetailedPlots) {
3188     // draw and save plots
3189     h_res_pt->Draw();
3190     rms = h_res_pt->GetRMS();
3191     sprintf(ctxt, "RMS = %.4f", rms);
3192     mySmallText(0.22, 0.82, 1, ctxt);
3193     c.SaveAs(DIR + type + "_res_pt.pdf");
3194 
3195     h_res_ptRel->Draw();
3196     rms = h_res_ptRel->GetRMS();
3197     sprintf(ctxt, "RMS = %.4f", rms);
3198     mySmallText(0.22, 0.82, 1, ctxt);
3199     c.SaveAs(DIR + type + "_res_ptRel.pdf");
3200 
3201     h_res_eta->Draw();
3202     rms = h_res_eta->GetRMS();
3203     sprintf(ctxt, "RMS = %.3e", rms);
3204     mySmallText(0.22, 0.82, 1, ctxt);
3205     c.SaveAs(DIR + type + "_res_eta.pdf");
3206 
3207     h_res_phi->Draw();
3208     rms = h_res_phi->GetRMS();
3209     sprintf(ctxt, "RMS = %.3e", rms);
3210     mySmallText(0.22, 0.82, 1, ctxt);
3211     c.SaveAs(DIR + type + "_res_phi.pdf");
3212 
3213     h_res_z0->Draw();
3214     rms = h_res_z0->GetRMS();
3215     sprintf(ctxt, "RMS = %.4f", rms);
3216     mySmallText(0.22, 0.82, 1, ctxt);
3217     c.SaveAs(DIR + type + "_res_z0.pdf");
3218 
3219     h_res_z0_C->Draw();
3220     rms = h_res_z0_C->GetRMS();
3221     sprintf(ctxt, "RMS = %.4f;", rms);
3222     mySmallText(0.22, 0.82, 1, ctxt);
3223     sprintf(ctxt, "|eta| < 0.8");
3224     mySmallText(0.22, 0.76, 1, ctxt);
3225     c.SaveAs(DIR + type + "_res_z0_C.pdf");
3226 
3227     h_res_z0_I->Draw();
3228     rms = h_res_z0_I->GetRMS();
3229     sprintf(ctxt, "RMS = %.4f;", rms);
3230     mySmallText(0.22, 0.82, 1, ctxt);
3231     sprintf(ctxt, "0.8 < |eta| < 1.6");
3232     mySmallText(0.22, 0.76, 1, ctxt);
3233     c.SaveAs(DIR + type + "_res_z0_I.pdf");
3234 
3235     h_res_z0_F->Draw();
3236     rms = h_res_z0_F->GetRMS();
3237     sprintf(ctxt, "RMS = %.4f;", rms);
3238     mySmallText(0.22, 0.82, 1, ctxt);
3239     sprintf(ctxt, "|eta| > 1.6");
3240     mySmallText(0.22, 0.76, 1, ctxt);
3241     c.SaveAs(DIR + type + "_res_z0_F.pdf");
3242 
3243     h_res_z0_C_L->Draw();
3244     h_res_z0_C_L->Write();
3245     rms = h_res_z0_C_L->GetRMS();
3246     sprintf(ctxt, "RMS = %.4f;", rms);
3247     mySmallText(0.22, 0.82, 1, ctxt);
3248     sprintf(ctxt, "|eta| < 0.8");
3249     mySmallText(0.22, 0.76, 1, ctxt);
3250     c.SaveAs(DIR + type + "_res_z0_C_L.pdf");
3251 
3252     h_res_z0_I_L->Draw();
3253     h_res_z0_I_L->Write();
3254     rms = h_res_z0_I_L->GetRMS();
3255     sprintf(ctxt, "RMS = %.4f;", rms);
3256     mySmallText(0.22, 0.82, 1, ctxt);
3257     sprintf(ctxt, "0.8 < |eta| < 1.6");
3258     mySmallText(0.22, 0.76, 1, ctxt);
3259     c.SaveAs(DIR + type + "_res_z0_I_L.pdf");
3260 
3261     h_res_z0_F_L->Draw();
3262     h_res_z0_F_L->Write();
3263     rms = h_res_z0_F_L->GetRMS();
3264     sprintf(ctxt, "RMS = %.4f;", rms);
3265     mySmallText(0.22, 0.82, 1, ctxt);
3266     sprintf(ctxt, "|eta| > 1.6");
3267     mySmallText(0.22, 0.76, 1, ctxt);
3268     c.SaveAs(DIR + type + "_res_z0_F_L.pdf");
3269 
3270     h_res_z0_C_H->Draw();
3271     h_res_z0_C_H->Write();
3272     rms = h_res_z0_C_H->GetRMS();
3273     sprintf(ctxt, "RMS = %.4f;", rms);
3274     mySmallText(0.22, 0.82, 1, ctxt);
3275     sprintf(ctxt, "|eta| < 0.8");
3276     mySmallText(0.22, 0.76, 1, ctxt);
3277     c.SaveAs(DIR + type + "_res_z0_C_H.pdf");
3278 
3279     h_res_z0_I_H->Draw();
3280     h_res_z0_I_H->Write();
3281     rms = h_res_z0_I_H->GetRMS();
3282     sprintf(ctxt, "RMS = %.4f;", rms);
3283     mySmallText(0.22, 0.82, 1, ctxt);
3284     sprintf(ctxt, "0.8 < |eta| < 1.6");
3285     mySmallText(0.22, 0.76, 1, ctxt);
3286     c.SaveAs(DIR + type + "_res_z0_I_H.pdf");
3287 
3288     h_res_z0_F_H->Draw();
3289     h_res_z0_F_H->Write();
3290     rms = h_res_z0_F_H->GetRMS();
3291     sprintf(ctxt, "RMS = %.4f;", rms);
3292     mySmallText(0.22, 0.82, 1, ctxt);
3293     sprintf(ctxt, "|eta| > 1.6");
3294     mySmallText(0.22, 0.76, 1, ctxt);
3295     c.SaveAs(DIR + type + "_res_z0_F_H.pdf");
3296 
3297     h_res_z0_L->Draw();
3298     h_res_z0_L->Write();
3299     rms = h_res_z0_L->GetRMS();
3300     sprintf(ctxt, "RMS = %.4f;", rms);
3301     mySmallText(0.22, 0.82, 1, ctxt);
3302     sprintf(ctxt, "p_{T} < 5 GeV");
3303     mySmallText(0.22, 0.76, 1, ctxt);
3304     c.SaveAs(DIR + type + "_res_z0_L.pdf");
3305 
3306     h_res_z0_H->Draw();
3307     h_res_z0_H->Write();
3308     rms = h_res_z0_H->GetRMS();
3309     sprintf(ctxt, "RMS = %.4f;", rms);
3310     mySmallText(0.22, 0.82, 1, ctxt);
3311     sprintf(ctxt, "p_{T} > 15 GeV");
3312     mySmallText(0.22, 0.76, 1, ctxt);
3313     c.SaveAs(DIR + type + "_res_z0_H.pdf");
3314 
3315     if (h_res_d0->GetEntries() > 0) {
3316       h_res_d0->Draw();
3317       rms = h_res_d0->GetRMS();
3318       sprintf(ctxt, "RMS = %.4f", rms);
3319       mySmallText(0.22, 0.82, 1, ctxt);
3320       c.SaveAs(DIR + type + "_res_d0.pdf");
3321 
3322       h_res_d0_C->Draw();
3323       h_res_d0_C->Write();
3324       rms = h_res_d0_C->GetRMS();
3325       sprintf(ctxt, "RMS = %.4f;", rms);
3326       mySmallText(0.22, 0.82, 1, ctxt);
3327       sprintf(ctxt, "|eta| < 0.8");
3328       mySmallText(0.22, 0.76, 1, ctxt);
3329       c.SaveAs(DIR + type + "_res_d0_C.pdf");
3330 
3331       h_res_d0_I->Draw();
3332       h_res_d0_I->Write();
3333       rms = h_res_d0_I->GetRMS();
3334       sprintf(ctxt, "RMS = %.4f;", rms);
3335       mySmallText(0.22, 0.82, 1, ctxt);
3336       sprintf(ctxt, "0.8 < |eta| < 1.6");
3337       mySmallText(0.22, 0.76, 1, ctxt);
3338       c.SaveAs(DIR + type + "_res_d0_I.pdf");
3339 
3340       h_res_d0_F->Draw();
3341       h_res_d0_F->Write();
3342       rms = h_res_d0_F->GetRMS();
3343       sprintf(ctxt, "RMS = %.4f;", rms);
3344       mySmallText(0.22, 0.82, 1, ctxt);
3345       sprintf(ctxt, "|eta| > 1.6");
3346       mySmallText(0.22, 0.76, 1, ctxt);
3347       c.SaveAs(DIR + type + "_res_d0_F.pdf");
3348 
3349       h_res_d0_C_L->Draw();
3350       h_res_d0_C_L->Write();
3351       rms = h_res_d0_C_L->GetRMS();
3352       sprintf(ctxt, "RMS = %.4f;", rms);
3353       mySmallText(0.22, 0.82, 1, ctxt);
3354       sprintf(ctxt, "|eta| < 0.8");
3355       mySmallText(0.22, 0.76, 1, ctxt);
3356       c.SaveAs(DIR + type + "_res_d0_C_L.pdf");
3357 
3358       h_res_d0_I_L->Draw();
3359       h_res_d0_I_L->Write();
3360       rms = h_res_d0_I_L->GetRMS();
3361       sprintf(ctxt, "RMS = %.4f;", rms);
3362       mySmallText(0.22, 0.82, 1, ctxt);
3363       sprintf(ctxt, "0.8 < |eta| < 1.6");
3364       mySmallText(0.22, 0.76, 1, ctxt);
3365       c.SaveAs(DIR + type + "_res_d0_I_L.pdf");
3366 
3367       h_res_d0_F_L->Draw();
3368       h_res_d0_F_L->Write();
3369       rms = h_res_d0_F_L->GetRMS();
3370       sprintf(ctxt, "RMS = %.4f;", rms);
3371       mySmallText(0.22, 0.82, 1, ctxt);
3372       sprintf(ctxt, "|eta| > 1.6");
3373       mySmallText(0.22, 0.76, 1, ctxt);
3374       c.SaveAs(DIR + type + "_res_d0_F_L.pdf");
3375 
3376       h_res_d0_C_H->Draw();
3377       h_res_d0_C_H->Write();
3378       rms = h_res_d0_C_H->GetRMS();
3379       sprintf(ctxt, "RMS = %.4f;", rms);
3380       mySmallText(0.22, 0.82, 1, ctxt);
3381       sprintf(ctxt, "|eta| < 0.8");
3382       mySmallText(0.22, 0.76, 1, ctxt);
3383       c.SaveAs(DIR + type + "_res_d0_C_H.pdf");
3384 
3385       h_res_d0_I_H->Draw();
3386       h_res_d0_I_H->Write();
3387       rms = h_res_d0_I_H->GetRMS();
3388       sprintf(ctxt, "RMS = %.4f;", rms);
3389       mySmallText(0.22, 0.82, 1, ctxt);
3390       sprintf(ctxt, "0.8 < |eta| < 1.6");
3391       mySmallText(0.22, 0.76, 1, ctxt);
3392       c.SaveAs(DIR + type + "_res_d0_I_H.pdf");
3393 
3394       h_res_d0_F_H->Draw();
3395       h_res_d0_F_H->Write();
3396       rms = h_res_d0_F_H->GetRMS();
3397       sprintf(ctxt, "RMS = %.4f;", rms);
3398       mySmallText(0.22, 0.82, 1, ctxt);
3399       sprintf(ctxt, "|eta| > 1.6");
3400       mySmallText(0.22, 0.76, 1, ctxt);
3401       c.SaveAs(DIR + type + "_res_d0_F_H.pdf");
3402 
3403       h_res_d0_L->Draw();
3404       h_res_d0_L->Write();
3405       rms = h_res_d0_L->GetRMS();
3406       sprintf(ctxt, "RMS = %.4f;", rms);
3407       mySmallText(0.22, 0.82, 1, ctxt);
3408       sprintf(ctxt, "p_{T} < 5 GeV");
3409       mySmallText(0.22, 0.76, 1, ctxt);
3410       c.SaveAs(DIR + type + "_res_d0_L.pdf");
3411 
3412       h_res_d0_H->Draw();
3413       h_res_d0_H->Write();
3414       rms = h_res_d0_H->GetRMS();
3415       sprintf(ctxt, "RMS = %.4f;", rms);
3416       mySmallText(0.22, 0.82, 1, ctxt);
3417       sprintf(ctxt, "p_{T} > 15 GeV");
3418       mySmallText(0.22, 0.76, 1, ctxt);
3419       c.SaveAs(DIR + type + "_res_d0_H.pdf");
3420     }
3421   }
3422 
3423   // ---------------------------------------------------------------------------------------------------------
3424   // "fake rates"
3425 
3426   h_trk_all_vspt->Sumw2();
3427   h_trk_loose_vspt->Sumw2();
3428   h_trk_genuine_vspt->Sumw2();
3429   h_trk_notloose_vspt->Sumw2();
3430   h_trk_notgenuine_vspt->Sumw2();
3431   h_trk_duplicate_vspt->Sumw2();
3432   h_tp_vspt->Sumw2();
3433 
3434   // fraction of not genuine tracks
3435   TH1F* h_notgenuine_pt = (TH1F*)h_trk_notgenuine_vspt->Clone();
3436   h_notgenuine_pt->SetName("notgenuine_pt");
3437   h_notgenuine_pt->GetYaxis()->SetTitle("Not genuine fraction");
3438   h_notgenuine_pt->Divide(h_trk_notgenuine_vspt, h_trk_all_vspt, 1.0, 1.0, "B");
3439 
3440   h_notgenuine_pt->Write();
3441   h_notgenuine_pt->Draw();
3442   c.SaveAs(DIR + type + "_notgenuine.pdf");
3443 
3444   // fraction of not loosely genuine tracks
3445   TH1F* h_notloose_pt = (TH1F*)h_trk_notloose_vspt->Clone();
3446   h_notloose_pt->SetName("notloose_pt");
3447   h_notloose_pt->GetYaxis()->SetTitle("Not loose fraction");
3448   h_notloose_pt->Divide(h_trk_notloose_vspt, h_trk_all_vspt, 1.0, 1.0, "B");
3449 
3450   h_notloose_pt->Write();
3451   h_notloose_pt->Draw();
3452   c.SaveAs(DIR + type + "_notloose.pdf");
3453 
3454   // fraction of DUPLICATE tracks (genuine and not matched)
3455   TH1F* h_duplicatefrac_pt = (TH1F*)h_trk_duplicate_vspt->Clone();
3456   h_duplicatefrac_pt->SetName("duplicatefrac_pt");
3457   h_duplicatefrac_pt->GetYaxis()->SetTitle("Duplicate fraction");
3458   h_duplicatefrac_pt->Divide(h_trk_duplicate_vspt, h_trk_all_vspt, 1.0, 1.0, "B");
3459 
3460   h_duplicatefrac_pt->Write();
3461   h_duplicatefrac_pt->Draw();
3462   c.SaveAs(DIR + type + "_duplicatefrac.pdf");
3463 
3464   // ---------------------------------------------------------------------------------------------------------
3465   // total track rates vs pt
3466 
3467   h_trk_all_vspt->Scale(1.0 / nevt);
3468   h_trk_loose_vspt->Scale(1.0 / nevt);
3469   h_trk_genuine_vspt->Scale(1.0 / nevt);
3470   h_trk_notloose_vspt->Scale(1.0 / nevt);
3471   h_trk_notgenuine_vspt->Scale(1.0 / nevt);
3472   h_trk_duplicate_vspt->Scale(1.0 / nevt);
3473   h_tp_vspt->Scale(1.0 / nevt);
3474 
3475   h_tp_vspt->GetYaxis()->SetTitle("Tracks / event");
3476   h_tp_vspt->GetXaxis()->SetTitle("Track p_{T} [GeV]");
3477   h_tp_vspt->SetLineColor(4);
3478   h_tp_vspt->SetLineStyle(2);
3479 
3480   h_trk_notgenuine_vspt->SetLineColor(2);
3481   h_trk_notgenuine_vspt->SetLineStyle(1);
3482 
3483   h_trk_duplicate_vspt->SetLineColor(8);
3484   h_trk_duplicate_vspt->SetLineStyle(2);
3485 
3486   float max = h_tp_vspt->GetMaximum();
3487   if (h_trk_all_vspt->GetMaximum() > max)
3488     max = h_trk_all_vspt->GetMaximum();
3489   h_tp_vspt->SetAxisRange(0.001, max * 1.5, "Y");
3490 
3491   h_tp_vspt->Draw("hist");
3492   h_trk_all_vspt->Draw("same,hist");
3493   h_tp_vspt->Draw("same,hist");
3494   h_trk_notgenuine_vspt->Draw("same,hist");
3495   //h_trk_duplicate_vspt->Draw("same,hist");
3496 
3497   h_trk_all_vspt->Write();
3498   h_trk_loose_vspt->Write();
3499   h_trk_genuine_vspt->Write();
3500   h_trk_notloose_vspt->Write();
3501   h_trk_notgenuine_vspt->Write();
3502   h_trk_duplicate_vspt->Write();
3503   h_tp_vspt->Write();
3504 
3505   char txt[500];
3506   sprintf(txt, "# tracks/event = %.1f", h_trk_all_vspt->GetSum());
3507   mySmallText(0.5, 0.85, 1, txt);
3508   char txt3[500];
3509   sprintf(txt3, "# TPs(stubs in #geq 4 layers)/");
3510   char txt2[500];
3511   sprintf(txt2, "event = %.1f", h_tp_vspt->GetSum());
3512   mySmallText(0.5, 0.79, 4, txt3);
3513   mySmallText(0.5, 0.74, 4, txt2);
3514 
3515   sprintf(txt, "# !genuine tracks/event = %.1f", h_trk_notgenuine_vspt->GetSum());
3516   mySmallText(0.5, 0.69, 2, txt);
3517   //sprintf(txt,"# duplicates/event = %.1f",h_trk_duplicate_vspt->GetSum());
3518   //mySmallText(0.5,0.64,8,txt);
3519 
3520   c.SaveAs(DIR + type + "_trackrate_vspt.pdf");
3521 
3522   gPad->SetLogy();
3523   c.SaveAs(DIR + type + "_trackrate_vspt_log.pdf");
3524   gPad->SetLogy(0);
3525 
3526   // ---------------------------------------------------------------------------------------------------------
3527   // sum track/ TP pt in jets
3528   /*
3529   if (TP_select_injet > 0) {
3530 
3531     TH1F* h_frac_sumpt_vspt = (TH1F*) h_jet_trk_sumpt_vspt->Clone();
3532     h_frac_sumpt_vspt->SetName("frac_sumpt_vspt");
3533     h_frac_sumpt_vspt->GetYaxis()->SetTitle("L1 sum(p_{T}) / TP sum(p_{T})");
3534     h_frac_sumpt_vspt->Divide(h_jet_trk_sumpt_vspt, h_jet_tp_sumpt_vspt, 1.0, 1.0, "B");
3535 
3536     TH1F* h_frac_sumpt_vseta = (TH1F*) h_jet_trk_sumpt_vseta->Clone();
3537     h_frac_sumpt_vseta->SetName("frac_sumpt_vseta");
3538     h_frac_sumpt_vseta->GetYaxis()->SetTitle("L1 sum(p_{T}) / TP sum(p_{T})");
3539     h_frac_sumpt_vseta->Divide(h_jet_trk_sumpt_vseta, h_jet_tp_sumpt_vseta, 1.0, 1.0, "B");
3540 
3541 
3542     TH1F* h_matchfrac_sumpt_vspt = (TH1F*) h_jet_matchtrk_sumpt_vspt->Clone();
3543     h_matchfrac_sumpt_vspt->SetName("matchfrac_sumpt_vspt");
3544     h_matchfrac_sumpt_vspt->GetYaxis()->SetTitle("Matched L1 sum(p_{T}) / TP sum(p_{T})");
3545     h_matchfrac_sumpt_vspt->Divide(h_jet_matchtrk_sumpt_vspt, h_jet_tp_sumpt_vspt, 1.0, 1.0, "B");
3546 
3547     TH1F* h_matchfrac_sumpt_vseta = (TH1F*) h_jet_matchtrk_sumpt_vseta->Clone();
3548     h_matchfrac_sumpt_vseta->SetName("matchfrac_sumpt_vseta");
3549     h_matchfrac_sumpt_vseta->GetYaxis()->SetTitle("Matched L1 sum(p_{T}) / TP sum(p_{T})");
3550     h_matchfrac_sumpt_vseta->Divide(h_jet_matchtrk_sumpt_vseta, h_jet_tp_sumpt_vseta, 1.0, 1.0, "B");
3551 
3552 
3553     h_frac_sumpt_vspt->Draw();
3554     c.SaveAs(DIR+type+"_sumpt_vspt.pdf");
3555 
3556     h_frac_sumpt_vseta->Draw();
3557     c.SaveAs(DIR+type+"_sumpt_vseta.pdf");
3558 
3559     h_matchfrac_sumpt_vspt->Draw();
3560     c.SaveAs(DIR+type+"_sumpt_match_vspt.pdf");
3561 
3562     h_matchfrac_sumpt_vseta->Draw();
3563     c.SaveAs(DIR+type+"_sumpt_match_vseta.pdf");
3564   }
3565   */
3566 
3567   // nbr tracks per event
3568 
3569   h_ntrk_pt2->Write();
3570   h_ntrk_pt3->Write();
3571   h_ntrk_pt10->Write();
3572 
3573   h_ntrkPerSector_all->Write();
3574   h_ntrkPerSector_pt2->Write();
3575   h_ntrkPerSector_pt3->Write();
3576   h_ntrkPerSector_pt4->Write();
3577 
3578   h_ntrkPerSector_all->Scale(1.0 / nevt);
3579   h_ntrkPerSector_pt2->Scale(1.0 / nevt);
3580   h_ntrkPerSector_pt3->Scale(1.0 / nevt);
3581   h_ntrkPerSector_pt4->Scale(1.0 / nevt);
3582 
3583   h_ntrkPerSector_all->GetYaxis()->SetTitle("Fraction of events");
3584   h_ntrkPerSector_all->GetXaxis()->SetTitle("Max number of transmitted tracks per #phi sector");
3585 
3586   h_ntrkPerSector_all->SetLineColor(1);
3587   h_ntrkPerSector_pt2->SetLineColor(4);
3588   h_ntrkPerSector_pt3->SetLineColor(2);
3589   h_ntrkPerSector_pt4->SetLineColor(8);
3590 
3591   max = h_ntrkPerSector_all->GetMaximum();
3592   h_ntrkPerSector_all->SetAxisRange(0.00001, max * 5, "Y");
3593   h_ntrkPerSector_all->SetAxisRange(0., 100, "X");
3594 
3595   h_ntrkPerSector_all->Draw("hist");
3596   h_ntrkPerSector_pt2->Draw("same,hist");
3597   h_ntrkPerSector_pt3->Draw("same,hist");
3598   h_ntrkPerSector_pt4->Draw("same,hist");
3599   gPad->SetLogy();
3600 
3601   TLegend* l = new TLegend(0.6, 0.55, 0.85, 0.85);
3602   l->SetFillStyle(0);
3603   l->SetBorderSize(0);
3604   l->SetTextSize(0.04);
3605   l->AddEntry(h_ntrkPerSector_all, "no p_{T}cut", "l");
3606   l->AddEntry(h_ntrkPerSector_pt2, "p_{T}^{track} > 2 GeV", "l");
3607   l->AddEntry(h_ntrkPerSector_pt3, "p_{T}^{track} > 3 GeV", "l");
3608   l->AddEntry(h_ntrkPerSector_pt4, "p_{T}^{track} > 4 GeV", "l");
3609   l->SetTextFont(42);
3610   l->Draw();
3611 
3612   c.SaveAs(DIR + type + "_trackRatePerPhiSector_log.pdf");
3613   gPad->SetLogy(0);
3614 
3615   h_ntrk_genuine_pt2->Write();
3616   h_ntrk_genuine_pt3->Write();
3617   h_ntrk_genuine_pt10->Write();
3618 
3619   if (doDetailedPlots) {
3620     h_ntrk_pt2->Draw();
3621     c.SaveAs(DIR + type + "_trackrate_pt2_perevt.pdf");
3622 
3623     h_ntrk_pt3->Draw();
3624     c.SaveAs(DIR + type + "_trackrate_pt3_perevt.pdf");
3625 
3626     h_ntrk_pt10->Draw();
3627     c.SaveAs(DIR + type + "_trackrate_pt10_perevt.pdf");
3628   }
3629 
3630   // number of tracks vs. eta, pT (trk_eta/trk_pt)
3631 
3632   if (doDetailedPlots) {
3633     h_trk_eta->Write();
3634     h_trk_pt->Write();
3635 
3636     h_trk_eta->Draw();
3637     c.SaveAs(DIR + type + "_trk_eta.pdf");
3638     h_trk_pt->Draw();
3639     c.SaveAs(DIR + type + "_trk_pt.pdf");
3640   }
3641 
3642   // Sample z0 resolution at a couple of rapidity points.
3643   float etaSample1 = h2_resVsEta_z0_68->GetXaxis()->GetBinCenter(1);
3644   float etaSample2 = h2_resVsEta_z0_68->GetXaxis()->GetBinCenter(0.8 * nETARANGE);
3645   float z0ResSample1 = h2_resVsEta_z0_68->GetBinContent(1);
3646   float z0ResSample2 = h2_resVsEta_z0_68->GetBinContent(0.8 * nETARANGE);
3647 
3648   fout->Close();
3649 
3650   // ---------------------------------------------------------------------------------------------------------
3651   //some printouts
3652 
3653   cout << "number of events = " << nevt << endl;
3654 
3655   float k = (float)n_match_eta1p0;
3656   float N = (float)n_all_eta1p0;
3657   if (std::abs(N) > 0)
3658     cout << endl
3659          << "efficiency for |eta| < 1.0 = " << k / N * 100.0 << " +- " << 1.0 / N * sqrt(k * (1.0 - k / N)) * 100.0
3660          << endl;
3661   k = (float)n_match_eta1p75;
3662   N = (float)n_all_eta1p75;
3663   if (std::abs(N) > 0)
3664     cout << "efficiency for 1.0 < |eta| < 1.75 = " << k / N * 100.0 << " +- "
3665          << 1.0 / N * sqrt(k * (1.0 - k / N)) * 100.0 << endl;
3666   k = (float)n_match_eta2p5;
3667   N = (float)n_all_eta2p5;
3668   if (std::abs(N) > 0)
3669     cout << "efficiency for 1.75 < |eta| < " << std::min(TP_maxEta, 2.5f) << " = " << k / N * 100.0 << " +- "
3670          << 1.0 / N * sqrt(k * (1.0 - k / N)) * 100.0 << endl;
3671   N = (float)n_all_eta1p0 + n_all_eta1p75 + n_all_eta2p5;
3672   k = (float)n_match_eta1p0 + n_match_eta1p75 + n_match_eta2p5;
3673   if (std::abs(N) > 0)
3674     cout << "combined efficiency for |eta| < " << std::min(TP_maxEta, 2.5f) << " = " << k / N * 100.0 << " +- "
3675          << 1.0 / N * sqrt(k * (1.0 - k / N)) * 100.0 << " = " << k << "/" << N << endl
3676          << endl;
3677 
3678   k = (float)n_match_ptg2;
3679   N = (float)n_all_ptg2;
3680   if (std::abs(N) > 0)
3681     cout << "efficiency for pt > " << std::max(TP_minPt, 2.0f) << " = " << k / N * 100.0 << " +- "
3682          << 1.0 / N * sqrt(k * (1.0 - k / N)) * 100.0 << endl;
3683   k = (float)n_match_pt2to8;
3684   N = (float)n_all_pt2to8;
3685   if (std::abs(N) > 0)
3686     cout << "efficiency for " << std::max(TP_minPt, 2.0f) << " < pt < 8.0 = " << k / N * 100.0 << " +- "
3687          << 1.0 / N * sqrt(k * (1.0 - k / N)) * 100.0 << endl;
3688   k = (float)n_match_ptg8;
3689   N = (float)n_all_ptg8;
3690   if (std::abs(N) > 0)
3691     cout << "efficiency for pt > 8.0 = " << k / N * 100.0 << " +- " << 1.0 / N * sqrt(k * (1.0 - k / N)) * 100.0
3692          << endl;
3693   k = (float)n_match_ptg40;
3694   N = (float)n_all_ptg40;
3695   if (std::abs(N) > 0)
3696     cout << "efficiency for pt > 40.0 = " << k / N * 100.0 << " +- " << 1.0 / N * sqrt(k * (1.0 - k / N)) * 100.0
3697          << endl
3698          << endl;
3699 
3700   // track rates
3701   cout << "# TP/event (pt > " << std::max(TP_minPt, 2.0f) << ") = " << (float)ntp_pt2 / nevt << endl;
3702   cout << "# TP/event (pt > 3.0) = " << (float)ntp_pt3 / nevt << endl;
3703   cout << "# TP/event (pt > 10.0) = " << (float)ntp_pt10 / nevt << endl;
3704 
3705   cout << "# tracks/event (no pt cut)= " << (float)ntrk / nevt << endl;
3706   cout << "# tracks/event (pt > " << std::max(TP_minPt, 2.0f) << ") = " << (float)ntrk_pt2 / nevt << endl;
3707   cout << "# tracks/event (pt > 3.0) = " << (float)ntrk_pt3 / nevt << endl;
3708   cout << "# tracks/event (pt > 10.0) = " << (float)ntrk_pt10 / nevt << endl << endl;
3709 
3710   // fake track rate
3711   if (ntrk_genuine_pt2 > 0) {
3712     cout << "Percentage fake tracks (pt > " << std::max(TP_minPt, 2.0f)
3713          << ") = " << 100. * (1. - float(ntrk_genuine_pt2) / float(ntrk_pt2)) << "%" << endl;
3714     cout << "Percentage duplicate tracks (pt > " << std::max(TP_minPt, 2.0f)
3715          << ")= " << 100. * float(ntp_nmatch) / float(ntrk_pt2) << "%" << endl
3716          << endl;
3717   }
3718 
3719   // z0 resolution
3720   cout << "z0 resolution = " << z0ResSample1 << "cm at |eta| = " << etaSample1 << endl;
3721   cout << "z0 resolution = " << z0ResSample2 << "cm at |eta| = " << etaSample2 << endl;
3722 }
3723 
3724 void SetPlotStyle() {
3725   // from ATLAS plot style macro
3726 
3727   // use plain black on white colors
3728   gStyle->SetFrameBorderMode(0);
3729   gStyle->SetFrameFillColor(0);
3730   gStyle->SetCanvasBorderMode(0);
3731   gStyle->SetCanvasColor(0);
3732   gStyle->SetPadBorderMode(0);
3733   gStyle->SetPadColor(0);
3734   gStyle->SetStatColor(0);
3735   gStyle->SetHistLineColor(1);
3736 
3737   gStyle->SetPalette(1);
3738 
3739   // set the paper & margin sizes
3740   gStyle->SetPaperSize(20, 26);
3741   gStyle->SetPadTopMargin(0.05);
3742   gStyle->SetPadRightMargin(0.05);
3743   gStyle->SetPadBottomMargin(0.16);
3744   gStyle->SetPadLeftMargin(0.16);
3745 
3746   // set title offsets (for axis label)
3747   gStyle->SetTitleXOffset(1.4);
3748   gStyle->SetTitleYOffset(1.4);
3749 
3750   // use large fonts
3751   gStyle->SetTextFont(42);
3752   gStyle->SetTextSize(0.05);
3753   gStyle->SetLabelFont(42, "x");
3754   gStyle->SetTitleFont(42, "x");
3755   gStyle->SetLabelFont(42, "y");
3756   gStyle->SetTitleFont(42, "y");
3757   gStyle->SetLabelFont(42, "z");
3758   gStyle->SetTitleFont(42, "z");
3759   gStyle->SetLabelSize(0.05, "x");
3760   gStyle->SetTitleSize(0.05, "x");
3761   gStyle->SetLabelSize(0.05, "y");
3762   gStyle->SetTitleSize(0.05, "y");
3763   gStyle->SetLabelSize(0.05, "z");
3764   gStyle->SetTitleSize(0.05, "z");
3765 
3766   // use bold lines and markers
3767   gStyle->SetMarkerStyle(20);
3768   gStyle->SetMarkerSize(1.2);
3769   gStyle->SetHistLineWidth(2.);
3770   gStyle->SetLineStyleString(2, "[12 12]");
3771 
3772   // get rid of error bar caps
3773   gStyle->SetEndErrorSize(0.);
3774 
3775   // do not display any of the standard histogram decorations
3776   gStyle->SetOptTitle(0);
3777   gStyle->SetOptStat(0);
3778   gStyle->SetOptFit(0);
3779 
3780   // put tick marks on top and RHS of plots
3781   gStyle->SetPadTickX(1);
3782   gStyle->SetPadTickY(1);
3783 }
3784 
3785 void mySmallText(Double_t x, Double_t y, Color_t color, char* text) {
3786   Double_t tsize = 0.044;
3787   TLatex l;
3788   l.SetTextSize(tsize);
3789   l.SetNDC();
3790   l.SetTextColor(color);
3791   l.DrawLatex(x, y, text);
3792 }
3793 
3794 double getIntervalContainingFractionOfEntries(TH1* absResidualHistogram, double quantileToCalculate, int minEntries) {
3795   double totalIntegral = absResidualHistogram->Integral(0, absResidualHistogram->GetNbinsX() + 1);
3796   double numEntries = absResidualHistogram->GetEntries();
3797 
3798   // Check that the interval is not somewhere in the overflow bin
3799   double maxAllowedEntriesInOverflow = totalIntegral * (1 - quantileToCalculate);
3800   double nEntriesInOverflow = absResidualHistogram->GetBinContent(absResidualHistogram->GetNbinsX() + 1);
3801   if (nEntriesInOverflow > maxAllowedEntriesInOverflow) {
3802     // cout << "WARNING : Cannot compute range corresponding to interval, as it is in the overflow bin" << endl;
3803     return absResidualHistogram->GetXaxis()->GetXmax() * 1.2;
3804   }
3805 
3806   // Calculate quantile for given interval
3807   double interval[1];
3808   double quantile[1] = {quantileToCalculate};
3809   if (totalIntegral > 0 && numEntries >= minEntries) {
3810     absResidualHistogram->GetQuantiles(1, interval, quantile);
3811   } else {
3812     cout << "WARNING: histo " << absResidualHistogram->GetName()
3813          << " empty or with too few entries, so can't calc quantiles." << endl;
3814     interval[0] = 0.;
3815   }
3816 
3817   return interval[0];
3818 }
3819 
3820 void makeResidualIntervalPlot(
3821     TString type, TString dir, TString variable, TH1F* h_68, TH1F* h_90, TH1F* h_99, double minY, double maxY) {
3822   TCanvas c;
3823 
3824   h_68->SetMinimum(minY);
3825   h_90->SetMinimum(minY);
3826   h_99->SetMinimum(minY);
3827 
3828   h_68->SetMaximum(maxY);
3829   h_90->SetMaximum(maxY);
3830   h_99->SetMaximum(maxY);
3831 
3832   h_68->SetMarkerStyle(20);
3833   h_90->SetMarkerStyle(26);
3834   h_99->SetMarkerStyle(24);
3835 
3836   h_68->Draw("P");
3837   h_68->Write();
3838   h_90->Draw("P same");
3839   h_90->Write();
3840   h_99->Draw("P same");
3841   h_99->Write();
3842 
3843   TLegend* l = new TLegend(0.65, 0.65, 0.85, 0.85);
3844   l->SetFillStyle(0);
3845   l->SetBorderSize(0);
3846   l->SetTextSize(0.04);
3847   l->AddEntry(h_99, "99%", "p");
3848   l->AddEntry(h_90, "90%", "p");
3849   l->AddEntry(h_68, "68%", "p");
3850   l->SetTextFont(42);
3851   l->Draw();
3852 
3853   c.SaveAs(dir + type + "_" + variable + "_interval.pdf");
3854 
3855   delete l;
3856 }