Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 11:58:54

0001 ///////////////////////////////////////////////////////////////////////////////
0002 //
0003 //   HBHEMuonOfflineSimAnalyzer h1(infile, outfile, mode, maxDHB, maxDHE);
0004 //   h1.Loop()
0005 //
0006 //      Offline analysis for MC files
0007 //
0008 //   infile     const char*  Name of the input file
0009 //   outfile    const char*  Name of the output file
0010 //                           (dyll_PU20_25_output.root)
0011 //   mode       int          Geometry file used 0:(defined by maxDHB/HE);
0012 //                           1 (Run 1; valid till 2016); 2 (Run 2; 2018);
0013 //                           3 (Run 3; post LS2); 4 (2017 Plan 1);
0014 //                           5 (Run 4; post LS3); default (3)
0015 //   maxDHB     int          Maximum number of depths for HB (4)
0016 //   maxDHE     int          Maximum number of depths for HE (7)
0017 //
0018 ///////////////////////////////////////////////////////////////////////////////
0019 
0020 #include <algorithm>
0021 #include <cmath>
0022 #include <iostream>
0023 #include <vector>
0024 #include <cstring>
0025 #include <TCanvas.h>
0026 #include <TChain.h>
0027 #include <TFile.h>
0028 #include <TH1.h>
0029 #include <TH1D.h>
0030 #include <TH2.h>
0031 #include <TH2D.h>
0032 #include <TProfile.h>
0033 #include <TROOT.h>
0034 #include <TStyle.h>
0035 #include <TTree.h>
0036 
0037 class HBHEMuonOfflineSimAnalyzer {
0038 private:
0039   TTree *fChain;  //!pointer to the analyzed TTree or TChain
0040   Int_t fCurrent;
0041 
0042   UInt_t Run_No;
0043   UInt_t Event_No;
0044   UInt_t LumiNumber;
0045   UInt_t BXNumber;
0046   double pt_of_muon;
0047   double eta_of_muon;
0048   double phi_of_muon;
0049   double p_of_muon;
0050   double ecal_3x3;
0051   unsigned int ecal_detID;
0052   double hcal_1x1;
0053   double matchedId;
0054   unsigned int hcal_detID;
0055   unsigned int hcal_cellHot;
0056   double activeLength;
0057   double hcal_edepth1;
0058   double hcal_edepth2;
0059   double hcal_edepth3;
0060   double hcal_edepth4;
0061   double hcal_activeL1;
0062   double hcal_activeL2;
0063   double hcal_activeL3;
0064   double hcal_activeL4;
0065   double activeLengthHot;
0066   double hcal_edepthHot1;
0067   double hcal_edepthHot2;
0068   double hcal_edepthHot3;
0069   double hcal_edepthHot4;
0070   double hcal_activeHotL1;
0071   double hcal_activeHotL2;
0072   double hcal_activeHotL3;
0073   double hcal_activeHotL4;
0074   double hcal_edepth5;
0075   double hcal_activeL5;
0076   double hcal_edepthHot5;
0077   double hcal_activeHotL5;
0078   double hcal_edepth6;
0079   double hcal_activeL6;
0080   double hcal_edepthHot6;
0081   double hcal_activeHotL6;
0082   double hcal_edepth7;
0083   double hcal_activeL7;
0084   double hcal_edepthHot7;
0085   double hcal_activeHotL7;
0086 
0087   TBranch *b_Run_No;            //!
0088   TBranch *b_Event_No;          //!
0089   TBranch *b_LumiNumber;        //!
0090   TBranch *b_BXNumber;          //!
0091   TBranch *b_pt_of_muon;        //!
0092   TBranch *b_eta_of_muon;       //!
0093   TBranch *b_phi_of_muon;       //!
0094   TBranch *b_p_of_muon;         //!
0095   TBranch *b_ecal_3x3;          //!
0096   TBranch *b_ecal_detID;        //!
0097   TBranch *b_hcal_1x1;          //!
0098   TBranch *b_hcal_detID;        //!
0099   TBranch *b_hcal_cellHot;      //!
0100   TBranch *b_activeLength;      //!
0101   TBranch *b_hcal_edepth1;      //!
0102   TBranch *b_hcal_edepth2;      //!
0103   TBranch *b_hcal_edepth3;      //!
0104   TBranch *b_hcal_edepth4;      //!
0105   TBranch *b_hcal_activeL1;     //!
0106   TBranch *b_hcal_activeL2;     //!
0107   TBranch *b_hcal_activeL3;     //!
0108   TBranch *b_hcal_activeL4;     //!
0109   TBranch *b_activeLengthHot;   //!
0110   TBranch *b_hcal_edepthHot1;   //!
0111   TBranch *b_hcal_edepthHot2;   //!
0112   TBranch *b_hcal_edepthHot3;   //!
0113   TBranch *b_hcal_edepthHot4;   //!
0114   TBranch *b_hcal_activeHotL1;  //!
0115   TBranch *b_hcal_activeHotL2;  //!
0116   TBranch *b_hcal_activeHotL3;  //!
0117   TBranch *b_hcal_activeHotL4;  //!
0118   TBranch *b_hcal_edepth5;      //!
0119   TBranch *b_hcal_activeL5;     //!
0120   TBranch *b_hcal_edepthHot5;   //!
0121   TBranch *b_hcal_activeHotL5;  //!
0122   TBranch *b_hcal_edepth6;      //!
0123   TBranch *b_hcal_activeL6;     //!
0124   TBranch *b_hcal_edepthHot6;   //!
0125   TBranch *b_hcal_activeHotL6;  //!
0126   TBranch *b_hcal_edepth7;      //!
0127   TBranch *b_hcal_activeL7;     //!
0128   TBranch *b_hcal_edepthHot7;   //!
0129   TBranch *b_hcal_activeHotL7;  //!
0130   TBranch *b_matchedId;         //!
0131 
0132 public:
0133   HBHEMuonOfflineSimAnalyzer(const char *infile,
0134                              const char *outfile = "dyll_PU20_25_output.root",
0135                              const int mode = 3,
0136                              const int maxDHB = 4,
0137                              const int maxDHE = 7);
0138   virtual ~HBHEMuonOfflineSimAnalyzer();
0139   virtual Int_t Cut(Long64_t entry);
0140   virtual Int_t GetEntry(Long64_t entry);
0141   virtual Long64_t LoadTree(Long64_t entry);
0142   virtual void Init(TTree *tree);
0143   virtual void Loop();
0144   virtual Bool_t Notify();
0145   virtual void Show(Long64_t entry = -1);
0146 
0147   std::vector<std::string> firedTriggers;
0148   void BookHistograms(const char *);
0149   void WriteHistograms();
0150   bool LooseMuon();
0151   bool tightMuon();
0152   bool SoftMuon();
0153   void etaPhiHcal(unsigned int detId, int &eta, int &phi, int &depth);
0154   void etaPhiEcal(unsigned int detId, int &type, int &zside, int &etaX, int &phiY, int &plane, int &strip);
0155   void calculateP(double pt, double eta, double &pM);
0156   void close();
0157   int NDepthBins(int ieta, int iphi);
0158   int NPhiBins(int ieta);
0159 
0160 private:
0161   static const bool debug_ = false;
0162   static const int maxDep = 7;
0163   static const int maxEta = 29;
0164   static const int maxPhi = 72;
0165   //3x16x72x2 + 5x4x72x2 + 5x9x36x2
0166   static const int maxHist = 20000;  //13032;
0167   int modeLHC_, maxDepthHB_, maxDepthHE_, maxDepth_;
0168   int nHist, nDepths[maxEta], nDepthsPhi[maxEta], indxEta[maxEta][maxDep][maxPhi];
0169   TFile *output_file;
0170 
0171   TH1D *h_Pt_Muon[3], *h_Eta_Muon[3], *h_Phi_Muon[3], *h_P_Muon[3];
0172   TH1D *h_PF_Muon[3], *h_GlobTrack_Chi[3], *h_Global_Muon_Hits[3];
0173   TH1D *h_MatchedStations[3], *h_Tight_TransImpactparameter[3];
0174   TH1D *h_Tight_LongitudinalImpactparameter[3], *h_InnerTrackPixelHits[3];
0175   TH1D *h_TrackerLayer[3], *h_IsolationR04[3], *h_Global_Muon[3];
0176   TH1D *h_LongImpactParameter[3], *h_LongImpactParameterBin1[3], *h_LongImpactParameterBin2[3];
0177 
0178   TH1D *h_TransImpactParameter[3], *h_TransImpactParameterBin1[3], *h_TransImpactParameterBin2[3];
0179   TH1D *h_Hot_MuonEnergy_hcal_ClosestCell[3][maxHist], *h_Hot_MuonEnergy_hcal_HotCell[3][maxHist],
0180       *h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[3][maxHist], *h_HotCell_MuonEnergy_phi[3][maxHist];
0181   TH2D *h_2D_Bin1[3], *h_2D_Bin2[3];
0182   TH1D *h_ecal_energy[3], *h_hcal_energy[3], *h_3x3_ecal[3], *h_1x1_hcal[3];
0183   TH1D *h_MuonHittingEcal[3], *h_HotCell[3], *h_MuonEnergy_hcal[3][maxHist];
0184   TH1D *h_Hot_MuonEnergy_hcal[3][maxHist];
0185   TH2D *hcal_ietaVsEnergy[3];
0186   TProfile *h_EtaX_hcal[3], *h_PhiY_hcal[3], *h_EtaX_ecal[3], *h_PhiY_ecal[3];
0187   TProfile *h_Eta_ecal[3], *h_Phi_ecal[3];
0188   TProfile *h_MuonEnergy_eta[3][maxDep], *h_MuonEnergy_phi[3][maxDep], *h_MuonEnergy_muon_eta[3][maxDep];
0189   TProfile *h_Hot_MuonEnergy_eta[3][maxDep], *h_Hot_MuonEnergy_phi[3][maxDep], *h_Hot_MuonEnergy_muon_eta[3][maxDep];
0190   TProfile *h_IsoHot_MuonEnergy_eta[3][maxDep], *h_IsoHot_MuonEnergy_phi[3][maxDep],
0191       *h_IsoHot_MuonEnergy_muon_eta[3][maxDep];
0192   TProfile *h_IsoWithoutHot_MuonEnergy_eta[3][maxDep], *h_IsoWithoutHot_MuonEnergy_phi[3][maxDep],
0193       *h_IsoWithoutHot_MuonEnergy_muon_eta[3][maxDep];
0194   TProfile *h_HotWithoutIso_MuonEnergy_eta[3][maxDep], *h_HotWithoutIso_MuonEnergy_phi[3][maxDep],
0195       *h_HotWithoutIso_MuonEnergy_muon_eta[3][maxDep];
0196 };
0197 
0198 HBHEMuonOfflineSimAnalyzer::HBHEMuonOfflineSimAnalyzer(
0199     const char *infile, const char *outFileName, const int mode, const int maxDHB, const int maxDHE) {
0200   modeLHC_ = mode;
0201   maxDepthHB_ = maxDHB;
0202   maxDepthHE_ = maxDHE;
0203   maxDepth_ = (maxDepthHB_ > maxDepthHE_) ? maxDepthHB_ : maxDepthHE_;
0204   // if parameter tree is not specified (or zero), connect the file
0205   // used to generate this class and read the Tree.
0206   TFile *f = new TFile(infile);
0207   TDirectory *dir = (TDirectory *)f->Get("HcalHBHEMuonAnalyzer");
0208   TTree *tree(0);
0209   dir->GetObject("TREE", tree);
0210   Init(tree);
0211 
0212   //Now book histograms
0213   BookHistograms(outFileName);
0214 }
0215 
0216 HBHEMuonOfflineSimAnalyzer::~HBHEMuonOfflineSimAnalyzer() {
0217   if (!fChain)
0218     return;
0219   delete fChain->GetCurrentFile();
0220 }
0221 
0222 Int_t HBHEMuonOfflineSimAnalyzer::Cut(Long64_t) {
0223   // This function may be called from Loop.
0224   // returns  1 if entry is accepted.
0225   // returns -1 otherwise.
0226   return 1;
0227 }
0228 
0229 Int_t HBHEMuonOfflineSimAnalyzer::GetEntry(Long64_t entry) {
0230   // Read contents of entry.
0231   if (!fChain)
0232     return 0;
0233   return fChain->GetEntry(entry);
0234 }
0235 
0236 Long64_t HBHEMuonOfflineSimAnalyzer::LoadTree(Long64_t entry) {
0237   // Set the environment to read one entry
0238   if (!fChain)
0239     return -5;
0240   Long64_t centry = fChain->LoadTree(entry);
0241   if (centry < 0)
0242     return centry;
0243   if (fChain->GetTreeNumber() != fCurrent) {
0244     fCurrent = fChain->GetTreeNumber();
0245     Notify();
0246   }
0247   return centry;
0248 }
0249 
0250 void HBHEMuonOfflineSimAnalyzer::Init(TTree *tree) {
0251   // The Init() function is called when the selector needs to initialize
0252   // a new tree or chain. Typically here the branch addresses and branch
0253   // pointers of the tree will be set.
0254   // It is normally not necessary to make changes to the generated
0255   // code, but the routine can be extended by the user if needed.
0256   // Init() will be called many times when running on PROOF
0257   // (once per file to be processed).
0258 
0259   // Set object pointer
0260   pt_of_muon = 0;
0261   eta_of_muon = 0;
0262   phi_of_muon = 0;
0263   p_of_muon = 0;
0264   ecal_3x3 = 0;
0265   ecal_detID = 0;
0266   hcal_1x1 = 0;
0267   hcal_detID = 0;
0268   hcal_cellHot = 0;
0269   activeLength = 0;
0270   hcal_edepth1 = 0;
0271   hcal_edepth2 = 0;
0272   hcal_edepth3 = 0;
0273   hcal_edepth4 = 0;
0274   hcal_activeL1 = 0;
0275   hcal_activeL2 = 0;
0276   hcal_activeL3 = 0;
0277   hcal_activeL4 = 0;
0278   activeLengthHot = 0;
0279   hcal_edepthHot1 = 0;
0280   hcal_edepthHot2 = 0;
0281   hcal_edepthHot3 = 0;
0282   hcal_edepthHot4 = 0;
0283   hcal_activeHotL1 = 0;
0284   hcal_activeHotL2 = 0;
0285   hcal_activeHotL3 = 0;
0286   hcal_activeHotL4 = 0;
0287   hcal_edepth5 = 0;
0288   hcal_activeL5 = 0;
0289   hcal_edepthHot5 = 0;
0290   hcal_activeHotL5 = 0;
0291   hcal_edepth6 = 0;
0292   hcal_activeL6 = 0;
0293   hcal_edepthHot6 = 0;
0294   hcal_activeHotL6 = 0;
0295   hcal_edepth7 = 0;
0296   hcal_activeL7 = 0;
0297   hcal_edepthHot7 = 0;
0298   hcal_activeHotL7 = 0;
0299   matchedId = 0;
0300   if (!tree)
0301     return;
0302   fChain = tree;
0303   fCurrent = -1;
0304   fChain->SetMakeClass(1);
0305 
0306   fChain->SetBranchAddress("Run_No", &Run_No, &b_Run_No);
0307   fChain->SetBranchAddress("Event_No", &Event_No, &b_Event_No);
0308   fChain->SetBranchAddress("LumiNumber", &LumiNumber, &b_LumiNumber);
0309   fChain->SetBranchAddress("BXNumber", &BXNumber, &b_BXNumber);
0310   fChain->SetBranchAddress("pt_of_muon", &pt_of_muon, &b_pt_of_muon);
0311   fChain->SetBranchAddress("eta_of_muon", &eta_of_muon, &b_eta_of_muon);
0312   fChain->SetBranchAddress("phi_of_muon", &phi_of_muon, &b_phi_of_muon);
0313   fChain->SetBranchAddress("p_of_muon", &p_of_muon, &b_p_of_muon);
0314   fChain->SetBranchAddress("ecal_3x3", &ecal_3x3, &b_ecal_3x3);
0315   fChain->SetBranchAddress("ecal_detID", &ecal_detID, &b_ecal_detID);
0316   fChain->SetBranchAddress("hcal_1x1", &hcal_1x1, &b_hcal_1x1);
0317   fChain->SetBranchAddress("matchedId", &matchedId, &b_matchedId);
0318   fChain->SetBranchAddress("hcal_detID", &hcal_detID, &b_hcal_detID);
0319   fChain->SetBranchAddress("hcal_cellHot", &hcal_cellHot, &b_hcal_cellHot);
0320   fChain->SetBranchAddress("activeLength", &activeLength, &b_activeLength);
0321   fChain->SetBranchAddress("hcal_edepth1", &hcal_edepth1, &b_hcal_edepth1);
0322   fChain->SetBranchAddress("hcal_edepth2", &hcal_edepth2, &b_hcal_edepth2);
0323   fChain->SetBranchAddress("hcal_edepth3", &hcal_edepth3, &b_hcal_edepth3);
0324   fChain->SetBranchAddress("hcal_edepth4", &hcal_edepth4, &b_hcal_edepth4);
0325   fChain->SetBranchAddress("hcal_edepth5", &hcal_edepth5, &b_hcal_edepth5);
0326   fChain->SetBranchAddress("hcal_edepth6", &hcal_edepth6, &b_hcal_edepth6);
0327   fChain->SetBranchAddress("hcal_edepth7", &hcal_edepth7, &b_hcal_edepth7);
0328   fChain->SetBranchAddress("hcal_activeL1", &hcal_activeL1, &b_hcal_activeL1);
0329   fChain->SetBranchAddress("hcal_activeL2", &hcal_activeL2, &b_hcal_activeL2);
0330   fChain->SetBranchAddress("hcal_activeL3", &hcal_activeL3, &b_hcal_activeL3);
0331   fChain->SetBranchAddress("hcal_activeL4", &hcal_activeL4, &b_hcal_activeL4);
0332   fChain->SetBranchAddress("hcal_activeL5", &hcal_activeL5, &b_hcal_activeL5);
0333   fChain->SetBranchAddress("hcal_activeL6", &hcal_activeL6, &b_hcal_activeL6);
0334   fChain->SetBranchAddress("hcal_activeL7", &hcal_activeL7, &b_hcal_activeL7);
0335   fChain->SetBranchAddress("activeLengthHot", &activeLengthHot, &b_activeLengthHot);
0336   fChain->SetBranchAddress("hcal_edepthHot1", &hcal_edepthHot1, &b_hcal_edepthHot1);
0337   fChain->SetBranchAddress("hcal_edepthHot2", &hcal_edepthHot2, &b_hcal_edepthHot2);
0338   fChain->SetBranchAddress("hcal_edepthHot3", &hcal_edepthHot3, &b_hcal_edepthHot3);
0339   fChain->SetBranchAddress("hcal_edepthHot4", &hcal_edepthHot4, &b_hcal_edepthHot4);
0340   fChain->SetBranchAddress("hcal_edepthHot5", &hcal_edepthHot5, &b_hcal_edepthHot5);
0341   fChain->SetBranchAddress("hcal_edepthHot6", &hcal_edepthHot6, &b_hcal_edepthHot6);
0342   fChain->SetBranchAddress("hcal_edepthHot7", &hcal_edepthHot7, &b_hcal_edepthHot7);
0343   fChain->SetBranchAddress("hcal_activeHotL1", &hcal_activeHotL1, &b_hcal_activeHotL1);
0344   fChain->SetBranchAddress("hcal_activeHotL2", &hcal_activeHotL2, &b_hcal_activeHotL2);
0345   fChain->SetBranchAddress("hcal_activeHotL3", &hcal_activeHotL3, &b_hcal_activeHotL3);
0346   fChain->SetBranchAddress("hcal_activeHotL4", &hcal_activeHotL4, &b_hcal_activeHotL4);
0347   fChain->SetBranchAddress("hcal_activeHotL5", &hcal_activeHotL5, &b_hcal_activeHotL5);
0348   fChain->SetBranchAddress("hcal_activeHotL6", &hcal_activeHotL6, &b_hcal_activeHotL6);
0349   fChain->SetBranchAddress("hcal_activeHotL7", &hcal_activeHotL7, &b_hcal_activeHotL7);
0350   Notify();
0351 }
0352 
0353 void HBHEMuonOfflineSimAnalyzer::Loop() {
0354   //declarations
0355   if (fChain == 0)
0356     return;
0357 
0358   Long64_t nentries = fChain->GetEntriesFast();
0359 
0360   if (debug_)
0361     std::cout << "nevent = " << nentries << std::endl;
0362 
0363   Long64_t nbytes = 0, nb = 0;
0364 
0365   for (Long64_t jentry = 0; jentry < nentries; jentry++) {
0366     Long64_t ientry = LoadTree(jentry);
0367     if (ientry < 0)
0368       break;
0369     nb = fChain->GetEntry(jentry);
0370     nbytes += nb;
0371 
0372     if (debug_) {
0373       std::cout << "ecal_det_id " << ecal_detID << std::endl;
0374       std::cout << "hcal_det_id " << std::hex << hcal_detID << std::dec;
0375     }
0376     int etaHcal, phiHcal, depthHcal;
0377     etaPhiHcal(hcal_detID, etaHcal, phiHcal, depthHcal);
0378 
0379     int eta = (etaHcal > 0) ? (etaHcal - 1) : -(1 + etaHcal);
0380     int nDepth = NDepthBins(eta + 1, phiHcal);
0381     int nPhi = NPhiBins(eta + 1);
0382 
0383     double phiYHcal = (phiHcal - 0.5);
0384     if (debug_)
0385       std::cout << "phiHcal" << phiHcal << " phiYHcal" << phiYHcal << std::endl;
0386 
0387     for (int cut = 0; cut < 3; ++cut) {
0388       bool select(false);
0389       if (cut == 0)
0390         select = tightMuon();
0391       else if (cut == 1)
0392         select = SoftMuon();
0393       else
0394         select = LooseMuon();
0395 
0396       if (select) {
0397         //    h_P_Muon[cut]->Fill(p_of_muon);
0398         h_P_Muon[cut]->Fill(p_of_muon);
0399         h_Pt_Muon[cut]->Fill(pt_of_muon);
0400         h_Eta_Muon[cut]->Fill(eta_of_muon);
0401 
0402         double energyFill;
0403         for (int dep = 0; dep < nDepth; ++dep) {
0404           if (debug_) {
0405             std::cout << "why on 15/2 only" << std::endl;
0406             std::cout << "dep:" << dep << std::endl;
0407           }
0408           int PHI = (nPhi > 36) ? (phiHcal - 1) : (phiHcal - 1) / 2;
0409           double en1(-9999), en2(-9999);
0410           if (dep == 0) {
0411             en1 = hcal_edepth1;
0412             en2 = hcal_edepthHot1;
0413             energyFill = (hcal_activeHotL1 > 0) ? hcal_activeHotL1 : 999;
0414           } else if (dep == 1) {
0415             en1 = hcal_edepth2;
0416             en2 = hcal_edepthHot2;
0417             energyFill = (hcal_activeHotL2 > 0) ? hcal_activeHotL2 : 999;
0418             if (debug_)
0419               std::cout << "problem here.. lets see if it got printed\n";
0420           } else if (dep == 2) {
0421             en1 = hcal_edepth3;
0422             en2 = hcal_edepthHot3;
0423             energyFill = (hcal_activeHotL3 > 0) ? hcal_activeHotL3 : 999;
0424           } else if (dep == 3) {
0425             en1 = hcal_edepth4;
0426             en2 = hcal_edepthHot4;
0427             if (debug_)
0428               std::cout << "Hello in 4" << std::endl;
0429             energyFill = (hcal_activeHotL4 > 0) ? hcal_activeHotL4 : 999;
0430           } else if (dep == 4) {
0431             en1 = hcal_edepth5;
0432             en2 = hcal_edepthHot5;
0433             energyFill = (hcal_activeHotL5 > 0) ? hcal_activeHotL5 : 999;
0434           } else if (dep == 5) {
0435             if (debug_)
0436               std::cout << "Energy in depth 6 " << maxDepth_ << ":" << hcal_edepth6 << ":" << hcal_edepthHot6
0437                         << std::endl;
0438             en1 = (maxDepth_ > 5) ? hcal_edepth6 : 0;
0439             en2 = (maxDepth_ > 5) ? hcal_edepthHot6 : 0;
0440             energyFill = (hcal_activeHotL6 > 0) ? hcal_activeHotL6 : 999;
0441           } else if (dep == 6) {
0442             if (debug_)
0443               std::cout << "Energy in depth 7 " << maxDepth_ << ":" << hcal_edepth7 << ":" << hcal_edepthHot7
0444                         << std::endl;
0445             en1 = (maxDepth_ > 6) ? hcal_edepth7 : 0;
0446             en2 = (maxDepth_ > 6) ? hcal_edepthHot7 : 0;
0447             energyFill = (hcal_activeHotL7 > 0) ? hcal_activeHotL7 : 999;
0448           }
0449 
0450           if (debug_) {
0451             std::cout << " Debug2" << std::endl;
0452             std::cout << "ok1" << en1 << std::endl;
0453             std::cout << "ok2" << en2 << std::endl;
0454           }
0455           bool ok1 = (en1 > -9999);
0456           bool ok2 = (en2 > -9999);
0457 
0458           if (debug_)
0459             std::cout << "Before Index" << std::endl;
0460 
0461           int ind = (etaHcal > 0) ? indxEta[eta][dep][PHI] : 1 + indxEta[eta][dep][PHI];
0462 
0463           if (debug_) {
0464             std::cout << "ieta " << eta << "depth " << dep << "indxEta[eta][dep]:" << indxEta[eta][dep][PHI]
0465                       << std::endl;
0466             std::cout << "index showing eta,depth:" << ind << std::endl;
0467             std::cout << "etaHcal: " << etaHcal << " eta " << eta << " dep " << dep << " indx " << ind << std::endl;
0468           }
0469           if (!(matchedId))
0470             continue;
0471           if (ok1) {
0472             if (debug_)
0473               std::cout << "enter ok1" << std::endl;
0474             if (hcal_cellHot == 1) {
0475               if (en2 > 0) {
0476                 h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[cut][ind]->Fill(en2 / energyFill);
0477               }
0478               if (debug_)
0479                 std::cout << "enter hot cell" << std::endl;
0480             }
0481           }
0482 
0483           if (ok2) {
0484             if (debug_)
0485               std::cout << "enter ok2" << std::endl;
0486             if (hcal_cellHot != 1) {
0487             }
0488           }
0489 
0490           if (debug_)
0491             std::cout << "ETA \t" << eta << "DEPTH \t" << dep << std::endl;
0492         }
0493       }
0494     }
0495   }
0496   close();
0497 }
0498 
0499 Bool_t HBHEMuonOfflineSimAnalyzer::Notify() {
0500   // The Notify() function is called when a new file is opened. This
0501   // can be either for a new TTree in a TChain or when when a new TTree
0502   // is started when using PROOF. It is normally not necessary to make changes
0503   // to the generated code, but the routine can be extended by the
0504   // user if needed. The return value is currently not used.
0505 
0506   return kTRUE;
0507 }
0508 
0509 void HBHEMuonOfflineSimAnalyzer::Show(Long64_t entry) {
0510   // Print contents of entry.
0511   // If entry is not specified, print current entry
0512   if (!fChain)
0513     return;
0514   fChain->Show(entry);
0515 }
0516 
0517 void HBHEMuonOfflineSimAnalyzer::BookHistograms(const char *fname) {
0518   output_file = TFile::Open(fname, "RECREATE");
0519   //output_file->cd();
0520   std::string type[] = {"tight", "soft", "loose"};
0521   char name[128], title[500];
0522 
0523   std::cout << "BookHistograms" << std::endl;
0524 
0525   nHist = 0;
0526   for (int eta = 0; eta < 29; ++eta) {
0527     int nDepth = NDepthBins(eta + 1, -1);
0528     int nPhi = NPhiBins(eta + 1);
0529     for (int depth = 0; depth < nDepth; depth++) {
0530       for (int PHI = 0; PHI < nPhi; ++PHI) {
0531         indxEta[eta][depth][PHI] = nHist;
0532         nHist += 2;
0533       }
0534     }
0535   }
0536 
0537   for (int i = 0; i < 3; ++i) {
0538     sprintf(name, "h_Pt_Muon_%s", type[i].c_str());
0539     sprintf(title, "p_{T} of %s muons (GeV)", type[i].c_str());
0540     h_Pt_Muon[i] = new TH1D(name, title, 100, 0, 200);
0541 
0542     sprintf(name, "h_Eta_Muon_%s", type[i].c_str());
0543     sprintf(title, "#eta of %s muons", type[i].c_str());
0544     h_Eta_Muon[i] = new TH1D(name, title, 50, -2.5, 2.5);
0545 
0546     sprintf(name, "h_Phi_Muon_%s", type[i].c_str());
0547     sprintf(title, "#phi of %s muons", type[i].c_str());
0548     h_Phi_Muon[i] = new TH1D(name, title, 100, -3.1415926, 3.1415926);
0549 
0550     sprintf(name, "h_P_Muon_%s", type[i].c_str());
0551     sprintf(title, "p of %s muons (GeV)", type[i].c_str());
0552     h_P_Muon[i] = new TH1D(name, title, 100, 0, 200);
0553 
0554     sprintf(name, "h_PF_Muon_%s", type[i].c_str());
0555     sprintf(title, "PF %s muons (GeV)", type[i].c_str());
0556     h_PF_Muon[i] = new TH1D(name, title, 2, 0, 2);
0557 
0558     sprintf(name, "h_Global_Muon_Chi2_%s", type[i].c_str());
0559     sprintf(title, "Chi2 Global %s muons (GeV)", type[i].c_str());
0560     h_GlobTrack_Chi[i] = new TH1D(name, title, 15, 0, 15);
0561 
0562     sprintf(name, "h_Global_Muon_Hits_%s", type[i].c_str());
0563     sprintf(title, "Global Hits %s muons (GeV)", type[i].c_str());
0564     h_Global_Muon_Hits[i] = new TH1D(name, title, 10, 0, 10);
0565 
0566     sprintf(name, "h_Matched_Stations_%s", type[i].c_str());
0567     sprintf(title, "Matched Stations %s muons (GeV)", type[i].c_str());
0568     h_MatchedStations[i] = new TH1D(name, title, 10, 0, 10);
0569 
0570     sprintf(name, "h_Transverse_ImpactParameter_%s", type[i].c_str());
0571     sprintf(title, "Transverse_ImpactParameter of %s muons (GeV)", type[i].c_str());
0572     h_Tight_TransImpactparameter[i] = new TH1D(name, title, 50, 0, 10);
0573 
0574     sprintf(name, "h_Longitudinal_ImpactParameter_%s", type[i].c_str());
0575     sprintf(title, "Longitudinal_ImpactParameter of %s muons (GeV)", type[i].c_str());
0576     h_Tight_LongitudinalImpactparameter[i] = new TH1D(name, title, 20, 0, 10);
0577 
0578     sprintf(name, "h_InnerTrack_PixelHits_%s", type[i].c_str());
0579     sprintf(title, "InnerTrack_PixelHits of %s muons (GeV)", type[i].c_str());
0580     h_InnerTrackPixelHits[i] = new TH1D(name, title, 20, 0, 20);
0581 
0582     sprintf(name, "h_TrackLayers_%s", type[i].c_str());
0583     sprintf(title, "No. of Tracker Layers of %s muons (GeV)", type[i].c_str());
0584     h_TrackerLayer[i] = new TH1D(name, title, 20, 0, 20);
0585     ;
0586 
0587     sprintf(name, "h_IsolationR04_%s", type[i].c_str());
0588     sprintf(title, "IsolationR04 %s muons (GeV)", type[i].c_str());
0589     h_IsolationR04[i] = new TH1D(name, title, 45, 0, 5);
0590     ;
0591 
0592     sprintf(name, "h_Global_Muon_%s", type[i].c_str());
0593     sprintf(title, "Global %s muons (GeV)", type[i].c_str());
0594     h_Global_Muon[i] = new TH1D(name, title, 2, 0, 2);
0595 
0596     sprintf(name, "h_TransImpactParameter_%s", type[i].c_str());
0597     sprintf(title, "TransImpactParameter of %s muons (GeV)", type[i].c_str());
0598     h_TransImpactParameter[i] = new TH1D(name, title, 100, 0, 0.5);
0599 
0600     sprintf(name, "h_TransImpactParameterBin1_%s", type[i].c_str());
0601     sprintf(title, "TransImpactParameter of %s muons (GeV) in -1.5 <= #phi <= 0.5", type[i].c_str());
0602     h_TransImpactParameterBin1[i] = new TH1D(name, title, 100, 0, 0.5);
0603 
0604     sprintf(name, "h_TransImpactParameterBin2_%s", type[i].c_str());
0605     sprintf(title, "TransImpactParameter of %s muons (GeV) in #phi> 0.5 and #phi< -1.5 ", type[i].c_str());
0606     h_TransImpactParameterBin2[i] = new TH1D(name, title, 100, 0, 0.5);
0607     //
0608     sprintf(name, "h_LongImpactParameter_%s", type[i].c_str());
0609     sprintf(title, "LongImpactParameter of %s muons (GeV)", type[i].c_str());
0610     h_LongImpactParameter[i] = new TH1D(name, title, 100, 0, 30);
0611 
0612     sprintf(name, "h_LongImpactParameterBin1_%s", type[i].c_str());
0613     sprintf(title, "LongImpactParameter of %s muons (GeV) in -1.5 <= #phi <= 0.5", type[i].c_str());
0614     h_LongImpactParameterBin1[i] = new TH1D(name, title, 100, 0, 30);
0615 
0616     sprintf(name, "h_LongImpactParameterBin2_%s", type[i].c_str());
0617     sprintf(title, "LongImpactParameter of %s muons (GeV) in #phi> 0.5 and #phi< -1.5 ", type[i].c_str());
0618     h_LongImpactParameterBin2[i] = new TH1D(name, title, 100, 0, 30);
0619 
0620     sprintf(name, "h_2D_Bin1_%s", type[i].c_str());
0621     sprintf(title, "Trans/Long ImpactParameter of %s muons (GeV) in -1.5 <= #phi< 0.5 ", type[i].c_str());
0622     h_2D_Bin1[i] = new TH2D(name, title, 100, 0, 0.5, 100, 0, 30);
0623 
0624     sprintf(name, "h_2D_Bin2_%s", type[i].c_str());
0625     sprintf(title, "Trans/Long ImpactParameter of %s muons (GeV) in #phi> 0.5 and #phi< -1.5 ", type[i].c_str());
0626     h_2D_Bin2[i] = new TH2D(name, title, 100, 0, 0.5, 100, 0, 30);
0627 
0628     sprintf(name, "h_ecal_energy_%s", type[i].c_str());
0629     sprintf(title, "ECAL energy for %s muons", type[i].c_str());
0630     h_ecal_energy[i] = new TH1D(name, title, 1000, -10.0, 90.0);
0631 
0632     sprintf(name, "h_hcal_energy_%s", type[i].c_str());
0633     sprintf(title, "HCAL energy for %s muons", type[i].c_str());
0634     h_hcal_energy[i] = new TH1D(name, title, 500, -10.0, 90.0);
0635 
0636     sprintf(name, "h_3x3_ecal_%s", type[i].c_str());
0637     sprintf(title, "ECAL energy in 3x3 for %s muons", type[i].c_str());
0638     h_3x3_ecal[i] = new TH1D(name, title, 1000, -10.0, 90.0);
0639 
0640     sprintf(name, "h_1x1_hcal_%s", type[i].c_str());
0641     sprintf(title, "HCAL energy in 1x1 for %s muons", type[i].c_str());
0642     h_1x1_hcal[i] = new TH1D(name, title, 500, -10.0, 90.0);
0643 
0644     sprintf(name, "h_EtaX_hcal_%s", type[i].c_str());
0645     sprintf(title, "HCAL energy as a function of i#eta for %s muons", type[i].c_str());
0646     h_EtaX_hcal[i] = new TProfile(name, title, 60, -30.0, 30.0);
0647 
0648     sprintf(name, "h_PhiY_hcal_%s", type[i].c_str());
0649     sprintf(title, "HCAL energy as a function of i#phi for %s muons", type[i].c_str());
0650     h_PhiY_hcal[i] = new TProfile(name, title, 72, 0, 72);
0651 
0652     sprintf(name, "h_EtaX_ecal_%s", type[i].c_str());
0653     sprintf(title, "EB energy as a function of i#eta for %s muons", type[i].c_str());
0654     h_EtaX_ecal[i] = new TProfile(name, title, 170, -85.0, 85.0);
0655 
0656     sprintf(name, "h_PhiY_ecal_%s", type[i].c_str());
0657     sprintf(title, "EB energy as a function of i#phi for %s muons", type[i].c_str());
0658     h_PhiY_ecal[i] = new TProfile(name, title, 360, 0, 360);
0659 
0660     sprintf(name, "h_Eta_ecal_%s", type[i].c_str());
0661     sprintf(title, "ECAL energy as a function of #eta for %s muons", type[i].c_str());
0662     h_Eta_ecal[i] = new TProfile(name, title, 100, -2.5, 2.5);
0663 
0664     sprintf(name, "h_Phi_ecal_%s", type[i].c_str());
0665     sprintf(title, "ECAL energy as a function of #phi for %s muons", type[i].c_str());
0666     h_Phi_ecal[i] = new TProfile(name, title, 100, -3.1415926, 3.1415926);
0667 
0668     sprintf(name, "h_MuonHittingEcal_%s", type[i].c_str());
0669     sprintf(title, "%s muons hitting ECAL", type[i].c_str());
0670     h_MuonHittingEcal[i] = new TH1D(name, title, 100, 0, 5.0);
0671 
0672     sprintf(name, "h_HotCell_%s", type[i].c_str());
0673     sprintf(title, "Hot cell for %s muons", type[i].c_str());
0674     h_HotCell[i] = new TH1D(name, title, 100, 0, 2);
0675 
0676     std::cout << "problem here" << std::endl;
0677     for (int eta = 0; eta < 29; ++eta) {
0678       int nDepth = NDepthBins(eta + 1, -1);
0679       int nPhi = NPhiBins(eta + 1);
0680       for (int depth = 0; depth < nDepth; ++depth) {
0681         for (int PHI = 0; PHI < nPhi; ++PHI) {
0682           int PHI0 = (nPhi == 72) ? PHI + 1 : 2 * PHI + 1;
0683           int ih = indxEta[eta][depth][PHI];
0684           std::cout << "eta:" << eta << " depth:" << depth << " PHI:" << PHI << ":" << PHI0 << " ih:" << ih
0685                     << std::endl;
0686 
0687           sprintf(name,
0688                   "h_Hot_MuonEnergy_hc_%d_%d_%d_%s_HotCell_ByActiveLength",
0689                   (eta + 1),
0690                   (depth + 1),
0691                   PHI0,
0692                   type[i].c_str());
0693           sprintf(title,
0694                   "HCAL energy in hot tower (i#eta=%d, depth=%d, i#phi = %d) for extrapolated %s muons (Hot Cell) "
0695                   "divided by Active Length",
0696                   (eta + 1),
0697                   (depth + 1),
0698                   PHI0,
0699                   type[i].c_str());
0700           h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[i][ih] = new TH1D(name, title, 4000, 0.0, 1.0);
0701           h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[i][ih]->Sumw2();
0702 
0703           ih++;
0704           sprintf(name,
0705                   "h_Hot_MuonEnergy_hc_%d_%d_%d_%s_HotCell_ByActiveLength",
0706                   -(eta + 1),
0707                   (depth + 1),
0708                   PHI0,
0709                   type[i].c_str());
0710           sprintf(title,
0711                   "HCAL energy in hot tower (i#eta=%d, depth=%d, i#phi=%d) for extrapolated %s muons (Hot Cell) "
0712                   "divided by Active Length",
0713                   -(eta + 1),
0714                   (depth + 1),
0715                   PHI0,
0716                   type[i].c_str());
0717           h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[i][ih] = new TH1D(name, title, 4000, 0.0, 1.0);
0718           h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[i][ih]->Sumw2();
0719         }
0720       }
0721       //output_file->cd();
0722     }
0723   }
0724   //output_file->cd();
0725 }
0726 
0727 bool HBHEMuonOfflineSimAnalyzer::LooseMuon() {
0728   if (pt_of_muon > 20.) {
0729     if (fabs(eta_of_muon) <= 5.0) {
0730       return true;
0731     }
0732   }
0733 
0734   return false;
0735 }
0736 
0737 bool HBHEMuonOfflineSimAnalyzer::SoftMuon() {
0738   if (pt_of_muon > 20.) {
0739     if (fabs(eta_of_muon) <= 5.0) {
0740       return true;
0741     }
0742   }
0743 
0744   return false;
0745 }
0746 
0747 bool HBHEMuonOfflineSimAnalyzer::tightMuon() {
0748   if (pt_of_muon > 20.) {
0749     if (fabs(eta_of_muon) <= 5.0) {
0750       return true;
0751     }
0752   }
0753 
0754   return false;
0755 }
0756 
0757 void HBHEMuonOfflineSimAnalyzer::etaPhiHcal(unsigned int detId, int &eta, int &phi, int &depth) {
0758   int zside, etaAbs;
0759   if ((detId & 0x1000000) == 0) {
0760     zside = (detId & 0x2000) ? (1) : (-1);
0761     etaAbs = (detId >> 7) & 0x3F;
0762     phi = detId & 0x7F;
0763     depth = (detId >> 14) & 0x1F;
0764   } else {
0765     zside = (detId & 0x80000) ? (1) : (-1);
0766     etaAbs = (detId >> 10) & 0x1FF;
0767     phi = detId & 0x3FF;
0768     depth = (detId >> 20) & 0xF;
0769   }
0770   eta = etaAbs * zside;
0771 }
0772 
0773 void HBHEMuonOfflineSimAnalyzer::etaPhiEcal(
0774     unsigned int detId, int &type, int &zside, int &etaX, int &phiY, int &plane, int &strip) {
0775   type = ((detId >> 25) & 0x7);
0776   plane = strip = 0;
0777   if (type == 1) {
0778     //Ecal Barrel
0779     zside = (detId & 0x10000) ? (1) : (-1);
0780     etaX = (detId >> 9) & 0x7F;
0781     phiY = detId & 0x1FF;
0782   } else if (type == 2) {
0783     zside = (detId & 0x4000) ? (1) : (-1);
0784     etaX = (detId >> 7) & 0x7F;
0785     phiY = (detId & 0x7F);
0786   } else if (type == 3) {
0787     zside = (detId & 0x80000) ? (1) : (-1);
0788     etaX = (detId >> 6) & 0x3F;
0789     /** get the sensor iy */
0790     phiY = (detId >> 12) & 0x3F;
0791     /** get the strip */
0792     plane = ((detId >> 18) & 0x1) + 1;
0793     strip = detId & 0x3F;
0794   } else {
0795     zside = etaX = phiY = 0;
0796   }
0797 }
0798 
0799 void HBHEMuonOfflineSimAnalyzer::calculateP(double pt, double eta, double &pM) {
0800   pM = (pt * cos(2 * (1 / atan(exp(eta)))));
0801 }
0802 
0803 void HBHEMuonOfflineSimAnalyzer::close() {
0804   output_file->cd();
0805   std::cout << "file yet to be Written" << std::endl;
0806   WriteHistograms();
0807   //    output_file->Write();
0808   std::cout << "file Written" << std::endl;
0809   output_file->Close();
0810   std::cout << "now doing return" << std::endl;
0811 }
0812 
0813 void HBHEMuonOfflineSimAnalyzer::WriteHistograms() {
0814   std::string type[] = {"tight", "soft", "loose"};
0815   char name[128];
0816 
0817   std::cout << "WriteHistograms" << std::endl;
0818   nHist = 0;
0819   for (int eta = 0; eta < 29; ++eta) {
0820     int nDepth = NDepthBins(eta + 1, -1);
0821     int nPhi = NPhiBins(eta + 1);
0822     if (debug_)
0823       std::cout << "Eta:" << eta << " nDepths " << nDepth << " nPhis " << nPhi << std::endl;
0824     for (int depth = 0; depth < nDepth; ++depth) {
0825       if (debug_)
0826         std::cout << "Eta:" << eta << "Depth:" << depth << std::endl;
0827       for (int PHI = 0; PHI < nPhi; ++PHI) {
0828         indxEta[eta][depth][PHI] = nHist;
0829         nHist += 2;
0830       }
0831     }
0832   }
0833 
0834   TDirectory *d_output_file[3][29];
0835   for (int i = 0; i < 3; ++i) {
0836     h_Pt_Muon[i]->Write();
0837     h_Eta_Muon[i]->Write();
0838     h_Phi_Muon[i]->Write();
0839     h_P_Muon[i]->Write();
0840     h_PF_Muon[i]->Write();
0841 
0842     h_GlobTrack_Chi[i]->Write();
0843     h_Global_Muon_Hits[i]->Write();
0844     h_MatchedStations[i]->Write();
0845 
0846     h_Tight_TransImpactparameter[i]->Write();
0847     h_Tight_LongitudinalImpactparameter[i]->Write();
0848 
0849     h_InnerTrackPixelHits[i]->Write();
0850     h_TrackerLayer[i]->Write();
0851     h_IsolationR04[i]->Write();
0852 
0853     h_Global_Muon[i]->Write();
0854     h_TransImpactParameter[i]->Write();
0855     ;
0856     h_TransImpactParameterBin1[i]->Write();
0857     h_TransImpactParameterBin2[i]->Write();
0858     //
0859     h_LongImpactParameter[i]->Write();
0860     h_LongImpactParameterBin1[i]->Write();
0861     h_LongImpactParameterBin2[i]->Write();
0862 
0863     h_ecal_energy[i]->Write();
0864     h_hcal_energy[i]->Write();
0865     ;
0866     h_3x3_ecal[i]->Write();
0867     h_1x1_hcal[i]->Write();
0868     ;
0869 
0870     h_EtaX_hcal[i]->Write();
0871     h_PhiY_hcal[i]->Write();
0872     ;
0873 
0874     h_EtaX_ecal[i]->Write();
0875     ;
0876     h_PhiY_ecal[i]->Write();
0877     ;
0878     h_Eta_ecal[i]->Write();
0879     ;
0880     h_Phi_ecal[i]->Write();
0881     ;
0882     h_MuonHittingEcal[i]->Write();
0883     ;
0884     h_HotCell[i]->Write();
0885     ;
0886 
0887     output_file->cd();
0888     for (int eta = 0; eta < 29; ++eta) {
0889       int nDepth = NDepthBins(eta + 1, -1);
0890       int nPhi = NPhiBins(eta + 1);
0891       sprintf(name, "Dir_muon_type_%s_ieta%d", type[i].c_str(), eta + 1);
0892       d_output_file[i][eta] = output_file->mkdir(name);
0893       d_output_file[i][eta]->cd();
0894       for (int depth = 0; depth < nDepth; ++depth) {
0895         for (int PHI = 0; PHI < nPhi; ++PHI) {
0896           int ih = indxEta[eta][depth][PHI];
0897           h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[i][ih]->Write();
0898           ih++;
0899           h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[i][ih]->Write();
0900         }
0901       }
0902       output_file->cd();
0903     }
0904   }
0905   output_file->cd();
0906 }
0907 
0908 int HBHEMuonOfflineSimAnalyzer::NDepthBins(int eta, int phi) {
0909   // Run 1 scenario
0910   int nDepthR1[29] = {1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 3, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 3, 3, 2};
0911   // Run 2 scenario from 2018
0912   int nDepthR2[29] = {1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 5, 6, 6, 6, 6, 6, 6, 6, 7, 7, 7, 3};
0913   // Run 3 scenario
0914   int nDepthR3[29] = {4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 3, 5, 6, 6, 6, 6, 6, 6, 6, 7, 7, 7, 3};
0915   // Run 4 scenario
0916   int nDepthR4[29] = {4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7};
0917   // for a test scenario with multi depth segmentation considered during Run 1
0918   //    int  nDepth[29]={3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,5,5,5,5,5,5,5,5,5,5,5,5,5};
0919   // modeLHC_ = 0 --> nbin defined maxDepthHB/HE
0920   //          = 1 -->      corresponds to Run 1 (valid till 2016)
0921   //          = 2 -->      corresponds to Run 2 (2018 geometry)
0922   //          = 3 -->      corresponds to Run 3 (post LS2)
0923   //          = 4 -->      corresponds to 2017 (Plan 1)
0924   //          = 5 -->      corresponds to Run 4 (post LS3)
0925   int nbin(0);
0926   if (modeLHC_ == 0) {
0927     if (eta <= 15) {
0928       nbin = maxDepthHB_;
0929     } else if (eta == 16) {
0930       nbin = 4;
0931     } else {
0932       nbin = maxDepthHE_;
0933     }
0934   } else if (modeLHC_ == 1) {
0935     nbin = nDepthR1[eta - 1];
0936   } else if (modeLHC_ == 2) {
0937     nbin = nDepthR2[eta - 1];
0938   } else if (modeLHC_ == 3) {
0939     nbin = nDepthR3[eta - 1];
0940   } else if (modeLHC_ == 4) {
0941     if (phi > 0) {
0942       if (eta >= 16 && phi >= 63 && phi <= 66) {
0943         nbin = nDepthR2[eta - 1];
0944       } else {
0945         nbin = nDepthR1[eta - 1];
0946       }
0947     } else {
0948       if (eta >= 16) {
0949         nbin = (nDepthR2[eta - 1] > nDepthR1[eta - 1]) ? nDepthR2[eta - 1] : nDepthR1[eta - 1];
0950       } else {
0951         nbin = nDepthR1[eta - 1];
0952       }
0953     }
0954   } else {
0955     if (eta > 0 && eta < 30) {
0956       nbin = nDepthR4[eta - 1];
0957     } else {
0958       nbin = nDepthR4[28];
0959     }
0960   }
0961   return nbin;
0962 }
0963 
0964 int HBHEMuonOfflineSimAnalyzer::NPhiBins(int eta) {
0965   int nphi = (eta <= 20) ? 72 : 36;
0966   return nphi;
0967 }