Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 14:14:02

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