Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-08-26 03:15:32

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