Back to home page

Project CMSSW displayed by LXR

 
 

    


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

0001 ///////////////////////////////////////////////////////////////////////////////
0002 //
0003 //  HBHEMuonOfflineAnalyzer h1(tree,outfile,rcorFile,flag,mode,maxDHB,maxDHE,
0004 //                             runLo,runHi,etaMin,etaMax,debug);
0005 //  HBHEMuonOfflineAnalyzer h1(infile,outfile,rcorFile,flag,mode,maxDHB,maxDHE,
0006 //                             runLo,runHi,etaMin,etaMax,debug);
0007 //   h1.Loop()
0008 //
0009 //   tree       TTree*       Pointer to the tree chain
0010 //   infile     const char*  Name of the input file
0011 //   outfile    const char*  Name of the output file
0012 //                           (dyll_PU20_25_output_10.root)
0013 //   rcorFile   consr char*  name of the text file having the correction factors
0014 //                           as a function of run numbers to be used for raddam
0015 //                           correction (default="", no corr.)
0016 //   flag       int          Flag of 2 digits ("to"): to where "o" decides if
0017 //                           corrected (1) or default (0) energy to be used;
0018 //                           "t" decides if all depths to be merged (1) or not
0019 //                           (0) (default is 0)
0020 //   mode       int          Geometry file used 0:(defined by maxDHB/HE);
0021 //                           1 (Run 1; valid till 2016); 2 (Run 2; 2018);
0022 //                           3 (Run 3; post LS2); 4 (2017 Plan 1);
0023 //                           5 (Run 4; post LS3); default (2)
0024 //   maxDHB     int          Maximum number of depths for HB (4)
0025 //   maxDHE     int          Maximum number of depths for HE (7)
0026 //   runLO      int          Minimum run number (1)
0027 //   runHI      int          Maximum run number (99999999)
0028 //   etaMin     int          Minimum (absolute) eta value (1)
0029 //   etaMax     int          Maximum (absolute) eta value (29)
0030 //
0031 ///////////////////////////////////////////////////////////////////////////////
0032 
0033 #include <algorithm>
0034 #include <cmath>
0035 #include <fstream>
0036 #include <iostream>
0037 #include <sstream>
0038 #include <vector>
0039 #include <string>
0040 #include <TCanvas.h>
0041 #include <TChain.h>
0042 #include <TFile.h>
0043 #include <TH1.h>
0044 #include <TH1D.h>
0045 #include <TH2.h>
0046 #include <TH2D.h>
0047 #include <TProfile.h>
0048 #include <TROOT.h>
0049 #include <TStyle.h>
0050 #include <TString.h>
0051 #include <TTree.h>
0052 
0053 class HBHEMuonOfflineAnalyzer {
0054 public:
0055   TChain *fChain;  //!pointer to the analyzed TTree/TChain
0056   Int_t fCurrent;  //!current Tree number in a TChain
0057 
0058   // Fixed size dimensions of array or collections stored in the TTree if any.
0059   // Declaration of leaf types
0060 
0061   UInt_t Event_No;
0062   UInt_t Run_No;
0063   UInt_t LumiNumber;
0064   UInt_t BXNumber;
0065   UInt_t GoodVertex;
0066   bool PF_Muon;
0067   bool Global_Muon;
0068   bool Tracker_muon;
0069   bool MuonIsTight;
0070   bool MuonIsMedium;
0071   double pt_of_muon;
0072   double eta_of_muon;
0073   double phi_of_muon;
0074   double energy_of_muon;
0075   double p_of_muon;
0076   float muon_trkKink;
0077   float muon_chi2LocalPosition;
0078   float muon_segComp;
0079   int TrackerLayer;
0080   int NumPixelLayers;
0081   int InnerTrackPixelHits;
0082   bool innerTrack;
0083   double chiTracker;
0084   double DxyTracker;
0085   double DzTracker;
0086   double innerTrackpt;
0087   double innerTracketa;
0088   double innerTrackphi;
0089   double tight_validFraction;
0090   bool OuterTrack;
0091   double OuterTrackPt;
0092   double OuterTrackEta;
0093   double OuterTrackPhi;
0094   double OuterTrackChi;
0095   int OuterTrackHits;
0096   int OuterTrackRHits;
0097   bool GlobalTrack;
0098   double GlobalTrckPt;
0099   double GlobalTrckEta;
0100   double GlobalTrckPhi;
0101   int Global_Muon_Hits;
0102   int MatchedStations;
0103   double GlobTrack_Chi;
0104   double Tight_LongitudinalImpactparameter;
0105   double Tight_TransImpactparameter;
0106   double IsolationR04;
0107   double IsolationR03;
0108   double ecal_3into3;
0109   double hcal_3into3;
0110   double tracker_3into3;
0111   bool matchedId;
0112   bool hcal_cellHot;
0113   double ecal_3x3;
0114   double hcal_1x1;
0115   unsigned int ecal_detID;
0116   unsigned int hcal_detID;
0117   unsigned int ehcal_detID;
0118   int hcal_ieta;
0119   int hcal_iphi;
0120   double hcal_edepth1;
0121   double hcal_activeL1;
0122   double hcal_edepthHot1;
0123   double hcal_activeHotL1;
0124   double hcal_edepthCorrect1;
0125   double hcal_edepthHotCorrect1;
0126   double hcal_cdepthHot1;
0127   double hcal_cdepthHotBG1;
0128   bool hcal_depthMatch1;
0129   bool hcal_depthMatchHot1;
0130   double hcal_edepth2;
0131   double hcal_activeL2;
0132   double hcal_edepthHot2;
0133   double hcal_activeHotL2;
0134   double hcal_edepthCorrect2;
0135   double hcal_edepthHotCorrect2;
0136   double hcal_cdepthHot2;
0137   double hcal_cdepthHotBG2;
0138   bool hcal_depthMatch2;
0139   bool hcal_depthMatchHot2;
0140   double hcal_edepth3;
0141   double hcal_activeL3;
0142   double hcal_edepthHot3;
0143   double hcal_activeHotL3;
0144   double hcal_edepthCorrect3;
0145   double hcal_edepthHotCorrect3;
0146   double hcal_cdepthHot3;
0147   double hcal_cdepthHotBG3;
0148   bool hcal_depthMatch3;
0149   bool hcal_depthMatchHot3;
0150   double hcal_edepth4;
0151   double hcal_activeL4;
0152   double hcal_edepthHot4;
0153   double hcal_activeHotL4;
0154   double hcal_edepthCorrect4;
0155   double hcal_edepthHotCorrect4;
0156   double hcal_cdepthHot4;
0157   double hcal_cdepthHotBG4;
0158   bool hcal_depthMatch4;
0159   bool hcal_depthMatchHot4;
0160   double hcal_edepth5;
0161   double hcal_activeL5;
0162   double hcal_edepthHot5;
0163   double hcal_activeHotL5;
0164   double hcal_edepthCorrect5;
0165   double hcal_edepthHotCorrect5;
0166   double hcal_cdepthHot5;
0167   double hcal_cdepthHotBG5;
0168   bool hcal_depthMatch5;
0169   bool hcal_depthMatchHot5;
0170   double hcal_edepth6;
0171   double hcal_activeL6;
0172   double hcal_edepthHot6;
0173   double hcal_activeHotL6;
0174   double hcal_edepthCorrect6;
0175   double hcal_edepthHotCorrect6;
0176   double hcal_cdepthHot6;
0177   double hcal_cdepthHotBG6;
0178   bool hcal_depthMatch6;
0179   bool hcal_depthMatchHot6;
0180   double hcal_edepth7;
0181   double hcal_activeL7;
0182   double hcal_edepthHot7;
0183   double hcal_activeHotL7;
0184   double hcal_edepthCorrect7;
0185   double hcal_edepthHotCorrect7;
0186   double hcal_cdepthHot7;
0187   double hcal_cdepthHotBG7;
0188   bool hcal_depthMatch7;
0189   bool hcal_depthMatchHot7;
0190   double activeLength;
0191   double activeLengthHot;
0192   std::vector<int> *hltresults;
0193   std::vector<std::string> *all_triggers;
0194 
0195   // List of branches
0196   TBranch *b_Event_No;                           //!
0197   TBranch *b_Run_No;                             //!
0198   TBranch *b_LumiNumber;                         //!
0199   TBranch *b_BXNumber;                           //!
0200   TBranch *b_GoodVertex;                         //!
0201   TBranch *b_PF_Muon;                            //!
0202   TBranch *b_Global_Muon;                        //!
0203   TBranch *b_Tracker_muon;                       //!
0204   TBranch *b_MuonIsTight;                        //!
0205   TBranch *b_MuonIsMedium;                       //!
0206   TBranch *b_pt_of_muon;                         //!
0207   TBranch *b_eta_of_muon;                        //!
0208   TBranch *b_phi_of_muon;                        //!
0209   TBranch *b_energy_of_muon;                     //!
0210   TBranch *b_p_of_muon;                          //!
0211   TBranch *b_muon_trkKink;                       //!
0212   TBranch *b_muon_chi2LocalPosition;             //!
0213   TBranch *b_muon_segComp;                       //!
0214   TBranch *b_TrackerLayer;                       //!
0215   TBranch *b_NumPixelLayers;                     //!
0216   TBranch *b_InnerTrackPixelHits;                //!
0217   TBranch *b_innerTrack;                         //!
0218   TBranch *b_chiTracker;                         //!
0219   TBranch *b_DxyTracker;                         //!
0220   TBranch *b_DzTracker;                          //!
0221   TBranch *b_innerTrackpt;                       //!
0222   TBranch *b_innerTracketa;                      //!
0223   TBranch *b_innerTrackphi;                      //!
0224   TBranch *b_tight_validFraction;                //!
0225   TBranch *b_OuterTrack;                         //!
0226   TBranch *b_OuterTrackPt;                       //!
0227   TBranch *b_OuterTrackEta;                      //!
0228   TBranch *b_OuterTrackPhi;                      //!
0229   TBranch *b_OuterTrackChi;                      //!
0230   TBranch *b_OuterTrackHits;                     //!
0231   TBranch *b_OuterTrackRHits;                    //!
0232   TBranch *b_GlobalTrack;                        //!
0233   TBranch *b_GlobalTrckPt;                       //!
0234   TBranch *b_GlobalTrckEta;                      //!
0235   TBranch *b_GlobalTrckPhi;                      //!
0236   TBranch *b_Global_Muon_Hits;                   //!
0237   TBranch *b_MatchedStations;                    //!
0238   TBranch *b_GlobTrack_Chi;                      //!
0239   TBranch *b_Tight_LongitudinalImpactparameter;  //!
0240   TBranch *b_Tight_TransImpactparameter;         //!
0241   TBranch *b_IsolationR04;                       //!
0242   TBranch *b_IsolationR03;                       //!
0243   TBranch *b_ecal_3into3;                        //!
0244   TBranch *b_hcal_3into3;                        //!
0245   TBranch *b_tracker_3into3;                     //!
0246   TBranch *b_matchedId;                          //!
0247   TBranch *b_hcal_cellHot;                       //!
0248   TBranch *b_ecal_3x3;                           //!
0249   TBranch *b_hcal_1x1;                           //!
0250   TBranch *b_hcal_detID;                         //!
0251   TBranch *b_ecal_detID;                         //!
0252   TBranch *b_ehcal_detID;                        //!
0253   TBranch *b_hcal_ieta;                          //!
0254   TBranch *b_hcal_iphi;                          //!
0255   TBranch *b_hcal_edepth1;                       //!
0256   TBranch *b_hcal_activeL1;                      //!
0257   TBranch *b_hcal_edepthHot1;                    //!
0258   TBranch *b_hcal_activeHotL1;                   //!
0259   TBranch *b_hcal_edepthCorrect1;                //!
0260   TBranch *b_hcal_edepthHotCorrect1;             //!
0261   TBranch *b_hcal_cdepthHot1;                    //!
0262   TBranch *b_hcal_cdepthHotBG1;                  //!
0263   TBranch *b_hcal_depthMatch1;                   //!
0264   TBranch *b_hcal_depthMatchHot1;                //!
0265   TBranch *b_hcal_edepth2;                       //!
0266   TBranch *b_hcal_activeL2;                      //!
0267   TBranch *b_hcal_edepthHot2;                    //!
0268   TBranch *b_hcal_activeHotL2;                   //!
0269   TBranch *b_hcal_edepthCorrect2;                //!
0270   TBranch *b_hcal_edepthHotCorrect2;             //!
0271   TBranch *b_hcal_cdepthHot2;                    //!
0272   TBranch *b_hcal_cdepthHotBG2;                  //!
0273   TBranch *b_hcal_depthMatch2;                   //!
0274   TBranch *b_hcal_depthMatchHot2;                //!
0275   TBranch *b_hcal_edepth3;                       //!
0276   TBranch *b_hcal_activeL3;                      //!
0277   TBranch *b_hcal_edepthHot3;                    //!
0278   TBranch *b_hcal_activeHotL3;                   //!
0279   TBranch *b_hcal_edepthCorrect3;                //!
0280   TBranch *b_hcal_edepthHotCorrect3;             //!
0281   TBranch *b_hcal_cdepthHot3;                    //!
0282   TBranch *b_hcal_cdepthHotBG3;                  //!
0283   TBranch *b_hcal_depthMatch3;                   //!
0284   TBranch *b_hcal_depthMatchHot3;                //!
0285   TBranch *b_hcal_edepth4;                       //!
0286   TBranch *b_hcal_activeL4;                      //!
0287   TBranch *b_hcal_edepthHot4;                    //!
0288   TBranch *b_hcal_activeHotL4;                   //!
0289   TBranch *b_hcal_edepthCorrect4;                //!
0290   TBranch *b_hcal_edepthHotCorrect4;             //!
0291   TBranch *b_hcal_cdepthHot4;                    //!
0292   TBranch *b_hcal_cdepthHotBG4;                  //!
0293   TBranch *b_hcal_depthMatch4;                   //!
0294   TBranch *b_hcal_depthMatchHot4;                //!
0295   TBranch *b_hcal_edepth5;                       //!
0296   TBranch *b_hcal_activeL5;                      //!
0297   TBranch *b_hcal_edepthHot5;                    //!
0298   TBranch *b_hcal_activeHotL5;                   //!
0299   TBranch *b_hcal_edepthCorrect5;                //!
0300   TBranch *b_hcal_edepthHotCorrect5;             //!
0301   TBranch *b_hcal_cdepthHot5;                    //!
0302   TBranch *b_hcal_cdepthHotBG5;                  //!
0303   TBranch *b_hcal_depthMatch5;                   //!
0304   TBranch *b_hcal_depthMatchHot5;                //!
0305   TBranch *b_hcal_edepth6;                       //!
0306   TBranch *b_hcal_activeL6;                      //!
0307   TBranch *b_hcal_edepthHot6;                    //!
0308   TBranch *b_hcal_activeHotL6;                   //!
0309   TBranch *b_hcal_edepthCorrect6;                //!
0310   TBranch *b_hcal_edepthHotCorrect6;             //!
0311   TBranch *b_hcal_cdepthHot6;                    //!
0312   TBranch *b_hcal_cdepthHotBG6;                  //!
0313   TBranch *b_hcal_depthMatch6;                   //!
0314   TBranch *b_hcal_depthMatchHot6;                //!
0315   TBranch *b_hcal_edepth7;                       //!
0316   TBranch *b_hcal_activeL7;                      //!
0317   TBranch *b_hcal_edepthHot7;                    //!
0318   TBranch *b_hcal_activeHotL7;                   //!
0319   TBranch *b_hcal_edepthCorrect7;                //!
0320   TBranch *b_hcal_edepthHotCorrect7;             //!
0321   TBranch *b_hcal_cdepthHot7;                    //!
0322   TBranch *b_hcal_cdepthHotBG7;                  //!
0323   TBranch *b_hcal_depthMatch7;                   //!
0324   TBranch *b_hcal_depthMatchHot7;                //!
0325   TBranch *b_activeLength;                       //!
0326   TBranch *b_activeLengthHot;                    //!
0327   TBranch *b_hltresults;                         //!
0328   TBranch *b_all_triggers;                       //!
0329 
0330   HBHEMuonOfflineAnalyzer(TChain *tree = 0,
0331                           const char *outfile = "dyll_PU20_25_output_10.root",
0332                           const char *rcorFileName = "",
0333                           int flag = 0,
0334                           int mode = 2,
0335                           int maxDHB = 4,
0336                           int maxDHE = 7,
0337                           int runLo = 1,
0338                           int runHi = 99999999,
0339                           int etaMin = 1,
0340                           int etaMax = 29,
0341                           bool debug = false);
0342   HBHEMuonOfflineAnalyzer(const char *infile,
0343                           const char *outfile = "dyll_PU20_25_output_10.root",
0344                           const char *rcorFileName = "",
0345                           int flag = 0,
0346                           int mode = 2,
0347                           int maxDHB = 4,
0348                           int maxDHE = 7,
0349                           int runLo = 1,
0350                           int runHi = 99999999,
0351                           int etaMin = 1,
0352                           int etaMax = 29,
0353                           bool debug = false);
0354   // mode of LHC is kept 1 for 2017 scenario as no change in depth segmentation
0355   // mode of LHC is 0 for 2021
0356   virtual ~HBHEMuonOfflineAnalyzer();
0357 
0358   virtual Int_t Cut(Long64_t entry);
0359   virtual Int_t GetEntry(Long64_t entry);
0360   virtual Long64_t LoadTree(Long64_t entry);
0361   virtual void Init(TChain *tree,
0362                     const char *rcorFileName,
0363                     int flag,
0364                     int mode,
0365                     int maxDHB,
0366                     int maxDHE,
0367                     int runLo,
0368                     int runHi,
0369                     int etaMin,
0370                     int etaMax);
0371   virtual void Loop();
0372   virtual Bool_t Notify();
0373   virtual void Show(Long64_t entry = -1);
0374 
0375   bool fillChain(TChain *chain, const char *inputFileList);
0376   bool readCorr(const char *rcorFileName);
0377   void bookHistograms(const char *);
0378   bool getEnergy(int dep, double &enb, double &enu, double &enh, double &enc, double &chgS, double &chgB, double &actL);
0379   void writeHistograms();
0380   bool looseMuon();
0381   bool tightMuon();
0382   bool softMuon();
0383   bool mediumMuon2016();
0384   void etaPhiHcal(unsigned int detId, int &eta, int &phi, int &depth);
0385   void etaPhiEcal(unsigned int detId, int &type, int &zside, int &etaX, int &phiY, int &plane, int &strip);
0386   void calculateP(double pt, double eta, double &pM);
0387   void close();
0388   int nDepthBins(int ieta, int iphi);
0389   int nPhiBins(int ieta);
0390   float getCorr(int run, unsigned int id);
0391   std::vector<std::string> splitString(const std::string &);
0392   unsigned int getDetIdHBHE(int ieta, int iphi, int depth);
0393   unsigned int getDetId(int subdet, int ieta, int iphi, int depth);
0394   unsigned int correctDetId(const unsigned int &detId);
0395   void unpackDetId(unsigned int detId, int &subdet, int &zside, int &ieta, int &iphi, int &depth);
0396 
0397 private:
0398   static const int maxDep_ = 7;
0399   static const int maxEta_ = 29;
0400   static const int maxPhi_ = 72;
0401   //3x16x72x2 + 5x4x72x2 + 5x9x36x2
0402   static const int maxHist_ = 20000;  //13032;
0403   static const int nCut_ = 1;
0404   static const unsigned int nmax_ = 10;
0405   const bool debug_;
0406   int modeLHC_, maxDepthHB_, maxDepthHE_, maxDepth_;
0407   int runLo_, runHi_, etaMin_, etaMax_;
0408   bool cFactor_, useCorrect_, mergeDepth_;
0409   int nHist, nDepths[maxEta_], nDepthsPhi[maxEta_];
0410   int indxEta[maxEta_][maxDep_][maxPhi_];
0411   TFile *output_file;
0412   std::map<unsigned int, float> corrFac_[nmax_];
0413   std::vector<int> runlow_;
0414 
0415   TTree *outtree_;
0416   int t_ieta, t_iphi, t_nvtx;
0417   double t_p, t_ediff;
0418   std::vector<double> t_ene, t_enec, t_actln, t_charge;
0419   std::vector<int> t_depth;
0420 
0421   TH1D *h_evtype, *h_Pt_Muon[3], *h_Eta_Muon[3], *h_Phi_Muon[3], *h_P_Muon[3];
0422   TH1D *h_PF_Muon[3], *h_GlobTrack_Chi[3], *h_Global_Muon_Hits[3];
0423   TH1D *h_MatchedStations[3], *h_Tight_TransImpactparameter[3];
0424   TH1D *h_Tight_LongitudinalImpactparameter[3], *h_InnerTrackPixelHits[3];
0425   TH1D *h_TrackerLayer[3], *h_IsolationR04[3], *h_Global_Muon[3];
0426   TH1D *h_LongImpactParameter[3], *h_LongImpactParameterBin1[3], *h_LongImpactParameterBin2[3];
0427 
0428   TH1D *h_TransImpactParameter[3], *h_TransImpactParameterBin1[3], *h_TransImpactParameterBin2[3];
0429   TH1D *h_Hot_MuonEnergy_hcal_ClosestCell[3][maxHist_], *h_Hot_MuonEnergy_hcal_HotCell[3][maxHist_],
0430       *h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[3][maxHist_], *h_HotCell_MuonEnergy_phi[3][maxHist_],
0431       *h_active_length_Fill[3][maxHist_], *h_p_muon_ineta[3][maxHist_], *h_charge_signal[3][maxHist_],
0432       *h_charge_bg[3][maxHist_];
0433   TH2D *h_2D_Bin1[3], *h_2D_Bin2[3];
0434   TH1D *h_ecal_energy[3], *h_hcal_energy[3], *h_3x3_ecal[3], *h_1x1_hcal[3];
0435   TH1D *h_MuonHittingEcal[3], *h_HotCell[3], *h_MuonEnergy_hcal[3][maxHist_];
0436   TH1D *h_Hot_MuonEnergy_hcal[3][maxHist_];
0437   TH2D *hcal_ietaVsEnergy[3];
0438   TProfile *h_EtaX_hcal[3], *h_PhiY_hcal[3], *h_EtaX_ecal[3], *h_PhiY_ecal[3];
0439   TProfile *h_Eta_ecal[3], *h_Phi_ecal[3];
0440   TProfile *h_MuonEnergy_eta[3][maxDep_], *h_MuonEnergy_phi[3][maxDep_], *h_MuonEnergy_muon_eta[3][maxDep_];
0441   TProfile *h_Hot_MuonEnergy_eta[3][maxDep_], *h_Hot_MuonEnergy_phi[3][maxDep_], *h_Hot_MuonEnergy_muon_eta[3][maxDep_];
0442   TProfile *h_IsoHot_MuonEnergy_eta[3][maxDep_], *h_IsoHot_MuonEnergy_phi[3][maxDep_],
0443       *h_IsoHot_MuonEnergy_muon_eta[3][maxDep_];
0444   TProfile *h_IsoWithoutHot_MuonEnergy_eta[3][maxDep_], *h_IsoWithoutHot_MuonEnergy_phi[3][maxDep_],
0445       *h_IsoWithoutHot_MuonEnergy_muon_eta[3][maxDep_];
0446   TProfile *h_HotWithoutIso_MuonEnergy_eta[3][maxDep_], *h_HotWithoutIso_MuonEnergy_phi[3][maxDep_],
0447       *h_HotWithoutIso_MuonEnergy_muon_eta[3][maxDep_];
0448 };
0449 
0450 HBHEMuonOfflineAnalyzer::HBHEMuonOfflineAnalyzer(TChain *tree,
0451                                                  const char *outFileName,
0452                                                  const char *rcorFileName,
0453                                                  int flag,
0454                                                  int mode,
0455                                                  int maxDHB,
0456                                                  int maxDHE,
0457                                                  int runLo,
0458                                                  int runHi,
0459                                                  int etaMin,
0460                                                  int etaMax,
0461                                                  bool deb)
0462     : debug_(deb), cFactor_(false) {
0463   Init(tree, rcorFileName, flag, mode, maxDHB, maxDHE, runLo, runHi, etaMin, etaMax);
0464 
0465   //Now book histograms
0466   bookHistograms(outFileName);
0467 }
0468 
0469 HBHEMuonOfflineAnalyzer::HBHEMuonOfflineAnalyzer(const char *infile,
0470                                                  const char *outFileName,
0471                                                  const char *rcorFileName,
0472                                                  int flag,
0473                                                  int mode,
0474                                                  int maxDHB,
0475                                                  int maxDHE,
0476                                                  int runLo,
0477                                                  int runHi,
0478                                                  int etaMin,
0479                                                  int etaMax,
0480                                                  bool deb)
0481     : debug_(deb), cFactor_(false) {
0482   TChain *chain = new TChain("hcalHBHEMuon/TREE");
0483   if (!fillChain(chain, infile)) {
0484     std::cout << "*****No valid tree chain can be obtained*****" << std::endl;
0485   } else {
0486     std::cout << "Proceed with a tree chain with " << chain->GetEntries() << " entries" << std::endl;
0487     Init(chain, rcorFileName, flag, mode, maxDHB, maxDHE, runLo, runHi, etaMin, etaMax);
0488 
0489     //Now book histograms
0490     bookHistograms(outFileName);
0491   }
0492 }
0493 
0494 HBHEMuonOfflineAnalyzer::~HBHEMuonOfflineAnalyzer() {
0495   if (!fChain)
0496     return;
0497   delete fChain->GetCurrentFile();
0498 }
0499 
0500 Int_t HBHEMuonOfflineAnalyzer::Cut(Long64_t) {
0501   // This function may be called from Loop.
0502   // returns  1 if entry is accepted.
0503   // returns -1 otherwise.
0504   return 1;
0505 }
0506 
0507 Int_t HBHEMuonOfflineAnalyzer::GetEntry(Long64_t entry) {
0508   // Read contents of entry.
0509   if (!fChain)
0510     return 0;
0511   return fChain->GetEntry(entry);
0512 }
0513 
0514 Long64_t HBHEMuonOfflineAnalyzer::LoadTree(Long64_t entry) {
0515   // Set the environment to read one entry
0516   if (!fChain)
0517     return -5;
0518   Long64_t centry = fChain->LoadTree(entry);
0519   if (centry < 0)
0520     return centry;
0521   if (fChain->GetTreeNumber() != fCurrent) {
0522     fCurrent = fChain->GetTreeNumber();
0523     Notify();
0524   }
0525   return centry;
0526 }
0527 
0528 void HBHEMuonOfflineAnalyzer::Init(TChain *tree,
0529                                    const char *rcorFileName,
0530                                    int flag,
0531                                    int mode,
0532                                    int maxDHB,
0533                                    int maxDHE,
0534                                    int runLo,
0535                                    int runHi,
0536                                    int etaMin,
0537                                    int etaMax) {
0538   modeLHC_ = mode;
0539   maxDepthHB_ = maxDHB;
0540   maxDepthHE_ = maxDHE;
0541   maxDepth_ = (maxDepthHB_ > maxDepthHE_) ? maxDepthHB_ : maxDepthHE_;
0542   runLo_ = runLo;
0543   runHi_ = runHi;
0544   etaMin_ = (etaMin > 0) ? etaMin : 1;
0545   etaMax_ = (etaMax <= 29) ? etaMax : 29;
0546   if (etaMax_ <= etaMin_) {
0547     if (etaMax_ == 29)
0548       etaMin_ = etaMax_ - 1;
0549     else
0550       etaMax_ = etaMin_ + 1;
0551   }
0552   useCorrect_ = ((flag % 10) > 0);
0553   mergeDepth_ = (((flag / 10) % 10) > 0);
0554   if (std::string(rcorFileName) != "")
0555     cFactor_ = readCorr(rcorFileName);
0556 
0557   // The Init() function is called when the selector needs to initialize
0558   // a new tree or chain. Typically here the branch addresses and branch
0559   // pointers of the tree will be set.
0560   // It is normally not necessary to make changes to the generated
0561   // code, but the routine can be extended by the user if needed.
0562   // Init() will be called many times when running on PROOF
0563   // (once per file to be processed).
0564 
0565   // Set object pointer
0566 
0567   PF_Muon = 0;
0568   Global_Muon = 0;
0569   Tracker_muon = 0;
0570   pt_of_muon = 0;
0571   eta_of_muon = 0;
0572   phi_of_muon = 0;
0573   energy_of_muon = 0;
0574   p_of_muon = 0;
0575   muon_trkKink = 0;
0576   muon_chi2LocalPosition = 0;
0577   muon_segComp = 0;
0578   TrackerLayer = 0;
0579   NumPixelLayers = 0;
0580   InnerTrackPixelHits = 0;
0581   innerTrack = 0;
0582   chiTracker = 0;
0583   DxyTracker = 0;
0584   DzTracker = 0;
0585   innerTrackpt = 0;
0586   innerTracketa = 0;
0587   innerTrackphi = 0;
0588   tight_validFraction = 0;
0589   OuterTrack = 0;
0590   OuterTrackPt = 0;
0591   OuterTrackEta = 0;
0592   OuterTrackPhi = 0;
0593   OuterTrackHits = 0;
0594   OuterTrackRHits = 0;
0595   OuterTrackChi = 0;
0596   GlobalTrack = 0;
0597   GlobalTrckPt = 0;
0598   GlobalTrckEta = 0;
0599   GlobalTrckPhi = 0;
0600   Global_Muon_Hits = 0;
0601   MatchedStations = 0;
0602   GlobTrack_Chi = 0;
0603   Tight_LongitudinalImpactparameter = 0;
0604   Tight_TransImpactparameter = 0;
0605   IsolationR04 = 0;
0606   IsolationR03 = 0;
0607   ecal_3into3 = 0;
0608   hcal_3into3 = 0;
0609   tracker_3into3 = 0;
0610   matchedId = 0;
0611   hcal_cellHot = 0;
0612   ecal_3x3 = 0;
0613   hcal_1x1 = 0;
0614   ecal_detID = 0;
0615   hcal_detID = 0;
0616   ehcal_detID = 0;
0617   hcal_edepth1 = 0;
0618   hcal_activeL1 = 0;
0619   hcal_edepthHot1 = 0;
0620   hcal_activeHotL1 = 0;
0621   hcal_edepthCorrect1 = 0;
0622   hcal_edepthHotCorrect1 = 0;
0623   hcal_cdepthHot1 = 0;
0624   hcal_cdepthHotBG1 = 0;
0625   hcal_depthMatch1 = 0;
0626   hcal_depthMatchHot1 = 0;
0627   hcal_edepth2 = 0;
0628   hcal_activeL2 = 0;
0629   hcal_edepthHot2 = 0;
0630   hcal_activeHotL2 = 0;
0631   hcal_edepthCorrect2 = 0;
0632   hcal_edepthHotCorrect2 = 0;
0633   hcal_cdepthHot2 = 0;
0634   hcal_cdepthHotBG2 = 0;
0635   hcal_depthMatch2 = 0;
0636   hcal_depthMatchHot2 = 0;
0637   hcal_edepth3 = 0;
0638   hcal_activeL3 = 0;
0639   hcal_edepthHot3 = 0;
0640   hcal_activeHotL3 = 0;
0641   hcal_edepthCorrect3 = 0;
0642   hcal_edepthHotCorrect3 = 0;
0643   hcal_cdepthHot3 = 0;
0644   hcal_cdepthHotBG3 = 0;
0645   hcal_depthMatch3 = 0;
0646   hcal_depthMatchHot3 = 0;
0647   hcal_edepth4 = 0;
0648   hcal_activeL4 = 0;
0649   hcal_edepthHot4 = 0;
0650   hcal_activeHotL4 = 0;
0651   hcal_edepthCorrect4 = 0;
0652   hcal_edepthHotCorrect4 = 0;
0653   hcal_cdepthHot4 = 0;
0654   hcal_cdepthHotBG4 = 0;
0655   hcal_depthMatch4 = 0;
0656   hcal_depthMatchHot4 = 0;
0657   hcal_edepth5 = 0;
0658   hcal_activeL5 = 0;
0659   hcal_edepthHot5 = 0;
0660   hcal_activeHotL5 = 0;
0661   hcal_edepthCorrect5 = 0;
0662   hcal_edepthHotCorrect5 = 0;
0663   hcal_cdepthHot5 = 0;
0664   hcal_cdepthHotBG5 = 0;
0665   hcal_depthMatch5 = 0;
0666   hcal_depthMatchHot5 = 0;
0667   hcal_edepth6 = 0;
0668   hcal_activeL6 = 0;
0669   hcal_edepthHot6 = 0;
0670   hcal_activeHotL6 = 0;
0671   hcal_edepthCorrect6 = 0;
0672   hcal_edepthHotCorrect6 = 0;
0673   hcal_cdepthHot6 = 0;
0674   hcal_cdepthHotBG6 = 0;
0675   hcal_depthMatch6 = 0;
0676   hcal_depthMatchHot6 = 0;
0677   hcal_edepth7 = 0;
0678   hcal_activeL7 = 0;
0679   hcal_edepthHot7 = 0;
0680   hcal_activeHotL7 = 0;
0681   hcal_edepthCorrect7 = 0;
0682   hcal_edepthHotCorrect7 = 0;
0683   hcal_cdepthHot7 = 0;
0684   hcal_cdepthHotBG7 = 0;
0685   hcal_depthMatch7 = 0;
0686   hcal_depthMatchHot7 = 0;
0687   activeLength = 0;
0688   activeLengthHot = 0;
0689   hltresults = 0;
0690   all_triggers = 0;
0691   // Set branch addresses and branch pointers
0692   if (!tree)
0693     return;
0694   fChain = tree;
0695   fCurrent = -1;
0696   fChain->SetMakeClass(1);
0697 
0698   fChain->SetBranchAddress("Event_No", &Event_No, &b_Event_No);
0699   fChain->SetBranchAddress("Run_No", &Run_No, &b_Run_No);
0700   fChain->SetBranchAddress("LumiNumber", &LumiNumber, &b_LumiNumber);
0701   fChain->SetBranchAddress("BXNumber", &BXNumber, &b_BXNumber);
0702   fChain->SetBranchAddress("GoodVertex", &GoodVertex, &b_GoodVertex);
0703   fChain->SetBranchAddress("PF_Muon", &PF_Muon, &b_PF_Muon);
0704   fChain->SetBranchAddress("Global_Muon", &Global_Muon, &b_Global_Muon);
0705   fChain->SetBranchAddress("Tracker_muon", &Tracker_muon, &b_Tracker_muon);
0706   fChain->SetBranchAddress("MuonIsTight", &MuonIsTight, &b_MuonIsTight);
0707   fChain->SetBranchAddress("MuonIsMedium", &MuonIsMedium, &b_MuonIsMedium);
0708   fChain->SetBranchAddress("pt_of_muon", &pt_of_muon, &b_pt_of_muon);
0709   fChain->SetBranchAddress("eta_of_muon", &eta_of_muon, &b_eta_of_muon);
0710   fChain->SetBranchAddress("phi_of_muon", &phi_of_muon, &b_phi_of_muon);
0711   fChain->SetBranchAddress("energy_of_muon", &energy_of_muon, &b_energy_of_muon);
0712   fChain->SetBranchAddress("p_of_muon", &p_of_muon, &b_p_of_muon);
0713   fChain->SetBranchAddress("muon_trkKink", &muon_trkKink, &b_muon_trkKink);
0714   fChain->SetBranchAddress("muon_chi2LocalPosition", &muon_chi2LocalPosition, &b_muon_chi2LocalPosition);
0715   fChain->SetBranchAddress("muon_segComp", &muon_segComp, &b_muon_segComp);
0716   fChain->SetBranchAddress("TrackerLayer", &TrackerLayer, &b_TrackerLayer);
0717   fChain->SetBranchAddress("NumPixelLayers", &NumPixelLayers, &b_NumPixelLayers);
0718   fChain->SetBranchAddress("InnerTrackPixelHits", &InnerTrackPixelHits, &b_InnerTrackPixelHits);
0719   fChain->SetBranchAddress("innerTrack", &innerTrack, &b_innerTrack);
0720   fChain->SetBranchAddress("chiTracker", &chiTracker, &b_chiTracker);
0721   fChain->SetBranchAddress("DxyTracker", &DxyTracker, &b_DxyTracker);
0722   fChain->SetBranchAddress("DzTracker", &DzTracker, &b_DzTracker);
0723   fChain->SetBranchAddress("innerTrackpt", &innerTrackpt, &b_innerTrackpt);
0724   fChain->SetBranchAddress("innerTracketa", &innerTracketa, &b_innerTracketa);
0725   fChain->SetBranchAddress("innerTrackphi", &innerTrackphi, &b_innerTrackphi);
0726   fChain->SetBranchAddress("tight_validFraction", &tight_validFraction, &b_tight_validFraction);
0727   fChain->SetBranchAddress("OuterTrack", &OuterTrack, &b_OuterTrack);
0728   fChain->SetBranchAddress("OuterTrackPt", &OuterTrackPt, &b_OuterTrackPt);
0729   fChain->SetBranchAddress("OuterTrackEta", &OuterTrackEta, &b_OuterTrackEta);
0730   fChain->SetBranchAddress("OuterTrackPhi", &OuterTrackPhi, &b_OuterTrackPhi);
0731   fChain->SetBranchAddress("OuterTrackChi", &OuterTrackChi, &b_OuterTrackChi);
0732   fChain->SetBranchAddress("OuterTrackHits", &OuterTrackHits, &b_OuterTrackHits);
0733   fChain->SetBranchAddress("OuterTrackRHits", &OuterTrackRHits, &b_OuterTrackRHits);
0734   fChain->SetBranchAddress("GlobalTrack", &GlobalTrack, &b_GlobalTrack);
0735   fChain->SetBranchAddress("GlobalTrckPt", &GlobalTrckPt, &b_GlobalTrckPt);
0736   fChain->SetBranchAddress("GlobalTrckEta", &GlobalTrckEta, &b_GlobalTrckEta);
0737   fChain->SetBranchAddress("GlobalTrckPhi", &GlobalTrckPhi, &b_GlobalTrckPhi);
0738   fChain->SetBranchAddress("Global_Muon_Hits", &Global_Muon_Hits, &b_Global_Muon_Hits);
0739   fChain->SetBranchAddress("MatchedStations", &MatchedStations, &b_MatchedStations);
0740   fChain->SetBranchAddress("GlobTrack_Chi", &GlobTrack_Chi, &b_GlobTrack_Chi);
0741   fChain->SetBranchAddress(
0742       "Tight_LongitudinalImpactparameter", &Tight_LongitudinalImpactparameter, &b_Tight_LongitudinalImpactparameter);
0743   fChain->SetBranchAddress("Tight_TransImpactparameter", &Tight_TransImpactparameter, &b_Tight_TransImpactparameter);
0744   fChain->SetBranchAddress("IsolationR04", &IsolationR04, &b_IsolationR04);
0745   fChain->SetBranchAddress("IsolationR03", &IsolationR03, &b_IsolationR03);
0746   fChain->SetBranchAddress("ecal_3into3", &ecal_3into3, &b_ecal_3into3);
0747   fChain->SetBranchAddress("hcal_3into3", &hcal_3into3, &b_hcal_3into3);
0748   fChain->SetBranchAddress("tracker_3into3", &tracker_3into3, &b_tracker_3into3);
0749   fChain->SetBranchAddress("matchedId", &matchedId, &b_matchedId);
0750   fChain->SetBranchAddress("hcal_cellHot", &hcal_cellHot, &b_hcal_cellHot);
0751   fChain->SetBranchAddress("ecal_3x3", &ecal_3x3, &b_ecal_3x3);
0752   fChain->SetBranchAddress("hcal_1x1", &hcal_1x1, &b_hcal_1x1);
0753   fChain->SetBranchAddress("ecal_detID", &ecal_detID, &b_ecal_detID);
0754   fChain->SetBranchAddress("hcal_detID", &hcal_detID, &b_hcal_detID);
0755   fChain->SetBranchAddress("ehcal_detID", &ehcal_detID, &b_ehcal_detID);
0756   fChain->SetBranchAddress("hcal_ieta", &hcal_ieta, &b_hcal_ieta);
0757   fChain->SetBranchAddress("hcal_iphi", &hcal_iphi, &b_hcal_iphi);
0758   fChain->SetBranchAddress("hcal_edepth1", &hcal_edepth1, &b_hcal_edepth1);
0759   fChain->SetBranchAddress("hcal_activeL1", &hcal_activeL1, &b_hcal_activeL1);
0760   fChain->SetBranchAddress("hcal_edepthHot1", &hcal_edepthHot1, &b_hcal_edepthHot1);
0761   fChain->SetBranchAddress("hcal_activeHotL1", &hcal_activeHotL1, &b_hcal_activeHotL1);
0762   fChain->SetBranchAddress("hcal_edepthCorrect1", &hcal_edepthCorrect1, &b_hcal_edepthCorrect1);
0763   fChain->SetBranchAddress("hcal_edepthHotCorrect1", &hcal_edepthHotCorrect1, &b_hcal_edepthHotCorrect1);
0764   fChain->SetBranchAddress("hcal_cdepthHot1", &hcal_cdepthHot1, &b_hcal_cdepthHot1);
0765   fChain->SetBranchAddress("hcal_cdepthHotBG1", &hcal_cdepthHotBG1, &b_hcal_cdepthHotBG1);
0766   fChain->SetBranchAddress("hcal_depthMatch1", &hcal_depthMatch1, &b_hcal_depthMatch1);
0767   fChain->SetBranchAddress("hcal_depthMatchHot1", &hcal_depthMatchHot1, &b_hcal_depthMatchHot1);
0768   fChain->SetBranchAddress("hcal_edepth2", &hcal_edepth2, &b_hcal_edepth2);
0769   fChain->SetBranchAddress("hcal_activeL2", &hcal_activeL2, &b_hcal_activeL2);
0770   fChain->SetBranchAddress("hcal_edepthHot2", &hcal_edepthHot2, &b_hcal_edepthHot2);
0771   fChain->SetBranchAddress("hcal_activeHotL2", &hcal_activeHotL2, &b_hcal_activeHotL2);
0772   fChain->SetBranchAddress("hcal_edepthCorrect2", &hcal_edepthCorrect2, &b_hcal_edepthCorrect2);
0773   fChain->SetBranchAddress("hcal_edepthHotCorrect2", &hcal_edepthHotCorrect2, &b_hcal_edepthHotCorrect2);
0774   fChain->SetBranchAddress("hcal_cdepthHot2", &hcal_cdepthHot2, &b_hcal_cdepthHot2);
0775   fChain->SetBranchAddress("hcal_cdepthHotBG2", &hcal_cdepthHotBG2, &b_hcal_cdepthHotBG2);
0776   fChain->SetBranchAddress("hcal_depthMatch2", &hcal_depthMatch2, &b_hcal_depthMatch2);
0777   fChain->SetBranchAddress("hcal_depthMatchHot2", &hcal_depthMatchHot2, &b_hcal_depthMatchHot2);
0778   fChain->SetBranchAddress("hcal_edepth3", &hcal_edepth3, &b_hcal_edepth3);
0779   fChain->SetBranchAddress("hcal_activeL3", &hcal_activeL3, &b_hcal_activeL3);
0780   fChain->SetBranchAddress("hcal_edepthHot3", &hcal_edepthHot3, &b_hcal_edepthHot3);
0781   fChain->SetBranchAddress("hcal_activeHotL3", &hcal_activeHotL3, &b_hcal_activeHotL3);
0782   fChain->SetBranchAddress("hcal_edepthCorrect3", &hcal_edepthCorrect3, &b_hcal_edepthCorrect3);
0783   fChain->SetBranchAddress("hcal_edepthHotCorrect3", &hcal_edepthHotCorrect3, &b_hcal_edepthHotCorrect3);
0784   fChain->SetBranchAddress("hcal_cdepthHot3", &hcal_cdepthHot3, &b_hcal_cdepthHot3);
0785   fChain->SetBranchAddress("hcal_cdepthHotBG3", &hcal_cdepthHotBG3, &b_hcal_cdepthHotBG3);
0786   fChain->SetBranchAddress("hcal_depthMatch3", &hcal_depthMatch3, &b_hcal_depthMatch3);
0787   fChain->SetBranchAddress("hcal_depthMatchHot3", &hcal_depthMatchHot3, &b_hcal_depthMatchHot3);
0788   fChain->SetBranchAddress("hcal_edepth4", &hcal_edepth4, &b_hcal_edepth4);
0789   fChain->SetBranchAddress("hcal_activeL4", &hcal_activeL4, &b_hcal_activeL4);
0790   fChain->SetBranchAddress("hcal_edepthHot4", &hcal_edepthHot4, &b_hcal_edepthHot4);
0791   fChain->SetBranchAddress("hcal_activeHotL4", &hcal_activeHotL4, &b_hcal_activeHotL4);
0792   fChain->SetBranchAddress("hcal_edepthCorrect4", &hcal_edepthCorrect4, &b_hcal_edepthCorrect4);
0793   fChain->SetBranchAddress("hcal_edepthHotCorrect4", &hcal_edepthHotCorrect4, &b_hcal_edepthHotCorrect4);
0794   fChain->SetBranchAddress("hcal_cdepthHot4", &hcal_cdepthHot4, &b_hcal_cdepthHot4);
0795   fChain->SetBranchAddress("hcal_cdepthHotBG4", &hcal_cdepthHotBG4, &b_hcal_cdepthHotBG4);
0796   fChain->SetBranchAddress("hcal_depthMatch4", &hcal_depthMatch4, &b_hcal_depthMatch4);
0797   fChain->SetBranchAddress("hcal_depthMatchHot4", &hcal_depthMatchHot4, &b_hcal_depthMatchHot4);
0798   fChain->SetBranchAddress("hcal_edepth5", &hcal_edepth5, &b_hcal_edepth5);
0799   fChain->SetBranchAddress("hcal_activeL5", &hcal_activeL5, &b_hcal_activeL5);
0800   fChain->SetBranchAddress("hcal_edepthHot5", &hcal_edepthHot5, &b_hcal_edepthHot5);
0801   fChain->SetBranchAddress("hcal_activeHotL5", &hcal_activeHotL5, &b_hcal_activeHotL5);
0802   fChain->SetBranchAddress("hcal_edepthCorrect5", &hcal_edepthCorrect5, &b_hcal_edepthCorrect5);
0803   fChain->SetBranchAddress("hcal_edepthHotCorrect5", &hcal_edepthHotCorrect5, &b_hcal_edepthHotCorrect5);
0804   fChain->SetBranchAddress("hcal_cdepthHot5", &hcal_cdepthHot5, &b_hcal_cdepthHot5);
0805   fChain->SetBranchAddress("hcal_cdepthHotBG5", &hcal_cdepthHotBG5, &b_hcal_cdepthHotBG5);
0806   fChain->SetBranchAddress("hcal_depthMatch5", &hcal_depthMatch5, &b_hcal_depthMatch5);
0807   fChain->SetBranchAddress("hcal_depthMatchHot5", &hcal_depthMatchHot5, &b_hcal_depthMatchHot5);
0808   fChain->SetBranchAddress("hcal_edepth6", &hcal_edepth6, &b_hcal_edepth6);
0809   fChain->SetBranchAddress("hcal_activeL6", &hcal_activeL6, &b_hcal_activeL6);
0810   fChain->SetBranchAddress("hcal_edepthHot6", &hcal_edepthHot6, &b_hcal_edepthHot6);
0811   fChain->SetBranchAddress("hcal_activeHotL6", &hcal_activeHotL6, &b_hcal_activeHotL6);
0812   fChain->SetBranchAddress("hcal_edepthCorrect6", &hcal_edepthCorrect6, &b_hcal_edepthCorrect6);
0813   fChain->SetBranchAddress("hcal_edepthHotCorrect6", &hcal_edepthHotCorrect6, &b_hcal_edepthHotCorrect6);
0814   fChain->SetBranchAddress("hcal_cdepthHot6", &hcal_cdepthHot6, &b_hcal_cdepthHot6);
0815   fChain->SetBranchAddress("hcal_cdepthHotBG6", &hcal_cdepthHotBG6, &b_hcal_cdepthHotBG6);
0816   fChain->SetBranchAddress("hcal_depthMatch6", &hcal_depthMatch6, &b_hcal_depthMatch6);
0817   fChain->SetBranchAddress("hcal_depthMatchHot6", &hcal_depthMatchHot6, &b_hcal_depthMatchHot6);
0818   fChain->SetBranchAddress("hcal_edepth7", &hcal_edepth7, &b_hcal_edepth7);
0819   fChain->SetBranchAddress("hcal_activeL7", &hcal_activeL7, &b_hcal_activeL7);
0820   fChain->SetBranchAddress("hcal_edepthHot7", &hcal_edepthHot7, &b_hcal_edepthHot7);
0821   fChain->SetBranchAddress("hcal_activeHotL7", &hcal_activeHotL7, &b_hcal_activeHotL7);
0822   fChain->SetBranchAddress("hcal_edepthCorrect7", &hcal_edepthCorrect7, &b_hcal_edepthCorrect7);
0823   fChain->SetBranchAddress("hcal_edepthHotCorrect7", &hcal_edepthHotCorrect7, &b_hcal_edepthHotCorrect7);
0824   fChain->SetBranchAddress("hcal_cdepthHot7", &hcal_cdepthHot7, &b_hcal_cdepthHot7);
0825   fChain->SetBranchAddress("hcal_cdepthHotBG7", &hcal_cdepthHotBG7, &b_hcal_cdepthHotBG7);
0826   fChain->SetBranchAddress("hcal_depthMatch7", &hcal_depthMatch7, &b_hcal_depthMatch7);
0827   fChain->SetBranchAddress("hcal_depthMatchHot7", &hcal_depthMatchHot7, &b_hcal_depthMatchHot7);
0828   fChain->SetBranchAddress("activeLength", &activeLength, &b_activeLength);
0829   fChain->SetBranchAddress("activeLengthHot", &activeLengthHot, &b_activeLengthHot);
0830   fChain->SetBranchAddress("hltresults", &hltresults, &b_hltresults);
0831   fChain->SetBranchAddress("all_triggers", &all_triggers, &b_all_triggers);
0832 
0833   Notify();
0834 }
0835 
0836 void HBHEMuonOfflineAnalyzer::Loop() {
0837   //declarations
0838   if (fChain == 0)
0839     return;
0840 
0841   Long64_t nentries = fChain->GetEntriesFast();
0842 
0843   if (debug_)
0844     std::cout << "nevent = " << nentries << std::endl;
0845 
0846   Long64_t nbytes = 0, nb = 0;
0847 
0848   for (Long64_t jentry = 0; jentry < nentries; jentry++) {
0849     Long64_t ientry = LoadTree(jentry);
0850     if (ientry < 0)
0851       break;
0852     nb = fChain->GetEntry(jentry);
0853     nbytes += nb;
0854     if ((int)(Run_No) < runLo_ || (int)(Run_No) > runHi_)
0855       continue;
0856     if (debug_)
0857       std::cout << "Run " << Run_No << " Event " << Event_No << " Muon pt " << pt_of_muon << std::endl;
0858 
0859     bool loose(false), soft(false), tight(false), pcut(false), ptcut(false);
0860     t_ene.clear();
0861     t_enec.clear();
0862     t_charge.clear();
0863     t_actln.clear();
0864     t_depth.clear();
0865 
0866     if (debug_)
0867       std::cout << "ecal_det_id " << ecal_detID << std::endl;
0868 
0869     int typeEcal, etaXEcal, phiYEcal, zsideEcal, planeEcal, stripEcal;
0870     etaPhiEcal(ecal_detID, typeEcal, zsideEcal, etaXEcal, phiYEcal, planeEcal, stripEcal);
0871     double etaEcal = (etaXEcal - 0.5) * zsideEcal;
0872     double phiEcal = phiYEcal - 0.5;
0873 
0874     if (debug_)
0875       std::cout << "hcal_det_id " << std::hex << hcal_detID << std::dec;
0876 
0877     int etaHcal, phiHcal, depthHcal;
0878     etaPhiHcal(hcal_detID, etaHcal, phiHcal, depthHcal);
0879 
0880     int eta = (etaHcal > 0) ? etaHcal - 1 : -(1 + etaHcal);
0881     double etaXHcal = (etaHcal > 0) ? etaHcal - 0.5 : etaHcal + 0.5;
0882     int nDepth = nDepthBins(eta + 1, phiHcal);
0883     int nPhi = nPhiBins(eta + 1);
0884     int PHI = (nPhi > 36) ? (phiHcal - 1) : (phiHcal - 1) / 2;
0885     double phiYHcal = (phiHcal - 0.5);
0886     t_ieta = etaHcal;
0887     t_iphi = PHI;
0888     t_p = p_of_muon;
0889     t_ediff = hcal_3into3 - hcal_1x1;
0890     t_nvtx = GoodVertex;
0891     if (p_of_muon > 20)
0892       pcut = true;
0893     if (pt_of_muon > 20)
0894       ptcut = true;
0895     if (looseMuon())
0896       loose = true;
0897     if (softMuon())
0898       soft = true;
0899     if (tightMuon())
0900       tight = true;
0901 
0902     if (debug_)
0903       std::cout << " etaHcal " << etaHcal << ":" << etaXHcal << " phiHcal " << phiHcal << ":" << phiYHcal << ":" << PHI
0904                 << " Depth " << nDepth << " Muon Pt " << pt_of_muon << " Isol " << IsolationR04 << std::endl;
0905 
0906     for (int cut = 0; cut < nCut_; ++cut) {
0907       bool select(false);
0908       if (cut == 0)
0909         select = tightMuon();
0910       else if (cut == 1)
0911         select = softMuon();
0912       else
0913         select = looseMuon();
0914 
0915       if (select && ((eta + 1) >= etaMin_) && ((eta + 1) <= etaMax_)) {
0916         h_Pt_Muon[cut]->Fill(pt_of_muon);
0917         h_Eta_Muon[cut]->Fill(eta_of_muon);
0918         h_Phi_Muon[cut]->Fill(phi_of_muon);
0919         h_PF_Muon[cut]->Fill(PF_Muon);
0920         h_GlobTrack_Chi[cut]->Fill(GlobTrack_Chi);
0921         h_Global_Muon_Hits[cut]->Fill(Global_Muon_Hits);
0922         h_MatchedStations[cut]->Fill(MatchedStations);
0923         h_Tight_TransImpactparameter[cut]->Fill(Tight_TransImpactparameter);
0924         h_Tight_LongitudinalImpactparameter[cut]->Fill(Tight_LongitudinalImpactparameter);
0925         h_InnerTrackPixelHits[cut]->Fill(InnerTrackPixelHits);
0926         h_TrackerLayer[cut]->Fill(TrackerLayer);
0927         h_IsolationR04[cut]->Fill(IsolationR04);
0928         h_Global_Muon[cut]->Fill(Global_Muon);
0929 
0930         h_TransImpactParameter[cut]->Fill(Tight_TransImpactparameter);
0931         h_LongImpactParameter[cut]->Fill(Tight_LongitudinalImpactparameter);
0932 
0933         //in Phi Bins
0934         if (((phi_of_muon) >= -1.5) || ((phi_of_muon) <= 0.5)) {
0935           h_TransImpactParameterBin1[cut]->Fill(Tight_TransImpactparameter);
0936           h_LongImpactParameterBin1[cut]->Fill(Tight_LongitudinalImpactparameter);
0937           h_2D_Bin1[cut]->Fill(Tight_TransImpactparameter, Tight_LongitudinalImpactparameter);
0938         }
0939 
0940         if ((phi_of_muon > 0.5) || (phi_of_muon < -1.5)) {
0941           h_TransImpactParameterBin2[cut]->Fill(Tight_TransImpactparameter);
0942           h_LongImpactParameterBin2[cut]->Fill(Tight_LongitudinalImpactparameter);
0943           h_2D_Bin2[cut]->Fill(Tight_TransImpactparameter, Tight_LongitudinalImpactparameter);
0944         }
0945 
0946         h_ecal_energy[cut]->Fill(ecal_3into3);
0947         h_3x3_ecal[cut]->Fill(ecal_3x3);
0948         h_Eta_ecal[cut]->Fill(eta_of_muon, ecal_3x3);
0949         h_Phi_ecal[cut]->Fill(phi_of_muon, ecal_3x3);
0950         h_MuonHittingEcal[cut]->Fill(typeEcal);
0951         if (typeEcal == 1) {
0952           h_EtaX_ecal[cut]->Fill(etaEcal, ecal_3x3);
0953           h_PhiY_ecal[cut]->Fill(phiEcal, ecal_3x3);
0954         }
0955 
0956         h_hcal_energy[cut]->Fill(hcal_3into3);
0957         h_1x1_hcal[cut]->Fill(hcal_1x1);
0958         h_EtaX_hcal[cut]->Fill(etaXHcal, hcal_1x1);
0959         h_PhiY_hcal[cut]->Fill(phiYHcal, hcal_1x1);
0960         h_HotCell[cut]->Fill(hcal_cellHot);
0961         if (mergeDepth_) {
0962           double en1(0), en2(0), energyFill(0), chargeS(0), chargeBG(0);
0963           double enh(0), enc(0);
0964           for (int dep = 0; dep < nDepth; ++dep) {
0965             double enb(0), enu(0), eh0(0), ec0(0), chgS(0), chgB(0), actL(0);
0966             getEnergy(dep, enb, enu, eh0, ec0, chgS, chgB, actL);
0967             en1 += ((useCorrect_) ? enu : enb);
0968             en2 += ((useCorrect_) ? ec0 : eh0);
0969             enh += (eh0);
0970             enc += (ec0);
0971             energyFill += (actL);
0972             chargeS += (chgS);
0973             chargeBG += (chgB);
0974           }
0975           int ind = (etaHcal > 0) ? indxEta[eta][0][PHI] : 1 + indxEta[eta][0][PHI];
0976           if (debug_)  // || eta==15 || eta==17)
0977             std::cout << "Matched Id " << matchedId << " Hot " << hcal_cellHot << " eta " << etaHcal << ":" << eta
0978                       << " phi " << phiHcal << ":" << PHI << " Index " << ind << " E " << en1 << ":" << en2 << ":"
0979                       << enh << ":" << enc << " L " << energyFill << " Charge " << chargeS << ":" << chargeBG
0980                       << std::endl;
0981           if (!(matchedId))
0982             continue;
0983           if (hcal_cellHot == 1) {
0984             if (energyFill > 0) {
0985               h_Hot_MuonEnergy_hcal_HotCell[cut][ind]->Fill(en2);
0986               h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[cut][ind]->Fill(en2 / energyFill);
0987               h_active_length_Fill[cut][ind]->Fill(energyFill);
0988               h_p_muon_ineta[cut][ind]->Fill(p_of_muon);
0989               h_charge_signal[cut][ind]->Fill(chargeS);
0990               h_charge_bg[cut][ind]->Fill(chargeBG);
0991 
0992               t_ene.push_back(enh);
0993               t_enec.push_back(enc);
0994               t_charge.push_back(chargeS);
0995               t_actln.push_back(energyFill);
0996               t_depth.push_back(0);
0997 
0998               outtree_->Fill();
0999             }
1000           }
1001         } else {
1002           bool fillTree(false);
1003           for (int dep = 0; dep < nDepth; ++dep) {
1004             if (debug_)
1005               std::cout << "dep:" << dep << std::endl;
1006 
1007             double energyFill(0), chargeS(-9999), chargeBG(-9999);
1008             double enh(-9999), enc(-9999), enb(0), enu(0);
1009             bool ok1 = getEnergy(dep, enb, enu, enh, enc, chargeS, chargeBG, energyFill);
1010             double en1 = ((useCorrect_) ? enu : enb);
1011             double en2 = ((useCorrect_) ? enc : enh);
1012             if (debug_)
1013               std::cout << "Hello in " << dep + 1 << " " << en1 << ":" << en2 << ":" << energyFill << std::endl;
1014 
1015             bool ok2 = ok1;
1016 
1017             if (debug_)
1018               std::cout << "Before Index " << ok1 << ":" << ok2 << std::endl;
1019 
1020             int ind = (etaHcal > 0) ? indxEta[eta][dep][PHI] : 1 + indxEta[eta][dep][PHI];
1021             if (debug_)  // || eta==15 || eta==17)
1022               std::cout << "Matched Id " << matchedId << " Hot " << hcal_cellHot << " eta " << etaHcal << ":" << eta
1023                         << " phi " << phiHcal << ":" << PHI << " depth " << dep << " Index " << ind << " E " << en1
1024                         << ":" << en2 << ":" << enh << ":" << enc << " L " << energyFill << " Charge " << chargeS << ":"
1025                         << chargeBG << std::endl;
1026             if (!(matchedId))
1027               continue;
1028             if (ok1) {
1029               if (debug_)
1030                 std::cout << "enter ok1" << std::endl;
1031 
1032               if (hcal_cellHot == 1) {
1033                 if (energyFill > 0) {
1034                   h_Hot_MuonEnergy_hcal_HotCell[cut][ind]->Fill(en2);
1035                   h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[cut][ind]->Fill(en2 / energyFill);
1036                   h_active_length_Fill[cut][ind]->Fill(energyFill);
1037                   h_p_muon_ineta[cut][ind]->Fill(p_of_muon);
1038                   h_charge_signal[cut][ind]->Fill(chargeS);
1039                   h_charge_bg[cut][ind]->Fill(chargeBG);
1040                   t_ene.push_back(enh);
1041                   t_enec.push_back(enc);
1042                   t_charge.push_back(chargeS);
1043                   t_actln.push_back(energyFill);
1044                   // added depth vector AmanKalsi
1045                   t_depth.push_back(dep);
1046                   fillTree = true;
1047                 } else {
1048                   t_ene.push_back(-999.0);
1049                   t_enec.push_back(-999.0);
1050                   t_charge.push_back(-999.0);
1051                   t_actln.push_back(-999.0);
1052                   t_depth.push_back(-999.0);
1053                 }
1054                 if (debug_)
1055                   std::cout << "enter hot cell" << std::endl;
1056               }
1057             }
1058 
1059             if (ok2) {
1060               if (debug_)
1061                 std::cout << "enter ok2" << std::endl;
1062               if (hcal_cellHot != 1) {
1063               }
1064             }
1065 
1066             if (debug_)
1067               std::cout << "ETA \t" << eta << "DEPTH \t" << dep << std::endl;
1068           }
1069           if (fillTree)
1070             outtree_->Fill();
1071         }
1072       }
1073     }
1074     int evtype(0);
1075     if (pcut)
1076       evtype += 1;
1077     if (ptcut)
1078       evtype += 2;
1079     if (loose)
1080       evtype += 4;
1081     if (soft)
1082       evtype += 8;
1083     if (tight)
1084       evtype += 16;
1085     h_evtype->Fill(evtype);
1086   }
1087   close();
1088 }
1089 
1090 Bool_t HBHEMuonOfflineAnalyzer::Notify() {
1091   // The Notify() function is called when a new file is opened. This
1092   // can be either for a new TTree in a TChain or when when a new TTree
1093   // is started when using PROOF. It is normally not necessary to make changes
1094   // to the generated code, but the routine can be extended by the
1095   // user if needed. The return value is currently not used.
1096 
1097   return kTRUE;
1098 }
1099 
1100 void HBHEMuonOfflineAnalyzer::Show(Long64_t entry) {
1101   // Print contents of entry.
1102   // If entry is not specified, print current entry
1103   if (!fChain)
1104     return;
1105   fChain->Show(entry);
1106 }
1107 
1108 bool HBHEMuonOfflineAnalyzer::fillChain(TChain *chain, const char *inputFileList) {
1109   std::string fname(inputFileList);
1110   if (fname.substr(fname.size() - 5, 5) == ".root") {
1111     chain->Add(fname.c_str());
1112   } else {
1113     std::ifstream infile(inputFileList);
1114     if (!infile.is_open()) {
1115       std::cout << "** ERROR: Can't open '" << inputFileList << "' for input" << std::endl;
1116       return false;
1117     }
1118     while (1) {
1119       infile >> fname;
1120       if (!infile.good())
1121         break;
1122       chain->Add(fname.c_str());
1123     }
1124     infile.close();
1125   }
1126   std::cout << "No. of Entries in this tree : " << chain->GetEntries() << std::endl;
1127   return true;
1128 }
1129 
1130 bool HBHEMuonOfflineAnalyzer::readCorr(const char *infile) {
1131   std::ifstream fInput(infile);
1132   unsigned int ncorr(0), all(0), good(0);
1133   if (!fInput.good()) {
1134     std::cout << "Cannot open file " << infile << std::endl;
1135   } else {
1136     char buffer[1024];
1137     while (fInput.getline(buffer, 1024)) {
1138       ++all;
1139       std::string bufferString(buffer);
1140       if (bufferString.substr(0, 5) == "#IOVs") {
1141         std::vector<std::string> items = splitString(bufferString.substr(6));
1142         ncorr = items.size() - 1;
1143         for (unsigned int n = 0; n < ncorr; ++n) {
1144           int run = std::atoi(items[n].c_str());
1145           runlow_.push_back(run);
1146         }
1147         std::cout << ncorr << ":" << runlow_.size() << " Run ranges" << std::endl;
1148         for (unsigned int n = 0; n < runlow_.size(); ++n)
1149           std::cout << " [" << n << "] " << runlow_[n];
1150         std::cout << std::endl;
1151       } else if (buffer[0] == '#') {
1152         continue;  //ignore other comments
1153       } else {
1154         std::vector<std::string> items = splitString(bufferString);
1155         if (items.size() != ncorr + 3) {
1156           std::cout << "Ignore  line: " << buffer << std::endl;
1157         } else {
1158           ++good;
1159           int ieta = std::atoi(items[0].c_str());
1160           int iphi = std::atoi(items[1].c_str());
1161           int depth = std::atoi(items[2].c_str());
1162           unsigned int id = getDetIdHBHE(ieta, iphi, depth);
1163           for (unsigned int n = 0; n < ncorr; ++n) {
1164             float corrf = std::atof(items[n + 3].c_str());
1165             if (n < nmax_)
1166               corrFac_[n][id] = corrf;
1167           }
1168           if (debug_) {
1169             std::cout << "ID " << std::hex << id << std::dec << ":" << id << " (eta " << ieta << " phi " << iphi
1170                       << " depth " << depth << ")";
1171             for (unsigned int n = 0; n < ncorr; ++n)
1172               std::cout << " " << corrFac_[n][id];
1173             std::cout << std::endl;
1174           }
1175         }
1176       }
1177     }
1178     fInput.close();
1179     std::cout << "Reads total of " << all << " and " << good << " good records" << std::endl;
1180   }
1181   return (good > 0);
1182 }
1183 
1184 void HBHEMuonOfflineAnalyzer::bookHistograms(const char *fname) {
1185   std::cout << "BookHistograms" << std::endl;
1186   output_file = TFile::Open(fname, "RECREATE");
1187   output_file->cd();
1188   outtree_ = new TTree("Lep_Tree", "Lep_Tree");
1189   outtree_->Branch("t_ieta", &t_ieta);
1190   outtree_->Branch("t_iphi", &t_iphi);
1191   outtree_->Branch("t_nvtx", &t_nvtx);
1192   outtree_->Branch("t_p", &t_p);
1193   outtree_->Branch("t_ediff", &t_ediff);
1194   outtree_->Branch("t_ene", &t_ene);
1195   outtree_->Branch("t_enec", &t_enec);
1196   outtree_->Branch("t_charge", &t_charge);
1197   outtree_->Branch("t_actln", &t_actln);
1198   outtree_->Branch("t_depth", &t_depth);
1199 
1200   std::string type[] = {"tight"};  //,"soft","loose"};
1201   char name[128], title[500];
1202 
1203   nHist = 0;
1204   for (int eta = etaMin_; eta <= etaMax_; ++eta) {
1205     int nDepth = nDepthBins(eta, -1);
1206     int nPhi = nPhiBins(eta);
1207     //std::cout<<"problem 2"<<std::endl;
1208     //std::cout<<"Eta: "<<eta<<" nDepth "<<nDepth<<" nPhi "<<nPhi<<std::endl;
1209     for (int depth = 0; depth < nDepth; depth++) {
1210       //std::cout<<"problem 3"<<std::endl;
1211       //std::cout<<"Eta:"<<eta<<"Depth:"<<depth<<std::endl;
1212       for (int PHI = 0; PHI < nPhi; ++PHI) {
1213         //std::cout<<"Eta:"<<eta<<"Depth:"<<depth<<"PHI:"<<PHI<<std::endl;
1214         indxEta[eta - 1][depth][PHI] = nHist;
1215         nHist += 2;
1216       }
1217     }
1218   }
1219   if (nHist >= maxHist_) {
1220     std::cout << "Problem here " << nHist << ":" << maxHist_ << std::endl;
1221   }
1222   //    TDirectory *d_output_file[nCut_][29];
1223   //output_file->cd();
1224   h_evtype = new TH1D("EvType", "Event Type", 100, 0, 100);
1225   for (int i = 0; i < nCut_; ++i) {
1226     sprintf(name, "h_Pt_Muon_%s", type[i].c_str());
1227     sprintf(title, "p_{T} of %s muons (GeV)", type[i].c_str());
1228     h_Pt_Muon[i] = new TH1D(name, title, 100, 0, 200);
1229 
1230     sprintf(name, "h_Eta_Muon_%s", type[i].c_str());
1231     sprintf(title, "#eta of %s muons", type[i].c_str());
1232     h_Eta_Muon[i] = new TH1D(name, title, 50, -2.5, 2.5);
1233 
1234     sprintf(name, "h_Phi_Muon_%s", type[i].c_str());
1235     sprintf(title, "#phi of %s muons", type[i].c_str());
1236     h_Phi_Muon[i] = new TH1D(name, title, 100, -3.1415926, 3.1415926);
1237 
1238     sprintf(name, "h_P_Muon_%s", type[i].c_str());
1239     sprintf(title, "p of %s muons (GeV)", type[i].c_str());
1240     h_P_Muon[i] = new TH1D(name, title, 100, 0, 200);
1241 
1242     sprintf(name, "h_PF_Muon_%s", type[i].c_str());
1243     sprintf(title, "PF %s muons (GeV)", type[i].c_str());
1244     h_PF_Muon[i] = new TH1D(name, title, 2, 0, 2);
1245 
1246     sprintf(name, "h_Global_Muon_Chi2_%s", type[i].c_str());
1247     sprintf(title, "Chi2 Global %s muons (GeV)", type[i].c_str());
1248     h_GlobTrack_Chi[i] = new TH1D(name, title, 15, 0, 15);
1249 
1250     sprintf(name, "h_Global_Muon_Hits_%s", type[i].c_str());
1251     sprintf(title, "Global Hits %s muons (GeV)", type[i].c_str());
1252     h_Global_Muon_Hits[i] = new TH1D(name, title, 10, 0, 10);
1253 
1254     sprintf(name, "h_Matched_Stations_%s", type[i].c_str());
1255     sprintf(title, "Matched Stations %s muons (GeV)", type[i].c_str());
1256     h_MatchedStations[i] = new TH1D(name, title, 10, 0, 10);
1257 
1258     sprintf(name, "h_Transverse_ImpactParameter_%s", type[i].c_str());
1259     sprintf(title, "Transverse_ImpactParameter of %s muons (GeV)", type[i].c_str());
1260     h_Tight_TransImpactparameter[i] = new TH1D(name, title, 50, 0, 10);
1261 
1262     sprintf(name, "h_Longitudinal_ImpactParameter_%s", type[i].c_str());
1263     sprintf(title, "Longitudinal_ImpactParameter of %s muons (GeV)", type[i].c_str());
1264     h_Tight_LongitudinalImpactparameter[i] = new TH1D(name, title, 20, 0, 10);
1265 
1266     sprintf(name, "h_InnerTrack_PixelHits_%s", type[i].c_str());
1267     sprintf(title, "InnerTrack_PixelHits of %s muons (GeV)", type[i].c_str());
1268     h_InnerTrackPixelHits[i] = new TH1D(name, title, 20, 0, 20);
1269 
1270     sprintf(name, "h_TrackLayers_%s", type[i].c_str());
1271     sprintf(title, "No. of Tracker Layers of %s muons (GeV)", type[i].c_str());
1272     h_TrackerLayer[i] = new TH1D(name, title, 20, 0, 20);
1273     ;
1274 
1275     sprintf(name, "h_IsolationR04_%s", type[i].c_str());
1276     sprintf(title, "IsolationR04 %s muons (GeV)", type[i].c_str());
1277     h_IsolationR04[i] = new TH1D(name, title, 45, 0, 5);
1278     ;
1279 
1280     sprintf(name, "h_Global_Muon_%s", type[i].c_str());
1281     sprintf(title, "Global %s muons (GeV)", type[i].c_str());
1282     h_Global_Muon[i] = new TH1D(name, title, 2, 0, 2);
1283 
1284     sprintf(name, "h_TransImpactParameter_%s", type[i].c_str());
1285     sprintf(title, "TransImpactParameter of %s muons (GeV)", type[i].c_str());
1286     h_TransImpactParameter[i] = new TH1D(name, title, 100, 0, 0.5);
1287 
1288     sprintf(name, "h_TransImpactParameterBin1_%s", type[i].c_str());
1289     sprintf(title, "TransImpactParameter of %s muons (GeV) in -1.5 <= #phi <= 0.5", type[i].c_str());
1290     h_TransImpactParameterBin1[i] = new TH1D(name, title, 100, 0, 0.5);
1291 
1292     sprintf(name, "h_TransImpactParameterBin2_%s", type[i].c_str());
1293     sprintf(title, "TransImpactParameter of %s muons (GeV) in #phi> 0.5 and #phi< -1.5 ", type[i].c_str());
1294     h_TransImpactParameterBin2[i] = new TH1D(name, title, 100, 0, 0.5);
1295     //
1296     sprintf(name, "h_LongImpactParameter_%s", type[i].c_str());
1297     sprintf(title, "LongImpactParameter of %s muons (GeV)", type[i].c_str());
1298     h_LongImpactParameter[i] = new TH1D(name, title, 100, 0, 30);
1299 
1300     sprintf(name, "h_LongImpactParameterBin1_%s", type[i].c_str());
1301     sprintf(title, "LongImpactParameter of %s muons (GeV) in -1.5 <= #phi <= 0.5", type[i].c_str());
1302     h_LongImpactParameterBin1[i] = new TH1D(name, title, 100, 0, 30);
1303 
1304     sprintf(name, "h_LongImpactParameterBin2_%s", type[i].c_str());
1305     sprintf(title, "LongImpactParameter of %s muons (GeV) in #phi> 0.5 and #phi< -1.5 ", type[i].c_str());
1306     h_LongImpactParameterBin2[i] = new TH1D(name, title, 100, 0, 30);
1307 
1308     sprintf(name, "h_2D_Bin1_%s", type[i].c_str());
1309     sprintf(title, "Trans/Long ImpactParameter of %s muons (GeV) in -1.5 <= #phi< 0.5 ", type[i].c_str());
1310     h_2D_Bin1[i] = new TH2D(name, title, 100, 0, 0.5, 100, 0, 30);
1311 
1312     sprintf(name, "h_2D_Bin2_%s", type[i].c_str());
1313     sprintf(title, "Trans/Long ImpactParameter of %s muons (GeV) in #phi> 0.5 and #phi< -1.5 ", type[i].c_str());
1314     h_2D_Bin2[i] = new TH2D(name, title, 100, 0, 0.5, 100, 0, 30);
1315 
1316     sprintf(name, "h_ecal_energy_%s", type[i].c_str());
1317     sprintf(title, "ECAL energy for %s muons", type[i].c_str());
1318     h_ecal_energy[i] = new TH1D(name, title, 1000, -10.0, 90.0);
1319 
1320     sprintf(name, "h_hcal_energy_%s", type[i].c_str());
1321     sprintf(title, "HCAL energy for %s muons", type[i].c_str());
1322     h_hcal_energy[i] = new TH1D(name, title, 500, -10.0, 90.0);
1323 
1324     sprintf(name, "h_3x3_ecal_%s", type[i].c_str());
1325     sprintf(title, "ECAL energy in 3x3 for %s muons", type[i].c_str());
1326     h_3x3_ecal[i] = new TH1D(name, title, 1000, -10.0, 90.0);
1327 
1328     sprintf(name, "h_1x1_hcal_%s", type[i].c_str());
1329     sprintf(title, "HCAL energy in 1x1 for %s muons", type[i].c_str());
1330     h_1x1_hcal[i] = new TH1D(name, title, 500, -10.0, 90.0);
1331 
1332     sprintf(name, "h_EtaX_hcal_%s", type[i].c_str());
1333     sprintf(title, "HCAL energy as a function of i#eta for %s muons", type[i].c_str());
1334     h_EtaX_hcal[i] = new TProfile(name, title, 60, -30.0, 30.0);
1335 
1336     sprintf(name, "h_PhiY_hcal_%s", type[i].c_str());
1337     sprintf(title, "HCAL energy as a function of i#phi for %s muons", type[i].c_str());
1338     h_PhiY_hcal[i] = new TProfile(name, title, 72, 0, 72);
1339 
1340     sprintf(name, "h_EtaX_ecal_%s", type[i].c_str());
1341     sprintf(title, "EB energy as a function of i#eta for %s muons", type[i].c_str());
1342     h_EtaX_ecal[i] = new TProfile(name, title, 170, -85.0, 85.0);
1343 
1344     sprintf(name, "h_PhiY_ecal_%s", type[i].c_str());
1345     sprintf(title, "EB energy as a function of i#phi for %s muons", type[i].c_str());
1346     h_PhiY_ecal[i] = new TProfile(name, title, 360, 0, 360);
1347 
1348     sprintf(name, "h_Eta_ecal_%s", type[i].c_str());
1349     sprintf(title, "ECAL energy as a function of #eta for %s muons", type[i].c_str());
1350     h_Eta_ecal[i] = new TProfile(name, title, 100, -2.5, 2.5);
1351 
1352     sprintf(name, "h_Phi_ecal_%s", type[i].c_str());
1353     sprintf(title, "ECAL energy as a function of #phi for %s muons", type[i].c_str());
1354     h_Phi_ecal[i] = new TProfile(name, title, 100, -3.1415926, 3.1415926);
1355 
1356     sprintf(name, "h_MuonHittingEcal_%s", type[i].c_str());
1357     sprintf(title, "%s muons hitting ECAL", type[i].c_str());
1358     h_MuonHittingEcal[i] = new TH1D(name, title, 100, 0, 5.0);
1359 
1360     sprintf(name, "h_HotCell_%s", type[i].c_str());
1361     sprintf(title, "Hot cell for %s muons", type[i].c_str());
1362     h_HotCell[i] = new TH1D(name, title, 100, 0, 2);
1363 
1364     //      output_file->cd();
1365     for (int eta = etaMin_; eta <= etaMax_; ++eta) {
1366       int nDepth = nDepthBins(eta, -1);
1367       int nPhi = nPhiBins(eta);
1368       //sprintf(name, "Dir_muon_type_%s_ieta%d",type[i].c_str(), eta);
1369       //d_output_file[i][eta]= output_file->mkdir(name);
1370       //output_file->cd(name);
1371       //d_output_file[i][eta]->cd();
1372       for (int depth = 0; depth < nDepth; ++depth) {
1373         for (int PHI = 0; PHI < nPhi; ++PHI) {
1374           int PHI0 = (nPhi == 72) ? PHI + 1 : 2 * PHI + 1;
1375           int ih = indxEta[eta - 1][depth][PHI];
1376           if (debug_)
1377             std::cout << "eta:" << eta << " depth:" << depth << " PHI:" << PHI << ":" << PHI0 << " ih:" << ih
1378                       << std::endl;
1379 
1380           sprintf(name, "h_Hot_MuonEnergy_hc_%d_%d_%d_%s_HotCell", eta, (depth + 1), PHI0, type[i].c_str());
1381           sprintf(title,
1382                   "HCAL energy in hot tower (i#eta=%d, depth=%d, i#phi = %d) for extrapolated %s muons (Hot Cell)",
1383                   eta,
1384                   (depth + 1),
1385                   PHI0,
1386                   type[i].c_str());
1387           h_Hot_MuonEnergy_hcal_HotCell[i][ih] = new TH1D(name, title, 4000, 0.0, 10.0);
1388           h_Hot_MuonEnergy_hcal_HotCell[i][ih]->Sumw2();
1389 
1390           sprintf(
1391               name, "h_Hot_MuonEnergy_hc_%d_%d_%d_%s_HotCell_ByActiveLength", eta, (depth + 1), PHI0, type[i].c_str());
1392           sprintf(title,
1393                   "HCAL energy in hot tower (i#eta=%d, depth=%d, i#phi = %d) for extrapolated %s muons (Hot Cell) "
1394                   "divided by Active Length",
1395                   eta,
1396                   (depth + 1),
1397                   PHI0,
1398                   type[i].c_str());
1399           h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[i][ih] = new TH1D(name, title, 4000, 0.0, 10.0);
1400           h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[i][ih]->Sumw2();
1401 
1402           sprintf(name, "h_active_length_Fill_%d_%d_%d_%s", eta, (depth + 1), PHI0, type[i].c_str());
1403           sprintf(title, "active_length%d_%d_%d_%s", eta, (depth + 1), PHI0, type[i].c_str());
1404           h_active_length_Fill[i][ih] = new TH1D(name, title, 20, 0, 20);
1405           h_active_length_Fill[i][ih]->Sumw2();
1406 
1407           sprintf(name, "h_p_muon_in_%d_%d_%d_%s", eta, (depth + 1), PHI0, type[i].c_str());
1408           sprintf(title, "p_muon_in%d_%d_%d_%s", eta, (depth + 1), PHI0, type[i].c_str());
1409           h_p_muon_ineta[i][ih] = new TH1D(name, title, 500, 0, 500);
1410           h_p_muon_ineta[i][ih]->Sumw2();
1411 
1412           sprintf(name, "h_charge_signal_in_%d_%d_%d_%s", eta, (depth + 1), PHI0, type[i].c_str());
1413           sprintf(name, "charge_signal_in_%d_%d_%d_%s", eta, (depth + 1), PHI0, type[i].c_str());
1414           h_charge_signal[i][ih] = new TH1D(name, title, 500, 0, 500);
1415           h_charge_signal[i][ih]->Sumw2();
1416 
1417           sprintf(name, "h_charge_bg_in_%d_%d_%d_%s", eta, (depth + 1), PHI0, type[i].c_str());
1418           sprintf(name, "charge_bg_in_%d_%d_%d_%s", eta, (depth + 1), PHI0, type[i].c_str());
1419           h_charge_bg[i][ih] = new TH1D(name, title, 500, 0, 500);
1420           h_charge_bg[i][ih]->Sumw2();
1421 
1422           ih++;
1423 
1424           sprintf(name, "h_Hot_MuonEnergy_hc_%d_%d_%d_%s_HotCell", -eta, (depth + 1), PHI0, type[i].c_str());
1425           sprintf(title,
1426                   "HCAL energy in hot tower (i#eta=%d, depth=%d, i#phi = %d) for extrapolated %s muons (Hot Cell)",
1427                   -eta,
1428                   (depth + 1),
1429                   PHI0,
1430                   type[i].c_str());
1431           h_Hot_MuonEnergy_hcal_HotCell[i][ih] = new TH1D(name, title, 4000, 0.0, 10.0);
1432           h_Hot_MuonEnergy_hcal_HotCell[i][ih]->Sumw2();
1433 
1434           sprintf(
1435               name, "h_Hot_MuonEnergy_hc_%d_%d_%d_%s_HotCell_ByActiveLength", -eta, (depth + 1), PHI0, type[i].c_str());
1436           sprintf(title,
1437                   "HCAL energy in hot tower (i#eta=%d, depth=%d, i#phi=%d) for extrapolated %s muons (Hot Cell) "
1438                   "divided by Active Length",
1439                   -eta,
1440                   (depth + 1),
1441                   PHI0,
1442                   type[i].c_str());
1443           h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[i][ih] = new TH1D(name, title, 4000, 0.0, 10.0);
1444           h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[i][ih]->Sumw2();
1445 
1446           sprintf(name, "h_active_length_Fill_%d_%d_%d_%s", -eta, (depth + 1), PHI0, type[i].c_str());
1447           sprintf(title, "active_length%d_%d_%d_%s", -eta, (depth + 1), PHI0, type[i].c_str());
1448           h_active_length_Fill[i][ih] = new TH1D(name, title, 20, 0, 20);
1449           h_active_length_Fill[i][ih]->Sumw2();
1450 
1451           sprintf(name, "h_p_muon_in_%d_%d_%d_%s", -eta, (depth + 1), PHI0, type[i].c_str());
1452           sprintf(title, "p_muon_in%d_%d_%d_%s", -eta, (depth + 1), PHI0, type[i].c_str());
1453           h_p_muon_ineta[i][ih] = new TH1D(name, title, 500, 0, 500);
1454           h_p_muon_ineta[i][ih]->Sumw2();
1455 
1456           sprintf(name, "h_charge_signal_in_%d_%d_%d_%s", -eta, (depth + 1), PHI0, type[i].c_str());
1457           sprintf(name, "charge_signal_in_%d_%d_%d_%s", -eta, (depth + 1), PHI0, type[i].c_str());
1458           h_charge_signal[i][ih] = new TH1D(name, title, 500, 0, 500);
1459           h_charge_signal[i][ih]->Sumw2();
1460 
1461           sprintf(name, "h_charge_bg_in_%d_%d_%d_%s", -eta, (depth + 1), PHI0, type[i].c_str());
1462           sprintf(name, "charge_bg_in_%d_%d_%d_%s", -eta, (depth + 1), PHI0, type[i].c_str());
1463           h_charge_bg[i][ih] = new TH1D(name, title, 500, 0, 500);
1464           h_charge_bg[i][ih]->Sumw2();
1465         }
1466       }
1467       //output_file->cd();
1468     }
1469   }
1470   //output_file->cd();
1471 }
1472 
1473 bool HBHEMuonOfflineAnalyzer::getEnergy(
1474     int dep, double &enb, double &enu, double &enh, double &enc, double &chgS, double &chgB, double &actL) {
1475   double cfac(1.0);
1476   bool flag(true);
1477   if (cFactor_) {
1478     int ieta = hcal_ieta;
1479     int iphi = hcal_iphi;
1480     unsigned int detId = getDetIdHBHE(ieta, iphi, dep + 1);
1481     cfac = getCorr(Run_No, detId);
1482   }
1483   if (dep == 0) {
1484     enb = cfac * hcal_edepth1;
1485     enu = cfac * hcal_edepthCorrect1;
1486     enh = cfac * hcal_edepthHot1;
1487     enc = cfac * hcal_edepthHotCorrect1;
1488     chgS = hcal_cdepthHot1;
1489     actL = hcal_activeHotL1;
1490     chgB = hcal_cdepthHotBG1;
1491   } else if (dep == 1) {
1492     enb = cfac * hcal_edepth2;
1493     enu = cfac * hcal_edepthCorrect2;
1494     enh = cfac * hcal_edepthHot2;
1495     enc = cfac * hcal_edepthHotCorrect2;
1496     chgS = hcal_cdepthHot2;
1497     actL = hcal_activeHotL2;
1498     chgB = hcal_cdepthHotBG2;
1499   } else if (dep == 2) {
1500     enb = cfac * hcal_edepth3;
1501     enu = cfac * hcal_edepthCorrect3;
1502     enh = cfac * hcal_edepthHot3;
1503     enc = cfac * hcal_edepthHotCorrect3;
1504     chgS = hcal_cdepthHot3;
1505     actL = hcal_activeHotL3;
1506     chgB = hcal_cdepthHotBG3;
1507   } else if (dep == 3) {
1508     enb = cfac * hcal_edepth4;
1509     enu = cfac * hcal_edepthCorrect4;
1510     enh = cfac * hcal_edepthHot4;
1511     enc = cfac * hcal_edepthHotCorrect4;
1512     chgS = hcal_cdepthHot4;
1513     actL = hcal_activeHotL4;
1514     chgB = hcal_cdepthHotBG4;
1515   } else if (dep == 4) {
1516     enb = cfac * hcal_edepth5;
1517     enu = cfac * hcal_edepthCorrect5;
1518     enh = cfac * hcal_edepthHot5;
1519     enc = cfac * hcal_edepthHotCorrect5;
1520     chgS = hcal_cdepthHot5;
1521     actL = hcal_activeHotL5;
1522     chgB = hcal_cdepthHotBG5;
1523   } else if (dep == 5) {
1524     if (dep <= maxDepth_) {
1525       enb = cfac * hcal_edepth6;
1526       enu = cfac * hcal_edepthCorrect6;
1527       enh = cfac * hcal_edepthHot6;
1528       enc = cfac * hcal_edepthHotCorrect6;
1529       chgS = hcal_cdepthHot6;
1530       actL = hcal_activeHotL6;
1531       chgB = hcal_cdepthHotBG6;
1532     } else {
1533       enb = enu = enh = enc = chgS = actL = chgB = 0;
1534       flag = false;
1535     }
1536   } else if (dep == 6) {
1537     if (dep <= maxDepth_) {
1538       enb = cfac * hcal_edepth7;
1539       enu = cfac * hcal_edepthCorrect7;
1540       enh = cfac * hcal_edepthHot7;
1541       enc = cfac * hcal_edepthHotCorrect7;
1542       chgS = hcal_cdepthHot7;
1543       actL = hcal_activeHotL7;
1544       chgB = hcal_cdepthHotBG7;
1545     } else {
1546       enb = enu = enh = enc = chgS = actL = chgB = 0;
1547       flag = false;
1548     }
1549   } else {
1550     enb = enu = enh = enc = chgS = actL = chgB = 0;
1551     flag = false;
1552   }
1553   return flag;
1554 }
1555 
1556 bool HBHEMuonOfflineAnalyzer::looseMuon() {
1557   if (pt_of_muon > 20.) {
1558     if (mediumMuon2016()) {
1559       if (IsolationR04 < 0.25) {
1560         return true;
1561       }
1562     }
1563   }
1564   return false;
1565 }
1566 
1567 bool HBHEMuonOfflineAnalyzer::softMuon() {
1568   if (pt_of_muon > 20.) {
1569     if (mediumMuon2016()) {
1570       if (IsolationR03 < 0.10) {
1571         return true;
1572       }
1573     }
1574   }
1575   return false;
1576 }
1577 
1578 bool HBHEMuonOfflineAnalyzer::tightMuon() {
1579   if (pt_of_muon > 20.) {
1580     if (mediumMuon2016()) {
1581       if (IsolationR04 < 0.15) {
1582         return true;
1583       }
1584     }
1585   }
1586   return false;
1587 }
1588 
1589 bool HBHEMuonOfflineAnalyzer::mediumMuon2016() {
1590   bool medium16 = (((PF_Muon) && (Global_Muon || Tracker_muon)) && (tight_validFraction > 0.49));
1591   if (!medium16)
1592     return medium16;
1593 
1594   bool goodGlob = (Global_Muon && GlobTrack_Chi < 3 && muon_chi2LocalPosition < 12 && muon_trkKink < 20);
1595   medium16 = muon_segComp > (goodGlob ? 0.303 : 0.451);
1596   return medium16;
1597 }
1598 
1599 void HBHEMuonOfflineAnalyzer::etaPhiHcal(unsigned int detId, int &eta, int &phi, int &depth) {
1600   int zside, etaAbs;
1601   if ((detId & 0x1000000) == 0) {
1602     zside = (detId & 0x2000) ? (1) : (-1);
1603     etaAbs = (detId >> 7) & 0x3F;
1604     phi = detId & 0x7F;
1605     depth = (detId >> 14) & 0x1F;
1606   } else {
1607     zside = (detId & 0x80000) ? (1) : (-1);
1608     etaAbs = (detId >> 10) & 0x1FF;
1609     phi = detId & 0x3FF;
1610     depth = (detId >> 20) & 0xF;
1611   }
1612   eta = etaAbs * zside;
1613 }
1614 
1615 void HBHEMuonOfflineAnalyzer::etaPhiEcal(
1616     unsigned int detId, int &type, int &zside, int &etaX, int &phiY, int &plane, int &strip) {
1617   type = ((detId >> 25) & 0x7);
1618   // std::cout<<"type"<<type<<std::endl;
1619   plane = strip = 0;
1620   if (type == 1) {
1621     //Ecal Barrel
1622     zside = (detId & 0x10000) ? (1) : (-1);
1623     etaX = (detId >> 9) & 0x7F;
1624     phiY = detId & 0x1FF;
1625   } else if (type == 2) {
1626     zside = (detId & 0x4000) ? (1) : (-1);
1627     etaX = (detId >> 7) & 0x7F;
1628     phiY = (detId & 0x7F);
1629   } else if (type == 3) {
1630     zside = (detId & 0x80000) ? (1) : (-1);
1631     etaX = (detId >> 6) & 0x3F;
1632     /** get the sensor iy */
1633     phiY = (detId >> 12) & 0x3F;
1634     /** get the strip */
1635     plane = ((detId >> 18) & 0x1) + 1;
1636     strip = detId & 0x3F;
1637   } else {
1638     zside = etaX = phiY = 0;
1639   }
1640 }
1641 
1642 void HBHEMuonOfflineAnalyzer::calculateP(double pt, double eta, double &pM) {
1643   pM = (pt * cos(2 * (1 / atan(exp(eta)))));
1644 }
1645 
1646 void HBHEMuonOfflineAnalyzer::close() {
1647   output_file->cd();
1648   std::cout << "file yet to be Written" << std::endl;
1649   outtree_->Write();
1650   writeHistograms();
1651   std::cout << "file Written" << std::endl;
1652   output_file->Close();
1653   std::cout << "now doing return" << std::endl;
1654 }
1655 
1656 void HBHEMuonOfflineAnalyzer::writeHistograms() {
1657   //output_file->cd();
1658   std::string type[] = {"tight"};  //,"soft","loose"};
1659   char name[128];
1660 
1661   std::cout << "WriteHistograms" << std::endl;
1662   nHist = 0;
1663   for (int eta = etaMin_; eta <= etaMax_; ++eta) {
1664     int nDepth = nDepthBins(eta, -1);
1665     int nPhi = nPhiBins(eta);
1666     if (debug_)
1667       std::cout << "Eta:" << eta << " nDepths " << nDepth << " nPhis " << nPhi << std::endl;
1668     for (int depth = 0; depth < nDepth; ++depth) {
1669       if (debug_)
1670         std::cout << "Eta:" << eta << "Depth:" << depth << std::endl;
1671       for (int PHI = 0; PHI < nPhi; ++PHI) {
1672         indxEta[eta - 1][depth][PHI] = nHist;
1673         nHist += 2;
1674       }
1675     }
1676   }
1677 
1678   TDirectory *d_output_file[nCut_][29];
1679   h_evtype->Write();
1680   //output_file->cd();
1681   for (int i = 0; i < nCut_; ++i) {
1682     h_Pt_Muon[i]->Write();
1683     h_Eta_Muon[i]->Write();
1684     h_Phi_Muon[i]->Write();
1685 
1686     h_P_Muon[i]->Write();
1687     h_PF_Muon[i]->Write();
1688     h_GlobTrack_Chi[i]->Write();
1689     h_Global_Muon_Hits[i]->Write();
1690 
1691     h_MatchedStations[i]->Write();
1692     h_Tight_TransImpactparameter[i]->Write();
1693     h_Tight_LongitudinalImpactparameter[i]->Write();
1694 
1695     h_InnerTrackPixelHits[i]->Write();
1696     h_TrackerLayer[i]->Write();
1697     h_IsolationR04[i]->Write();
1698 
1699     h_Global_Muon[i]->Write();
1700     h_TransImpactParameter[i]->Write();
1701     ;
1702     h_TransImpactParameterBin1[i]->Write();
1703     h_TransImpactParameterBin2[i]->Write();
1704     //
1705     h_LongImpactParameter[i]->Write();
1706     h_LongImpactParameterBin1[i]->Write();
1707     h_LongImpactParameterBin2[i]->Write();
1708 
1709     h_ecal_energy[i]->Write();
1710     h_hcal_energy[i]->Write();
1711     ;
1712     h_3x3_ecal[i]->Write();
1713     h_1x1_hcal[i]->Write();
1714     ;
1715 
1716     h_EtaX_hcal[i]->Write();
1717     h_PhiY_hcal[i]->Write();
1718     ;
1719 
1720     h_EtaX_ecal[i]->Write();
1721     ;
1722     h_PhiY_ecal[i]->Write();
1723     ;
1724     h_Eta_ecal[i]->Write();
1725     ;
1726     h_Phi_ecal[i]->Write();
1727     ;
1728     h_MuonHittingEcal[i]->Write();
1729     ;
1730     h_HotCell[i]->Write();
1731     ;
1732 
1733     output_file->cd();
1734     for (int eta = etaMin_; eta <= etaMax_; ++eta) {
1735       int nDepth = nDepthBins(eta, -1);
1736       int nPhi = nPhiBins(eta);
1737       sprintf(name, "Dir_muon_type_%s_ieta%d", type[i].c_str(), eta);
1738       d_output_file[i][eta - 1] = output_file->mkdir(name);
1739       //output_file->cd(name);
1740       d_output_file[i][eta - 1]->cd();
1741       for (int depth = 0; depth < nDepth; ++depth) {
1742         for (int PHI = 0; PHI < nPhi; ++PHI) {
1743           //    std::cout<<"eta:"<<eta<<"depth:"<<depth<<"PHI:"<<PHI<<std::endl;
1744           int ih = indxEta[eta - 1][depth][PHI];
1745           //    std::cout<<"ih:"<<ih<<std::endl;
1746           h_Hot_MuonEnergy_hcal_HotCell[i][ih]->Write();
1747 
1748           h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[i][ih]->Write();
1749           h_active_length_Fill[i][ih]->Write();
1750           h_p_muon_ineta[i][ih]->Write();
1751           h_charge_signal[i][ih]->Write();
1752           h_charge_bg[i][ih]->Write();
1753           ih++;
1754           h_Hot_MuonEnergy_hcal_HotCell[i][ih]->Write();
1755           h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[i][ih]->Write();
1756           h_active_length_Fill[i][ih]->Write();
1757           h_p_muon_ineta[i][ih]->Write();
1758           h_charge_signal[i][ih]->Write();
1759           h_charge_bg[i][ih]->Write();
1760         }
1761       }
1762       output_file->cd();
1763     }
1764   }
1765   output_file->cd();
1766 }
1767 
1768 int HBHEMuonOfflineAnalyzer::nDepthBins(int eta, int phi) {
1769   // Run 1 scenario
1770   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};
1771   // Run 2 scenario from 2018
1772   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};
1773   // Run 3 scenario
1774   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};
1775   // Run 4 scenario
1776   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};
1777   // for a test scenario with multi depth segmentation considered during Run 1
1778   //    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};
1779   // modeLHC_ = 0 --> nbin defined maxDepthHB/HE
1780   //          = 1 -->      corresponds to Run 1 (valid till 2016)
1781   //          = 2 -->      corresponds to Run 2 (2018 geometry)
1782   //          = 3 -->      corresponds to Run 3 (post LS2)
1783   //          = 4 -->      corresponds to 2017 (Plan 1)
1784   //          = 5 -->      corresponds to Run 4 (post LS3)
1785   int nbin(0);
1786   if (modeLHC_ == 0) {
1787     if (eta <= 15) {
1788       nbin = maxDepthHB_;
1789     } else if (eta == 16) {
1790       nbin = 4;
1791     } else {
1792       nbin = maxDepthHE_;
1793     }
1794   } else if (modeLHC_ == 1) {
1795     nbin = nDepthR1[eta - 1];
1796   } else if (modeLHC_ == 2) {
1797     nbin = nDepthR2[eta - 1];
1798   } else if (modeLHC_ == 3) {
1799     nbin = nDepthR3[eta - 1];
1800   } else if (modeLHC_ == 4) {
1801     if (phi > 0) {
1802       if (eta >= 16 && phi >= 63 && phi <= 66) {
1803         nbin = nDepthR2[eta - 1];
1804       } else {
1805         nbin = nDepthR1[eta - 1];
1806       }
1807     } else {
1808       if (eta >= 16) {
1809         nbin = (nDepthR2[eta - 1] > nDepthR1[eta - 1]) ? nDepthR2[eta - 1] : nDepthR1[eta - 1];
1810       } else {
1811         nbin = nDepthR1[eta - 1];
1812       }
1813     }
1814   } else {
1815     if (eta > 0 && eta < 30) {
1816       nbin = nDepthR4[eta - 1];
1817     } else {
1818       nbin = nDepthR4[28];
1819     }
1820   }
1821   return nbin;
1822 }
1823 
1824 int HBHEMuonOfflineAnalyzer::nPhiBins(int eta) {
1825   int nphi = (eta <= 20) ? 72 : 36;
1826   if (modeLHC_ == 5 && eta > 16)
1827     nphi = 360;
1828   return nphi;
1829 }
1830 
1831 float HBHEMuonOfflineAnalyzer::getCorr(int run, unsigned int id) {
1832   float cfac(1.0);
1833   int ip(-1);
1834   for (unsigned int k = 0; k < runlow_.size(); ++k) {
1835     unsigned int i = runlow_.size() - k - 1;
1836     if (run >= runlow_[i]) {
1837       ip = (int)(i);
1838       break;
1839     }
1840   }
1841   if (debug_) {
1842     std::cout << "Run " << run << " Perdiod " << ip << std::endl;
1843   }
1844   unsigned idx = correctDetId(id);
1845   if (ip >= 0) {
1846     std::map<unsigned int, float>::iterator itr = corrFac_[ip].find(idx);
1847     if (itr != corrFac_[ip].end())
1848       cfac = itr->second;
1849   }
1850   if (debug_) {
1851     int subdet, zside, ieta, iphi, depth;
1852     unpackDetId(idx, subdet, zside, ieta, iphi, depth);
1853     std::cout << "ID " << std::hex << id << std::dec << " (Sub " << subdet << " eta " << zside * ieta << " phi " << iphi
1854               << " depth " << depth << ")  Factor " << cfac << std::endl;
1855   }
1856   return cfac;
1857 }
1858 
1859 std::vector<std::string> HBHEMuonOfflineAnalyzer::splitString(const std::string &fLine) {
1860   std::vector<std::string> result;
1861   int start = 0;
1862   bool empty = true;
1863   for (unsigned i = 0; i <= fLine.size(); i++) {
1864     if (fLine[i] == ' ' || i == fLine.size()) {
1865       if (!empty) {
1866         std::string item(fLine, start, i - start);
1867         result.push_back(item);
1868         empty = true;
1869       }
1870       start = i + 1;
1871     } else {
1872       if (empty)
1873         empty = false;
1874     }
1875   }
1876   return result;
1877 }
1878 
1879 unsigned int HBHEMuonOfflineAnalyzer::getDetIdHBHE(int ieta, int iphi, int depth) {
1880   int eta = std::abs(ieta);
1881   int subdet(0);
1882   if (eta > 16)
1883     subdet = 2;
1884   else if (eta == 16 && depth > 2)
1885     subdet = 2;
1886   else
1887     subdet = 1;
1888   return getDetId(subdet, ieta, iphi, depth);
1889 }
1890 
1891 unsigned int HBHEMuonOfflineAnalyzer::getDetId(int subdet, int ieta, int iphi, int depth) {
1892   // All numbers used here are described as masks/offsets in
1893   // DataFormats/HcalDetId/interface/HcalDetId.h
1894   unsigned int id = ((4 << 28) | ((subdet & 0x7) << 25));
1895   id |= ((0x1000000) | ((depth & 0xF) << 20) | ((ieta > 0) ? (0x80000 | (ieta << 10)) : ((-ieta) << 10)) |
1896          (iphi & 0x3FF));
1897   return id;
1898 }
1899 
1900 unsigned int HBHEMuonOfflineAnalyzer::correctDetId(const unsigned int &detId) {
1901   int subdet, ieta, zside, depth, iphi;
1902   unpackDetId(detId, subdet, zside, ieta, iphi, depth);
1903   if (subdet == 0) {
1904     if (ieta > 16)
1905       subdet = 2;
1906     else if (ieta == 16 && depth > 2)
1907       subdet = 2;
1908     else
1909       subdet = 1;
1910   }
1911   unsigned int id = getDetId(subdet, ieta * zside, iphi, depth);
1912   if ((id != detId) && debug_) {
1913     std::cout << "Correct Id " << std::hex << detId << " to " << id << std::dec << "(Sub " << subdet << " eta "
1914               << ieta * zside << " phi " << iphi << " depth " << depth << ")" << std::endl;
1915   }
1916   return id;
1917 }
1918 
1919 void HBHEMuonOfflineAnalyzer::unpackDetId(
1920     unsigned int detId, int &subdet, int &zside, int &ieta, int &iphi, int &depth) {
1921   // The maskings are defined in DataFormats/DetId/interface/DetId.h
1922   //                      and in DataFormats/HcalDetId/interface/HcalDetId.h
1923   // The macro does not invoke the classes there and use them
1924   subdet = ((detId >> 25) & (0x7));
1925   if ((detId & 0x1000000) == 0) {
1926     ieta = ((detId >> 7) & 0x3F);
1927     zside = (detId & 0x2000) ? (1) : (-1);
1928     depth = ((detId >> 14) & 0x1F);
1929     iphi = (detId & 0x3F);
1930   } else {
1931     ieta = ((detId >> 10) & 0x1FF);
1932     zside = (detId & 0x80000) ? (1) : (-1);
1933     depth = ((detId >> 20) & 0xF);
1934     iphi = (detId & 0x3FF);
1935   }
1936 }