Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 11:59:00

0001 #include <fstream>
0002 #include <vector>
0003 #include <TTree.h>
0004 
0005 // user include files
0006 #include "FWCore/Framework/interface/Frameworkfwd.h"
0007 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0008 #include "FWCore/Framework/interface/Event.h"
0009 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0010 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0011 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0012 #include "FWCore/ServiceRegistry/interface/Service.h"
0013 #include "FWCore/Common/interface/TriggerNames.h"
0014 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0015 
0016 #include "DataFormats/HcalDetId/interface/HcalDetId.h"
0017 #include "DataFormats/MuonReco/interface/Muon.h"
0018 #include "DataFormats/MuonReco/interface/MuonFwd.h"
0019 #include "DataFormats/MuonReco/interface/MuonSelectors.h"
0020 #include "DataFormats/VertexReco/interface/Vertex.h"
0021 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0022 #include "DataFormats/TrackReco/interface/Track.h"
0023 
0024 //////////////trigger info////////////////////////////////////
0025 
0026 #include "DataFormats/Common/interface/TriggerResults.h"
0027 #include "DataFormats/HLTReco/interface/TriggerObject.h"
0028 #include "DataFormats/HLTReco/interface/TriggerEvent.h"
0029 #include "DataFormats/EcalDetId/interface/EcalSubdetector.h"
0030 #include "DataFormats/EcalDetId/interface/EBDetId.h"
0031 #include "DataFormats/EcalDetId/interface/EEDetId.h"
0032 #include "DataFormats/HcalDetId/interface/HcalDetId.h"
0033 #include "DataFormats/EcalRecHit/interface/EcalRecHitCollections.h"
0034 #include "DataFormats/HcalRecHit/interface/HcalRecHitCollections.h"
0035 
0036 #include "HLTrigger/HLTcore/interface/HLTConfigProvider.h"
0037 #include "HLTrigger/HLTcore/interface/HLTConfigData.h"
0038 
0039 #include "CondFormats/HcalObjects/interface/HcalRespCorrs.h"
0040 #include "CondFormats/DataRecord/interface/HcalRespCorrsRcd.h"
0041 
0042 #include "CondFormats/DataRecord/interface/EcalPFRecHitThresholdsRcd.h"
0043 #include "CondFormats/EcalObjects/interface/EcalPFRecHitThresholds.h"
0044 #include "CondFormats/DataRecord/interface/EcalChannelStatusRcd.h"
0045 #include "RecoLocalCalo/EcalRecAlgos/interface/EcalSeverityLevelAlgo.h"
0046 #include "RecoLocalCalo/EcalRecAlgos/interface/EcalSeverityLevelAlgoRcd.h"
0047 
0048 #include "CalibFormats/HcalObjects/interface/HcalCalibrations.h"
0049 #include "CalibFormats/HcalObjects/interface/HcalDbService.h"
0050 #include "CalibFormats/HcalObjects/interface/HcalDbRecord.h"
0051 
0052 #include "Calibration/IsolatedParticles/interface/CaloPropagateTrack.h"
0053 #include "Calibration/IsolatedParticles/interface/eECALMatrix.h"
0054 #include "Calibration/IsolatedParticles/interface/eHCALMatrix.h"
0055 
0056 #include "Geometry/CaloGeometry/interface/CaloGeometry.h"
0057 #include "Geometry/Records/interface/CaloGeometryRecord.h"
0058 #include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h"
0059 #include "Geometry/Records/interface/CaloTopologyRecord.h"
0060 #include "Geometry/CaloTopology/interface/CaloSubdetectorTopology.h"
0061 #include "Geometry/CaloTopology/interface/HcalTopology.h"
0062 #include "Geometry/CaloTopology/interface/CaloTopology.h"
0063 #include "Geometry/HcalCommonData/interface/HcalDDDRecConstants.h"
0064 #include "Geometry/HcalTowerAlgo/interface/HcalGeometry.h"
0065 #include "MagneticField/Engine/interface/MagneticField.h"
0066 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0067 
0068 //#define EDM_ML_DEBUG
0069 
0070 class HcalHBHEMuonAnalyzer : public edm::one::EDAnalyzer<edm::one::WatchRuns, edm::one::SharedResources> {
0071 public:
0072   explicit HcalHBHEMuonAnalyzer(const edm::ParameterSet&);
0073 
0074   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0075 
0076 private:
0077   void beginJob() override;
0078   void analyze(edm::Event const&, edm::EventSetup const&) override;
0079   void beginRun(edm::Run const&, edm::EventSetup const&) override;
0080   void endRun(edm::Run const&, edm::EventSetup const&) override {}
0081   void clearVectors();
0082   int matchId(const HcalDetId&, const HcalDetId&);
0083   double activeLength(const DetId&, bool);
0084   bool isGoodVertex(const reco::Vertex& vtx);
0085   double respCorr(const DetId& id);
0086   double gainFactor(const HcalDbService* dbserv, const HcalDetId& id);
0087   int depth16HE(int ieta, int iphi);
0088   bool goodCell(const HcalDetId& hcid, const reco::Track* pTrack, const CaloGeometry* geo, const MagneticField* bField);
0089 
0090   // ----------member data ---------------------------
0091   HLTConfigProvider hltConfig_;
0092   const edm::InputTag hlTriggerResults_;
0093   const edm::InputTag labelEBRecHit_, labelEERecHit_, labelHBHERecHit_;
0094   const std::string labelVtx_, labelMuon_, fileInCorr_;
0095   const std::vector<std::string> triggers_;
0096   const double pMinMuon_;
0097   const int verbosity_, useRaw_;
0098   const bool unCorrect_, collapseDepth_, isItPlan1_;
0099   const bool ignoreHECorr_, isItPreRecHit_;
0100   const bool getCharge_, writeRespCorr_;
0101   const int maxDepth_;
0102   const std::string modnam_, procnm_;
0103   const bool usePFThresh_;
0104 
0105   const edm::EDGetTokenT<edm::TriggerResults> tok_trigRes_;
0106   const edm::EDGetTokenT<reco::VertexCollection> tok_Vtx_;
0107   const edm::EDGetTokenT<EcalRecHitCollection> tok_EB_;
0108   const edm::EDGetTokenT<EcalRecHitCollection> tok_EE_;
0109   const edm::EDGetTokenT<HBHERecHitCollection> tok_HBHE_;
0110   const edm::EDGetTokenT<reco::MuonCollection> tok_Muon_;
0111 
0112   const edm::ESGetToken<HcalDDDRecConstants, HcalRecNumberingRecord> tok_ddrec_;
0113   const edm::ESGetToken<HcalTopology, HcalRecNumberingRecord> tok_htopo_;
0114   const edm::ESGetToken<HcalRespCorrs, HcalRespCorrsRcd> tok_respcorr_;
0115   const edm::ESGetToken<CaloGeometry, CaloGeometryRecord> tok_geom_;
0116   const edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> tok_magField_;
0117   const edm::ESGetToken<EcalChannelStatus, EcalChannelStatusRcd> tok_chan_;
0118   const edm::ESGetToken<EcalSeverityLevelAlgo, EcalSeverityLevelAlgoRcd> tok_sevlv_;
0119   const edm::ESGetToken<CaloTopology, CaloTopologyRecord> tok_topo_;
0120   const edm::ESGetToken<HcalDbService, HcalDbRecord> tok_dbservice_;
0121   const edm::ESGetToken<EcalPFRecHitThresholds, EcalPFRecHitThresholdsRcd> tok_ecalPFRecHitThresholds_;
0122 
0123   const HcalDDDRecConstants* hdc_;
0124   const HcalTopology* theHBHETopology_;
0125   const CaloGeometry* geo_;
0126   HcalRespCorrs* respCorrs_;
0127 
0128   bool mergedDepth_, useMyCorr_;
0129   int kount_;
0130 
0131   //////////////////////////////////////////////////////
0132   static const int depthMax_ = 7;
0133   TTree* tree_;
0134   unsigned int runNumber_, eventNumber_, lumiNumber_, bxNumber_;
0135   unsigned int goodVertex_;
0136   bool muon_is_good_, muon_global_, muon_tracker_;
0137   bool muon_is_tight_, muon_is_medium_;
0138   double ptGlob_, etaGlob_, phiGlob_, energyMuon_, pMuon_;
0139   float muon_trkKink_, muon_chi2LocalPosition_, muon_segComp_;
0140   int trackerLayer_, numPixelLayers_, tight_PixelHits_;
0141   bool innerTrack_, outerTrack_, globalTrack_;
0142   double chiTracker_, dxyTracker_, dzTracker_;
0143   double innerTrackpt_, innerTracketa_, innerTrackphi_;
0144   double tight_validFraction_, outerTrackChi_;
0145   double outerTrackPt_, outerTrackEta_, outerTrackPhi_;
0146   int outerTrackHits_, outerTrackRHits_;
0147   double globalTrckPt_, globalTrckEta_, globalTrckPhi_;
0148   int globalMuonHits_, matchedStat_;
0149   double chiGlobal_, tight_LongPara_, tight_TransImpara_;
0150   double isolationR04_, isolationR03_;
0151   double ecalEnergy_, hcalEnergy_, hoEnergy_;
0152   bool matchedId_, hcalHot_;
0153   double ecal3x3Energy_, hcal1x1Energy_;
0154   unsigned int ecalDetId_, hcalDetId_, ehcalDetId_;
0155   int hcal_ieta_, hcal_iphi_;
0156   double hcalDepthEnergy_[depthMax_];
0157   double hcalDepthActiveLength_[depthMax_];
0158   double hcalDepthEnergyHot_[depthMax_];
0159   double hcalDepthActiveLengthHot_[depthMax_];
0160   double hcalDepthChargeHot_[depthMax_];
0161   double hcalDepthChargeHotBG_[depthMax_];
0162   double hcalDepthEnergyCorr_[depthMax_];
0163   double hcalDepthEnergyHotCorr_[depthMax_];
0164   bool hcalDepthMatch_[depthMax_];
0165   bool hcalDepthMatchHot_[depthMax_];
0166   double hcalActiveLength_, hcalActiveLengthHot_;
0167   std::vector<std::string> all_triggers_;
0168   std::vector<int> hltresults_;
0169 
0170   std::vector<HcalDDDRecConstants::HcalActiveLength> actHB, actHE;
0171   std::map<DetId, double> corrValue_;
0172   ////////////////////////////////////////////////////////////
0173 };
0174 
0175 HcalHBHEMuonAnalyzer::HcalHBHEMuonAnalyzer(const edm::ParameterSet& iConfig)
0176     : hlTriggerResults_(iConfig.getParameter<edm::InputTag>("hlTriggerResults")),
0177       labelEBRecHit_(iConfig.getParameter<edm::InputTag>("labelEBRecHit")),
0178       labelEERecHit_(iConfig.getParameter<edm::InputTag>("labelEERecHit")),
0179       labelHBHERecHit_(iConfig.getParameter<edm::InputTag>("labelHBHERecHit")),
0180       labelVtx_(iConfig.getParameter<std::string>("labelVertex")),
0181       labelMuon_(iConfig.getParameter<std::string>("labelMuon")),
0182       fileInCorr_(iConfig.getUntrackedParameter<std::string>("fileInCorr", "")),
0183       triggers_(iConfig.getParameter<std::vector<std::string>>("triggers")),
0184       pMinMuon_(iConfig.getParameter<double>("pMinMuon")),
0185       verbosity_(iConfig.getUntrackedParameter<int>("verbosity", 0)),
0186       useRaw_(iConfig.getParameter<int>("useRaw")),
0187       unCorrect_(iConfig.getParameter<bool>("unCorrect")),
0188       collapseDepth_(iConfig.getParameter<bool>("collapseDepth")),
0189       isItPlan1_(iConfig.getParameter<bool>("isItPlan1")),
0190       ignoreHECorr_(iConfig.getUntrackedParameter<bool>("ignoreHECorr", false)),
0191       isItPreRecHit_(iConfig.getUntrackedParameter<bool>("isItPreRecHit", false)),
0192       getCharge_(iConfig.getParameter<bool>("getCharge")),
0193       writeRespCorr_(iConfig.getUntrackedParameter<bool>("writeRespCorr", false)),
0194       maxDepth_(iConfig.getUntrackedParameter<int>("maxDepth", 7)),
0195       modnam_(iConfig.getUntrackedParameter<std::string>("moduleName", "")),
0196       procnm_(iConfig.getUntrackedParameter<std::string>("processName", "")),
0197       usePFThresh_(iConfig.getParameter<bool>("usePFThreshold")),
0198       tok_trigRes_(consumes<edm::TriggerResults>(hlTriggerResults_)),
0199       tok_Vtx_((modnam_.empty()) ? consumes<reco::VertexCollection>(labelVtx_)
0200                                  : consumes<reco::VertexCollection>(edm::InputTag(modnam_, labelVtx_, procnm_))),
0201       tok_EB_(consumes<EcalRecHitCollection>(labelEBRecHit_)),
0202       tok_EE_(consumes<EcalRecHitCollection>(labelEERecHit_)),
0203       tok_HBHE_(consumes<HBHERecHitCollection>(labelHBHERecHit_)),
0204       tok_Muon_((modnam_.empty()) ? consumes<reco::MuonCollection>(labelMuon_)
0205                                   : consumes<reco::MuonCollection>(edm::InputTag(modnam_, labelMuon_, procnm_))),
0206       tok_ddrec_(esConsumes<HcalDDDRecConstants, HcalRecNumberingRecord, edm::Transition::BeginRun>()),
0207       tok_htopo_(esConsumes<HcalTopology, HcalRecNumberingRecord, edm::Transition::BeginRun>()),
0208       tok_respcorr_(esConsumes<HcalRespCorrs, HcalRespCorrsRcd, edm::Transition::BeginRun>()),
0209       tok_geom_(esConsumes<CaloGeometry, CaloGeometryRecord, edm::Transition::BeginRun>()),
0210       tok_magField_(esConsumes<MagneticField, IdealMagneticFieldRecord>()),
0211       tok_chan_(esConsumes<EcalChannelStatus, EcalChannelStatusRcd>()),
0212       tok_sevlv_(esConsumes<EcalSeverityLevelAlgo, EcalSeverityLevelAlgoRcd>()),
0213       tok_topo_(esConsumes<CaloTopology, CaloTopologyRecord>()),
0214       tok_dbservice_(esConsumes<HcalDbService, HcalDbRecord>()),
0215       tok_ecalPFRecHitThresholds_(esConsumes<EcalPFRecHitThresholds, EcalPFRecHitThresholdsRcd>()),
0216       hdc_(nullptr),
0217       theHBHETopology_(nullptr),
0218       respCorrs_(nullptr) {
0219   usesResource(TFileService::kSharedResource);
0220   //now do what ever initialization is needed
0221   kount_ = 0;
0222   mergedDepth_ = (!isItPreRecHit_) || (collapseDepth_);
0223 
0224   if (modnam_.empty()) {
0225     edm::LogVerbatim("HBHEMuon") << "Labels used: Trig " << hlTriggerResults_ << " Vtx " << labelVtx_ << " EB "
0226                                  << labelEBRecHit_ << " EE " << labelEERecHit_ << " HBHE " << labelHBHERecHit_ << " MU "
0227                                  << labelMuon_;
0228   } else {
0229     edm::LogVerbatim("HBHEMuon") << "Labels used Trig " << hlTriggerResults_ << "\n  Vtx  "
0230                                  << edm::InputTag(modnam_, labelVtx_, procnm_) << "\n  EB   " << labelEBRecHit_
0231                                  << "\n  EE   " << labelEERecHit_ << "\n  HBHE " << labelHBHERecHit_ << "\n  MU   "
0232                                  << edm::InputTag(modnam_, labelMuon_, procnm_);
0233   }
0234 
0235   if (!fileInCorr_.empty()) {
0236     std::ifstream infile(fileInCorr_.c_str());
0237     if (infile.is_open()) {
0238       while (true) {
0239         unsigned int id;
0240         double cfac;
0241         infile >> id >> cfac;
0242         if (!infile.good())
0243           break;
0244         corrValue_[DetId(id)] = cfac;
0245       }
0246       infile.close();
0247     }
0248   }
0249   useMyCorr_ = (!corrValue_.empty());
0250   edm::LogVerbatim("HBHEMuon") << "Flags used: UseRaw " << useRaw_ << " GetCharge " << getCharge_ << " UnCorrect "
0251                                << unCorrect_ << " IgnoreHECorr " << ignoreHECorr_ << " CollapseDepth " << collapseDepth_
0252                                << ":" << mergedDepth_ << " IsItPlan1 " << isItPlan1_ << " IsItPreRecHit "
0253                                << isItPreRecHit_ << " UseMyCorr " << useMyCorr_ << " pMinMuon " << pMinMuon_
0254                                << " usePFThresh " << usePFThresh_;
0255 }
0256 
0257 //
0258 // member functions
0259 //
0260 
0261 // ------------ method called for each event  ------------
0262 void HcalHBHEMuonAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0263   ++kount_;
0264   clearVectors();
0265   std::vector<bool> muon_is_good, muon_global, muon_tracker;
0266   std::vector<bool> muon_is_tight, muon_is_medium;
0267   std::vector<double> ptGlob, etaGlob, phiGlob, energyMuon, pMuon;
0268   std::vector<float> muon_trkKink, muon_chi2LocalPosition, muon_segComp;
0269   std::vector<int> trackerLayer, numPixelLayers, tight_PixelHits;
0270   std::vector<bool> innerTrack, outerTrack, globalTrack;
0271   std::vector<double> chiTracker, dxyTracker, dzTracker;
0272   std::vector<double> innerTrackpt, innerTracketa, innerTrackphi;
0273   std::vector<double> tight_validFraction, outerTrackChi;
0274   std::vector<double> outerTrackPt, outerTrackEta, outerTrackPhi;
0275   std::vector<int> outerTrackHits, outerTrackRHits;
0276   std::vector<double> globalTrckPt, globalTrckEta, globalTrckPhi;
0277   std::vector<int> globalMuonHits, matchedStat;
0278   std::vector<double> chiGlobal, tight_LongPara, tight_TransImpara;
0279   std::vector<double> isolationR04, isolationR03;
0280   std::vector<double> ecalEnergy, hcalEnergy, hoEnergy;
0281   std::vector<bool> matchedId, hcalHot;
0282   std::vector<double> ecal3x3Energy, hcal1x1Energy;
0283   std::vector<unsigned int> ecalDetId, hcalDetId, ehcalDetId;
0284   std::vector<int> hcal_ieta, hcal_iphi;
0285   std::vector<double> hcalDepthEnergy[depthMax_];
0286   std::vector<double> hcalDepthActiveLength[depthMax_];
0287   std::vector<double> hcalDepthEnergyHot[depthMax_];
0288   std::vector<double> hcalDepthActiveLengthHot[depthMax_];
0289   std::vector<double> hcalDepthChargeHot[depthMax_];
0290   std::vector<double> hcalDepthChargeHotBG[depthMax_];
0291   std::vector<double> hcalDepthEnergyCorr[depthMax_];
0292   std::vector<double> hcalDepthEnergyHotCorr[depthMax_];
0293   std::vector<bool> hcalDepthMatch[depthMax_];
0294   std::vector<bool> hcalDepthMatchHot[depthMax_];
0295   std::vector<double> hcalActiveLength, hcalActiveLengthHot;
0296   runNumber_ = iEvent.id().run();
0297   eventNumber_ = iEvent.id().event();
0298   lumiNumber_ = iEvent.id().luminosityBlock();
0299   bxNumber_ = iEvent.bunchCrossing();
0300 #ifdef EDM_ML_DEBUG
0301   edm::LogVerbatim("HBHEMuon") << "Run " << runNumber_ << " Event " << eventNumber_ << " Lumi " << lumiNumber_ << " BX "
0302                                << bxNumber_ << std::endl;
0303 #endif
0304   const edm::Handle<edm::TriggerResults>& _Triggers = iEvent.getHandle(tok_trigRes_);
0305 #ifdef EDM_ML_DEBUG
0306   if ((verbosity_ / 10000) % 10 > 0)
0307     edm::LogVerbatim("HBHEMuon") << "Size of all triggers " << all_triggers_.size();
0308 #endif
0309   int Ntriggers = static_cast<int>(all_triggers_.size());
0310 #ifdef EDM_ML_DEBUG
0311   if ((verbosity_ / 10000) % 10 > 0)
0312     edm::LogVerbatim("HBHEMuon") << "Size of HLT MENU: " << _Triggers->size();
0313 #endif
0314   if (_Triggers.isValid()) {
0315     const edm::TriggerNames& triggerNames_ = iEvent.triggerNames(*_Triggers);
0316     std::vector<int> index;
0317     for (int i = 0; i < Ntriggers; i++) {
0318       index.push_back(triggerNames_.triggerIndex(all_triggers_[i]));
0319       int triggerSize = static_cast<int>(_Triggers->size());
0320 #ifdef EDM_ML_DEBUG
0321       if ((verbosity_ / 10000) % 10 > 0)
0322         edm::LogVerbatim("HBHEMuon") << "outside loop " << index[i] << "\ntriggerSize " << triggerSize;
0323 #endif
0324       if (index[i] < triggerSize) {
0325         hltresults_.push_back(_Triggers->accept(index[i]));
0326 #ifdef EDM_ML_DEBUG
0327         if ((verbosity_ / 10000) % 10 > 0)
0328           edm::LogVerbatim("HBHEMuon") << "Trigger_info " << triggerSize << " triggerSize " << index[i]
0329                                        << " trigger_index " << hltresults_.at(i) << " hltresult";
0330 #endif
0331       } else {
0332         if ((verbosity_ / 10000) % 10 > 0)
0333           edm::LogVerbatim("HBHEMuon") << "Requested HLT path \""
0334                                        << "\" does not exist";
0335       }
0336     }
0337   }
0338 
0339   // get handles to calogeometry and calotopology
0340   const MagneticField* bField = &iSetup.getData(tok_magField_);
0341   const EcalChannelStatus* theEcalChStatus = &iSetup.getData(tok_chan_);
0342   const EcalSeverityLevelAlgo* sevlv = &iSetup.getData(tok_sevlv_);
0343   const CaloTopology* caloTopology = &iSetup.getData(tok_topo_);
0344   const HcalDbService* conditions = &iSetup.getData(tok_dbservice_);
0345   const EcalPFRecHitThresholds* eThresholds = &iSetup.getData(tok_ecalPFRecHitThresholds_);
0346 
0347   // Relevant blocks from iEvent
0348   const edm::Handle<reco::VertexCollection>& vtx = iEvent.getHandle(tok_Vtx_);
0349 
0350   edm::Handle<EcalRecHitCollection> barrelRecHitsHandle = iEvent.getHandle(tok_EB_);
0351   edm::Handle<EcalRecHitCollection> endcapRecHitsHandle = iEvent.getHandle(tok_EE_);
0352 
0353   edm::Handle<HBHERecHitCollection> hbhe = iEvent.getHandle(tok_HBHE_);
0354 
0355   const edm::Handle<reco::MuonCollection>& _Muon = iEvent.getHandle(tok_Muon_);
0356 
0357   // require a good vertex
0358   math::XYZPoint pvx;
0359   goodVertex_ = 0;
0360   if (!vtx.isValid()) {
0361 #ifdef EDM_ML_DEBUG
0362     edm::LogVerbatim("HBHEMuon") << "No Good Vertex found == Reject";
0363 #endif
0364     return;
0365   }
0366   reco::VertexCollection::const_iterator firstGoodVertex = vtx->end();
0367   for (reco::VertexCollection::const_iterator it = vtx->begin(); it != vtx->end(); it++) {
0368     if (isGoodVertex(*it)) {
0369       if (firstGoodVertex == vtx->end())
0370         firstGoodVertex = it;
0371       ++goodVertex_;
0372     }
0373   }
0374   if (firstGoodVertex != vtx->end())
0375     pvx = firstGoodVertex->position();
0376 
0377   bool accept(false);
0378   if (_Muon.isValid() && barrelRecHitsHandle.isValid() && endcapRecHitsHandle.isValid() && hbhe.isValid()) {
0379     for (const auto& RecMuon : (*(_Muon.product()))) {
0380       muon_is_good.push_back(RecMuon.isPFMuon());
0381       muon_global.push_back(RecMuon.isGlobalMuon());
0382       muon_tracker.push_back(RecMuon.isTrackerMuon());
0383       ptGlob.push_back(RecMuon.pt());
0384       etaGlob.push_back(RecMuon.eta());
0385       phiGlob.push_back(RecMuon.phi());
0386       energyMuon.push_back(RecMuon.energy());
0387       pMuon.push_back(RecMuon.p());
0388 #ifdef EDM_ML_DEBUG
0389       edm::LogVerbatim("HBHEMuon") << "Energy:" << RecMuon.energy() << " P:" << RecMuon.p();
0390 #endif
0391       muon_is_tight.push_back(muon::isTightMuon(RecMuon, *firstGoodVertex));
0392       muon_is_medium.push_back(muon::isMediumMuon(RecMuon));
0393       muon_trkKink.push_back(RecMuon.combinedQuality().trkKink);
0394       muon_chi2LocalPosition.push_back(RecMuon.combinedQuality().chi2LocalPosition);
0395       muon_segComp.push_back(muon::segmentCompatibility(RecMuon));
0396       // acessing tracker hits info
0397       if (RecMuon.track().isNonnull()) {
0398         trackerLayer.push_back(RecMuon.track()->hitPattern().trackerLayersWithMeasurement());
0399       } else {
0400         trackerLayer.push_back(-1);
0401       }
0402       if (RecMuon.innerTrack().isNonnull()) {
0403         innerTrack.push_back(true);
0404         numPixelLayers.push_back(RecMuon.innerTrack()->hitPattern().pixelLayersWithMeasurement());
0405         chiTracker.push_back(RecMuon.innerTrack()->normalizedChi2());
0406         dxyTracker.push_back(fabs(RecMuon.innerTrack()->dxy(pvx)));
0407         dzTracker.push_back(fabs(RecMuon.innerTrack()->dz(pvx)));
0408         innerTrackpt.push_back(RecMuon.innerTrack()->pt());
0409         innerTracketa.push_back(RecMuon.innerTrack()->eta());
0410         innerTrackphi.push_back(RecMuon.innerTrack()->phi());
0411         tight_PixelHits.push_back(RecMuon.innerTrack()->hitPattern().numberOfValidPixelHits());
0412         tight_validFraction.push_back(RecMuon.innerTrack()->validFraction());
0413       } else {
0414         innerTrack.push_back(false);
0415         numPixelLayers.push_back(0);
0416         chiTracker.push_back(0);
0417         dxyTracker.push_back(0);
0418         dzTracker.push_back(0);
0419         innerTrackpt.push_back(0);
0420         innerTracketa.push_back(0);
0421         innerTrackphi.push_back(0);
0422         tight_PixelHits.push_back(0);
0423         tight_validFraction.push_back(-99);
0424       }
0425       // outer track info
0426       if (RecMuon.outerTrack().isNonnull()) {
0427         outerTrack.push_back(true);
0428         outerTrackPt.push_back(RecMuon.outerTrack()->pt());
0429         outerTrackEta.push_back(RecMuon.outerTrack()->eta());
0430         outerTrackPhi.push_back(RecMuon.outerTrack()->phi());
0431         outerTrackChi.push_back(RecMuon.outerTrack()->normalizedChi2());
0432         outerTrackHits.push_back(RecMuon.outerTrack()->numberOfValidHits());
0433         outerTrackRHits.push_back(RecMuon.outerTrack()->recHitsSize());
0434       } else {
0435         outerTrack.push_back(false);
0436         outerTrackPt.push_back(0);
0437         outerTrackEta.push_back(0);
0438         outerTrackPhi.push_back(0);
0439         outerTrackChi.push_back(0);
0440         outerTrackHits.push_back(0);
0441         outerTrackRHits.push_back(0);
0442       }
0443       // Tight Muon cuts
0444       if (RecMuon.globalTrack().isNonnull()) {
0445         globalTrack.push_back(true);
0446         chiGlobal.push_back(RecMuon.globalTrack()->normalizedChi2());
0447         globalMuonHits.push_back(RecMuon.globalTrack()->hitPattern().numberOfValidMuonHits());
0448         matchedStat.push_back(RecMuon.numberOfMatchedStations());
0449         globalTrckPt.push_back(RecMuon.globalTrack()->pt());
0450         globalTrckEta.push_back(RecMuon.globalTrack()->eta());
0451         globalTrckPhi.push_back(RecMuon.globalTrack()->phi());
0452         tight_TransImpara.push_back(fabs(RecMuon.muonBestTrack()->dxy(pvx)));
0453         tight_LongPara.push_back(fabs(RecMuon.muonBestTrack()->dz(pvx)));
0454       } else {
0455         globalTrack.push_back(false);
0456         chiGlobal.push_back(0);
0457         globalMuonHits.push_back(0);
0458         matchedStat.push_back(0);
0459         globalTrckPt.push_back(0);
0460         globalTrckEta.push_back(0);
0461         globalTrckPhi.push_back(0);
0462         tight_TransImpara.push_back(0);
0463         tight_LongPara.push_back(0);
0464       }
0465 
0466       isolationR04.push_back(
0467           ((RecMuon.pfIsolationR04().sumChargedHadronPt +
0468             std::max(0.,
0469                      RecMuon.pfIsolationR04().sumNeutralHadronEt + RecMuon.pfIsolationR04().sumPhotonEt -
0470                          (0.5 * RecMuon.pfIsolationR04().sumPUPt))) /
0471            RecMuon.pt()));
0472 
0473       isolationR03.push_back(
0474           ((RecMuon.pfIsolationR03().sumChargedHadronPt +
0475             std::max(0.,
0476                      RecMuon.pfIsolationR03().sumNeutralHadronEt + RecMuon.pfIsolationR03().sumPhotonEt -
0477                          (0.5 * RecMuon.pfIsolationR03().sumPUPt))) /
0478            RecMuon.pt()));
0479 
0480       ecalEnergy.push_back(RecMuon.calEnergy().emS9);
0481       hcalEnergy.push_back(RecMuon.calEnergy().hadS9);
0482       hoEnergy.push_back(RecMuon.calEnergy().hoS9);
0483 
0484       double eEcal(0), eHcal(0), activeLengthTot(0), activeLengthHotTot(0);
0485       double eHcalDepth[depthMax_], eHcalDepthHot[depthMax_];
0486       double eHcalDepthC[depthMax_], eHcalDepthHotC[depthMax_];
0487       double cHcalDepthHot[depthMax_], cHcalDepthHotBG[depthMax_];
0488       double activeL[depthMax_], activeHotL[depthMax_];
0489       bool matchDepth[depthMax_], matchDepthHot[depthMax_];
0490       HcalDetId eHcalDetId[depthMax_];
0491       unsigned int isHot(0);
0492       bool tmpmatch(false);
0493       int ieta(-1000), iphi(-1000);
0494       for (int i = 0; i < depthMax_; ++i) {
0495         eHcalDepth[i] = eHcalDepthHot[i] = 0;
0496         eHcalDepthC[i] = eHcalDepthHotC[i] = 0;
0497         cHcalDepthHot[i] = cHcalDepthHotBG[i] = 0;
0498         activeL[i] = activeHotL[i] = 0;
0499         matchDepth[i] = matchDepthHot[i] = true;
0500       }
0501       if (RecMuon.innerTrack().isNonnull()) {
0502         const reco::Track* pTrack = (RecMuon.innerTrack()).get();
0503         spr::propagatedTrackID trackID = spr::propagateCALO(pTrack, geo_, bField, (((verbosity_ / 100) % 10 > 0)));
0504         if ((RecMuon.p() > pMinMuon_) && (trackID.okHCAL))
0505           accept = true;
0506 
0507         ecalDetId.push_back((trackID.detIdECAL)());
0508         hcalDetId.push_back((trackID.detIdHCAL)());
0509         ehcalDetId.push_back((trackID.detIdEHCAL)());
0510 
0511         HcalDetId check;
0512         std::pair<bool, HcalDetId> info = spr::propagateHCALBack(pTrack, geo_, bField, (((verbosity_ / 100) % 10 > 0)));
0513         if (info.first) {
0514           check = info.second;
0515         }
0516 
0517         bool okE = trackID.okECAL;
0518         if (okE) {
0519           const DetId isoCell(trackID.detIdECAL);
0520           std::pair<double, bool> e3x3 = (usePFThresh_ ? spr::eECALmatrix(isoCell,
0521                                                                           barrelRecHitsHandle,
0522                                                                           endcapRecHitsHandle,
0523                                                                           *theEcalChStatus,
0524                                                                           geo_,
0525                                                                           caloTopology,
0526                                                                           sevlv,
0527                                                                           eThresholds,
0528                                                                           1,
0529                                                                           1,
0530                                                                           false)
0531                                                        : spr::eECALmatrix(isoCell,
0532                                                                           barrelRecHitsHandle,
0533                                                                           endcapRecHitsHandle,
0534                                                                           *theEcalChStatus,
0535                                                                           geo_,
0536                                                                           caloTopology,
0537                                                                           sevlv,
0538                                                                           1,
0539                                                                           1,
0540                                                                           -100.0,
0541                                                                           -100.0,
0542                                                                           -500.0,
0543                                                                           500.0,
0544                                                                           false));
0545           eEcal = e3x3.first;
0546           okE = e3x3.second;
0547         }
0548 #ifdef EDM_ML_DEBUG
0549         edm::LogVerbatim("HBHEMuon") << "Propagate Track to ECAL: " << okE << ":" << trackID.okECAL << " E " << eEcal;
0550 #endif
0551 
0552         if (trackID.okHCAL) {
0553           DetId closestCell(trackID.detIdHCAL);
0554           HcalDetId hcidt(closestCell.rawId());
0555           if ((hcidt.ieta() == check.ieta()) && (hcidt.iphi() == check.iphi()))
0556             tmpmatch = true;
0557 #ifdef EDM_ML_DEBUG
0558           edm::LogVerbatim("HBHEMuon") << "Front " << hcidt << " Back " << info.first << ":" << check << " Match "
0559                                        << tmpmatch;
0560 #endif
0561 
0562           HcalSubdetector subdet = hcidt.subdet();
0563           ieta = hcidt.ieta();
0564           iphi = hcidt.iphi();
0565           bool hborhe = (std::abs(ieta) == 16);
0566 
0567           eHcal = spr::eHCALmatrix(theHBHETopology_,
0568                                    closestCell,
0569                                    hbhe,
0570                                    0,
0571                                    0,
0572                                    false,
0573                                    true,
0574                                    -100.0,
0575                                    -100.0,
0576                                    -100.0,
0577                                    -100.0,
0578                                    -500.,
0579                                    500.,
0580                                    useRaw_);
0581           std::vector<std::pair<double, int>> ehdepth;
0582           spr::energyHCALCell((HcalDetId)closestCell,
0583                               hbhe,
0584                               ehdepth,
0585                               maxDepth_,
0586                               -100.0,
0587                               -100.0,
0588                               -100.0,
0589                               -100.0,
0590                               -500.0,
0591                               500.0,
0592                               useRaw_,
0593                               depth16HE(ieta, iphi),
0594                               (((verbosity_ / 1000) % 10) > 0));
0595           for (int i = 0; i < depthMax_; ++i)
0596             eHcalDetId[i] = HcalDetId();
0597           for (unsigned int i = 0; i < ehdepth.size(); ++i) {
0598             HcalSubdetector subdet0 =
0599                 (hborhe) ? ((ehdepth[i].second >= depth16HE(ieta, iphi)) ? HcalEndcap : HcalBarrel) : subdet;
0600             HcalDetId hcid0(subdet0, ieta, iphi, ehdepth[i].second);
0601             double actL = activeLength(DetId(hcid0), ((verbosity_ / 100000) % 10));
0602             double ene = ehdepth[i].first;
0603 #ifdef EDM_ML_DEBUG
0604             if ((verbosity_ / 100000) % 10)
0605               edm::LogVerbatim("HBHEMuon")
0606                   << "DetId " << subdet0 << ":" << ieta << ":" << iphi << ":" << ehdepth[i].second << " E " << ene;
0607 #endif
0608             bool tmpC(false);
0609             if (ene > 0.0) {
0610               if (!(theHBHETopology_->validHcal(hcid0))) {
0611                 edm::LogWarning("HBHEMuon") << "(1) Invalid ID " << hcid0 << " with E = " << ene;
0612                 edm::LogWarning("HBHEMuon") << HcalDetId(closestCell) << " with " << ehdepth.size() << " depths:";
0613                 for (const auto& ehd : ehdepth)
0614                   edm::LogWarning("HBHEMuon") << " " << ehd.second << ":" << ehd.first;
0615               } else {
0616                 tmpC = goodCell(hcid0, pTrack, geo_, bField);
0617                 double enec(ene);
0618                 if (unCorrect_) {
0619                   double corr = (ignoreHECorr_ && (subdet0 == HcalEndcap)) ? 1.0 : respCorr(DetId(hcid0));
0620                   if (corr != 0)
0621                     ene /= corr;
0622 #ifdef EDM_ML_DEBUG
0623                   HcalDetId id = (isItPlan1_ && isItPreRecHit_) ? hdc_->mergedDepthDetId(hcid0) : hcid0;
0624                   edm::LogVerbatim("HBHEMuon") << hcid0 << ":" << id << " Corr " << corr;
0625 #endif
0626                 }
0627                 int depth = ehdepth[i].second - 1;
0628                 if (collapseDepth_) {
0629                   HcalDetId id = hdc_->mergedDepthDetId(hcid0);
0630                   depth = id.depth() - 1;
0631                 }
0632                 eHcalDepth[depth] += ene;
0633                 eHcalDepthC[depth] += enec;
0634                 activeL[depth] += actL;
0635                 activeLengthTot += actL;
0636                 matchDepth[depth] = (matchDepth[depth] && tmpC);
0637 #ifdef EDM_ML_DEBUG
0638                 if ((verbosity_ % 10) > 0)
0639                   edm::LogVerbatim("HBHEMuon")
0640                       << hcid0 << " E " << ene << ":" << enec << " L " << actL << " Match " << tmpC;
0641 #endif
0642               }
0643             }
0644           }
0645 #ifdef EDM_ML_DEBUG
0646           if ((verbosity_ % 10) > 0) {
0647             edm::LogVerbatim("HBHEMuon") << hcidt << " Match " << tmpmatch << " Depths " << ehdepth.size();
0648             for (unsigned int k = 0; k < ehdepth.size(); ++k)
0649               edm::LogVerbatim("HBHEMuon") << " [" << k << ":" << ehdepth[k].second << "] Match " << matchDepth[k]
0650                                            << " E " << eHcalDepth[k] << ":" << eHcalDepthC[k] << " L " << activeL[k];
0651           }
0652 #endif
0653           HcalDetId hotCell;
0654           spr::eHCALmatrix(geo_, theHBHETopology_, closestCell, hbhe, 1, 1, hotCell, false, useRaw_, false);
0655           isHot = matchId(closestCell, hotCell);
0656           if (hotCell != HcalDetId()) {
0657             subdet = HcalDetId(hotCell).subdet();
0658             ieta = HcalDetId(hotCell).ieta();
0659             iphi = HcalDetId(hotCell).iphi();
0660             hborhe = (std::abs(ieta) == 16);
0661             std::vector<std::pair<double, int>> ehdepth;
0662             spr::energyHCALCell(hotCell,
0663                                 hbhe,
0664                                 ehdepth,
0665                                 maxDepth_,
0666                                 -100.0,
0667                                 -100.0,
0668                                 -100.0,
0669                                 -100.0,
0670                                 -500.0,
0671                                 500.0,
0672                                 useRaw_,
0673                                 depth16HE(ieta, iphi),
0674                                 false);  //(((verbosity_/1000)%10)>0    ));
0675             for (int i = 0; i < depthMax_; ++i)
0676               eHcalDetId[i] = HcalDetId();
0677             for (unsigned int i = 0; i < ehdepth.size(); ++i) {
0678               HcalSubdetector subdet0 =
0679                   (hborhe) ? ((ehdepth[i].second >= depth16HE(ieta, iphi)) ? HcalEndcap : HcalBarrel) : subdet;
0680               HcalDetId hcid0(subdet0, ieta, iphi, ehdepth[i].second);
0681               double actL = activeLength(DetId(hcid0), ((verbosity_ / 100000) % 10));
0682               double ene = ehdepth[i].first;
0683 #ifdef EDM_ML_DEBUG
0684               if ((verbosity_ / 100000) % 10)
0685                 edm::LogVerbatim("HBHEMuon")
0686                     << "DetId " << subdet0 << ":" << ieta << ":" << iphi << ":" << ehdepth[i].second << " E " << ene;
0687 #endif
0688               bool tmpC(false);
0689               if (ene > 0.0) {
0690                 if (!(theHBHETopology_->validHcal(hcid0))) {
0691                   edm::LogWarning("HBHEMuon") << "(2) Invalid ID " << hcid0 << " with E = " << ene;
0692                   edm::LogWarning("HBHEMuon") << HcalDetId(hotCell) << " with " << ehdepth.size() << " depths:";
0693                   for (const auto& ehd : ehdepth)
0694                     edm::LogWarning("HBHEMuon") << " " << ehd.second << ":" << ehd.first;
0695                 } else {
0696                   tmpC = goodCell(hcid0, pTrack, geo_, bField);
0697                   double chg(ene), enec(ene);
0698                   if (unCorrect_) {
0699                     double corr = (ignoreHECorr_ && (subdet0 == HcalEndcap)) ? 1.0 : respCorr(DetId(hcid0));
0700                     if (corr != 0)
0701                       ene /= corr;
0702 #ifdef EDM_ML_DEBUG
0703                     HcalDetId id = (isItPlan1_ && isItPreRecHit_) ? hdc_->mergedDepthDetId(hcid0) : hcid0;
0704                     edm::LogVerbatim("HBHEMuon")
0705                         << hcid0 << ":" << id << " Corr " << corr << " E " << ene << ":" << enec;
0706 #endif
0707                   }
0708                   if (getCharge_) {
0709                     double gain = gainFactor(conditions, hcid0);
0710                     if (gain != 0)
0711                       chg /= gain;
0712 #ifdef EDM_ML_DEBUG
0713                     edm::LogVerbatim("HBHEMuon") << hcid0 << " Gain " << gain << " C " << chg;
0714 #endif
0715                   }
0716                   int depth = ehdepth[i].second - 1;
0717                   if (collapseDepth_) {
0718                     HcalDetId id = hdc_->mergedDepthDetId(hcid0);
0719                     depth = id.depth() - 1;
0720                   }
0721                   eHcalDepthHot[depth] += ene;
0722                   eHcalDepthHotC[depth] += enec;
0723                   cHcalDepthHot[depth] += chg;
0724                   activeHotL[depth] += actL;
0725                   activeLengthHotTot += actL;
0726                   matchDepthHot[depth] = (matchDepthHot[depth] && tmpC);
0727 #ifdef EDM_ML_DEBUG
0728                   if ((verbosity_ % 10) > 0)
0729                     edm::LogVerbatim("HBHEMuon") << hcid0 << " depth " << depth << " E " << ene << ":" << enec << " C "
0730                                                  << chg << " L " << actL << " Match " << tmpC;
0731 #endif
0732                 }
0733               }
0734             }
0735 
0736             HcalDetId oppCell(subdet, -ieta, iphi, HcalDetId(hotCell).depth());
0737             std::vector<std::pair<double, int>> ehdeptho;
0738             spr::energyHCALCell(oppCell,
0739                                 hbhe,
0740                                 ehdeptho,
0741                                 maxDepth_,
0742                                 -100.0,
0743                                 -100.0,
0744                                 -100.0,
0745                                 -100.0,
0746                                 -500.0,
0747                                 500.0,
0748                                 useRaw_,
0749                                 depth16HE(-ieta, iphi),
0750                                 false);  //(((verbosity_/1000)%10)>0));
0751             for (unsigned int i = 0; i < ehdeptho.size(); ++i) {
0752               HcalSubdetector subdet0 =
0753                   (hborhe) ? ((ehdeptho[i].second >= depth16HE(-ieta, iphi)) ? HcalEndcap : HcalBarrel) : subdet;
0754               HcalDetId hcid0(subdet0, -ieta, iphi, ehdeptho[i].second);
0755               double ene = ehdeptho[i].first;
0756               if (ene > 0.0) {
0757                 if (!(theHBHETopology_->validHcal(hcid0))) {
0758                   edm::LogWarning("HBHEMuon") << "(3) Invalid ID " << hcid0 << " with E = " << ene;
0759                   edm::LogWarning("HBHEMuon") << oppCell << " with " << ehdeptho.size() << " depths:";
0760                   for (const auto& ehd : ehdeptho)
0761                     edm::LogWarning("HBHEMuon") << " " << ehd.second << ":" << ehd.first;
0762                 } else {
0763                   double chg(ene);
0764                   if (unCorrect_) {
0765                     double corr = (ignoreHECorr_ && (subdet0 == HcalEndcap)) ? 1.0 : respCorr(DetId(hcid0));
0766                     if (corr != 0)
0767                       ene /= corr;
0768 #ifdef EDM_ML_DEBUG
0769                     HcalDetId id = (isItPlan1_ && isItPreRecHit_) ? hdc_->mergedDepthDetId(hcid0) : hcid0;
0770                     edm::LogVerbatim("HBHEMuon")
0771                         << hcid0 << ":" << id << " Corr " << corr << " E " << ene << ":" << ehdeptho[i].first;
0772 #endif
0773                   }
0774                   if (getCharge_) {
0775                     double gain = gainFactor(conditions, hcid0);
0776                     if (gain != 0)
0777                       chg /= gain;
0778 #ifdef EDM_ML_DEBUG
0779                     edm::LogVerbatim("HBHEMuon") << hcid0 << " Gain " << gain << " C " << chg;
0780 #endif
0781                   }
0782                   int depth = ehdeptho[i].second - 1;
0783                   if (collapseDepth_) {
0784                     HcalDetId id = hdc_->mergedDepthDetId(hcid0);
0785                     depth = id.depth() - 1;
0786                   }
0787                   cHcalDepthHotBG[depth] += chg;
0788 #ifdef EDM_ML_DEBUG
0789                   if ((verbosity_ % 10) > 0)
0790                     edm::LogVerbatim("HBHEMuon") << hcid0 << " Depth " << depth << " E " << ene << " C " << chg;
0791 #endif
0792                 }
0793               }
0794             }
0795           }
0796         }
0797 #ifdef EDM_ML_DEBUG
0798         edm::LogVerbatim("HBHEMuon") << "Propagate Track to HCAL: " << trackID.okHCAL << " Match " << tmpmatch
0799                                      << " Hot " << isHot << " Energy " << eHcal << std::endl;
0800 #endif
0801 
0802       } else {
0803         ecalDetId.push_back(0);
0804         hcalDetId.push_back(0);
0805         ehcalDetId.push_back(0);
0806       }
0807 
0808       matchedId.push_back(tmpmatch);
0809       ecal3x3Energy.push_back(eEcal);
0810       hcal1x1Energy.push_back(eHcal);
0811       hcal_ieta.push_back(ieta);
0812       hcal_iphi.push_back(iphi);
0813       for (int i = 0; i < depthMax_; ++i) {
0814         hcalDepthEnergy[i].push_back(eHcalDepth[i]);
0815         hcalDepthActiveLength[i].push_back(activeL[i]);
0816         hcalDepthEnergyHot[i].push_back(eHcalDepthHot[i]);
0817         hcalDepthActiveLengthHot[i].push_back(activeHotL[i]);
0818         hcalDepthEnergyCorr[i].push_back(eHcalDepthC[i]);
0819         hcalDepthEnergyHotCorr[i].push_back(eHcalDepthHotC[i]);
0820         hcalDepthChargeHot[i].push_back(cHcalDepthHot[i]);
0821         hcalDepthChargeHotBG[i].push_back(cHcalDepthHotBG[i]);
0822         hcalDepthMatch[i].push_back(matchDepth[i]);
0823         hcalDepthMatchHot[i].push_back(matchDepthHot[i]);
0824       }
0825       hcalActiveLength.push_back(activeLengthTot);
0826       hcalHot.push_back(isHot);
0827       hcalActiveLengthHot.push_back(activeLengthHotTot);
0828     }
0829   }
0830   if (accept) {
0831 #ifdef EDM_ML_DEBUG
0832     for (unsigned int i = 0; i < hcal_ieta.size(); ++i)
0833       edm::LogVerbatim("HBHEMuon") << "[" << i << "] ieta/iphi for entry to "
0834                                    << "HCAL has value of " << hcal_ieta[i] << ":" << hcal_iphi[i];
0835 #endif
0836     for (unsigned int k = 0; k < muon_is_good.size(); ++k) {
0837       muon_is_good_ = muon_is_good[k];
0838       muon_global_ = muon_global[k];
0839       muon_tracker_ = muon_tracker[k];
0840       muon_is_tight_ = muon_is_tight[k];
0841       muon_is_medium_ = muon_is_medium[k];
0842       ptGlob_ = ptGlob[k];
0843       etaGlob_ = etaGlob[k];
0844       phiGlob_ = phiGlob[k];
0845       energyMuon_ = energyMuon[k];
0846       pMuon_ = pMuon[k];
0847       muon_trkKink_ = muon_trkKink[k];
0848       muon_chi2LocalPosition_ = muon_chi2LocalPosition[k];
0849       muon_segComp_ = muon_segComp[k];
0850       trackerLayer_ = trackerLayer[k];
0851       numPixelLayers_ = numPixelLayers[k];
0852       tight_PixelHits_ = tight_PixelHits[k];
0853       innerTrack_ = innerTrack[k];
0854       outerTrack_ = outerTrack[k];
0855       globalTrack_ = globalTrack[k];
0856       chiTracker_ = chiTracker[k];
0857       dxyTracker_ = dxyTracker[k];
0858       dzTracker_ = dzTracker[k];
0859       innerTrackpt_ = innerTrackpt[k];
0860       innerTracketa_ = innerTracketa[k];
0861       innerTrackphi_ = innerTrackphi[k];
0862       tight_validFraction_ = tight_validFraction[k];
0863       outerTrackChi_ = outerTrackChi[k];
0864       outerTrackPt_ = outerTrackPt[k];
0865       outerTrackEta_ = outerTrackEta[k];
0866       outerTrackPhi_ = outerTrackPhi[k];
0867       outerTrackHits_ = outerTrackHits[k];
0868       outerTrackRHits_ = outerTrackRHits[k];
0869       globalTrckPt_ = globalTrckPt[k];
0870       globalTrckEta_ = globalTrckEta[k];
0871       globalTrckPhi_ = globalTrckPhi[k];
0872       globalMuonHits_ = globalMuonHits[k];
0873       matchedStat_ = matchedStat[k];
0874       chiGlobal_ = chiGlobal[k];
0875       tight_LongPara_ = tight_LongPara[k];
0876       tight_TransImpara_ = tight_TransImpara[k];
0877       isolationR04_ = isolationR04[k];
0878       isolationR03_ = isolationR03[k];
0879       ecalEnergy_ = ecalEnergy[k];
0880       hcalEnergy_ = hcalEnergy[k];
0881       hoEnergy_ = hoEnergy[k];
0882       matchedId_ = matchedId[k];
0883       hcalHot_ = hcalHot[k];
0884       ecal3x3Energy_ = ecal3x3Energy[k];
0885       hcal1x1Energy_ = hcal1x1Energy[k];
0886       ecalDetId_ = ecalDetId[k];
0887       hcalDetId_ = hcalDetId[k];
0888       ehcalDetId_ = ehcalDetId[k];
0889       hcal_ieta_ = hcal_ieta[k];
0890       hcal_iphi_ = hcal_iphi[k];
0891       for (int i = 0; i < depthMax_; ++i) {
0892         hcalDepthEnergy_[i] = hcalDepthEnergy[i][k];
0893         hcalDepthActiveLength_[i] = hcalDepthActiveLength[i][k];
0894         hcalDepthEnergyHot_[i] = hcalDepthEnergyHot[i][k];
0895         hcalDepthActiveLengthHot_[i] = hcalDepthActiveLengthHot[i][k];
0896         hcalDepthChargeHot_[i] = hcalDepthChargeHot[i][k];
0897         hcalDepthChargeHotBG_[i] = hcalDepthChargeHotBG[i][k];
0898         hcalDepthEnergyCorr_[i] = hcalDepthEnergyCorr[i][k];
0899         hcalDepthEnergyHotCorr_[i] = hcalDepthEnergyHotCorr[i][k];
0900         hcalDepthMatch_[i] = hcalDepthMatch[i][k];
0901         hcalDepthMatchHot_[i] = hcalDepthMatchHot[i][k];
0902       }
0903       hcalActiveLength_ = hcalActiveLength[k];
0904       hcalActiveLengthHot_ = hcalActiveLengthHot[k];
0905       tree_->Fill();
0906     }
0907   }
0908 }
0909 
0910 // ------------ method called once each job just before starting event loop  ------------
0911 void HcalHBHEMuonAnalyzer::beginJob() {
0912   edm::Service<TFileService> fs;
0913   tree_ = fs->make<TTree>("TREE", "TREE");
0914   tree_->Branch("Event_No", &eventNumber_);
0915   tree_->Branch("Run_No", &runNumber_);
0916   tree_->Branch("LumiNumber", &lumiNumber_);
0917   tree_->Branch("BXNumber", &bxNumber_);
0918   tree_->Branch("GoodVertex", &goodVertex_);
0919   tree_->Branch("PF_Muon", &muon_is_good_);
0920   tree_->Branch("Global_Muon", &muon_global_);
0921   tree_->Branch("Tracker_muon", &muon_tracker_);
0922   tree_->Branch("MuonIsTight", &muon_is_tight_);
0923   tree_->Branch("MuonIsMedium", &muon_is_medium_);
0924   tree_->Branch("pt_of_muon", &ptGlob_);
0925   tree_->Branch("eta_of_muon", &etaGlob_);
0926   tree_->Branch("phi_of_muon", &phiGlob_);
0927   tree_->Branch("energy_of_muon", &energyMuon_);
0928   tree_->Branch("p_of_muon", &pMuon_);
0929   tree_->Branch("muon_trkKink", &muon_trkKink_);
0930   tree_->Branch("muon_chi2LocalPosition", &muon_chi2LocalPosition_);
0931   tree_->Branch("muon_segComp", &muon_segComp_);
0932 
0933   tree_->Branch("TrackerLayer", &trackerLayer_);
0934   tree_->Branch("NumPixelLayers", &numPixelLayers_);
0935   tree_->Branch("InnerTrackPixelHits", &tight_PixelHits_);
0936   tree_->Branch("innerTrack", &innerTrack_);
0937   tree_->Branch("chiTracker", &chiTracker_);
0938   tree_->Branch("DxyTracker", &dxyTracker_);
0939   tree_->Branch("DzTracker", &dzTracker_);
0940   tree_->Branch("innerTrackpt", &innerTrackpt_);
0941   tree_->Branch("innerTracketa", &innerTracketa_);
0942   tree_->Branch("innerTrackphi", &innerTrackphi_);
0943   tree_->Branch("tight_validFraction", &tight_validFraction_);
0944 
0945   tree_->Branch("OuterTrack", &outerTrack_);
0946   tree_->Branch("OuterTrackChi", &outerTrackChi_);
0947   tree_->Branch("OuterTrackPt", &outerTrackPt_);
0948   tree_->Branch("OuterTrackEta", &outerTrackEta_);
0949   tree_->Branch("OuterTrackPhi", &outerTrackPhi_);
0950   tree_->Branch("OuterTrackHits", &outerTrackHits_);
0951   tree_->Branch("OuterTrackRHits", &outerTrackRHits_);
0952 
0953   tree_->Branch("GlobalTrack", &globalTrack_);
0954   tree_->Branch("GlobalTrckPt", &globalTrckPt_);
0955   tree_->Branch("GlobalTrckEta", &globalTrckEta_);
0956   tree_->Branch("GlobalTrckPhi", &globalTrckPhi_);
0957   tree_->Branch("Global_Muon_Hits", &globalMuonHits_);
0958   tree_->Branch("MatchedStations", &matchedStat_);
0959   tree_->Branch("GlobTrack_Chi", &chiGlobal_);
0960   tree_->Branch("Tight_LongitudinalImpactparameter", &tight_LongPara_);
0961   tree_->Branch("Tight_TransImpactparameter", &tight_TransImpara_);
0962 
0963   tree_->Branch("IsolationR04", &isolationR04_);
0964   tree_->Branch("IsolationR03", &isolationR03_);
0965   tree_->Branch("ecal_3into3", &ecalEnergy_);
0966   tree_->Branch("hcal_3into3", &hcalEnergy_);
0967   tree_->Branch("tracker_3into3", &hoEnergy_);
0968 
0969   tree_->Branch("matchedId", &matchedId_);
0970   tree_->Branch("hcal_cellHot", &hcalHot_);
0971 
0972   tree_->Branch("ecal_3x3", &ecal3x3Energy_);
0973   tree_->Branch("hcal_1x1", &hcal1x1Energy_);
0974   tree_->Branch("ecal_detID", &ecalDetId_);
0975   tree_->Branch("hcal_detID", &hcalDetId_);
0976   tree_->Branch("ehcal_detID", &ehcalDetId_);
0977   tree_->Branch("hcal_ieta", &hcal_ieta_);
0978   tree_->Branch("hcal_iphi", &hcal_iphi_);
0979 
0980   char name[100];
0981   for (int k = 0; k < maxDepth_; ++k) {
0982     sprintf(name, "hcal_edepth%d", (k + 1));
0983     tree_->Branch(name, &hcalDepthEnergy_[k]);
0984     sprintf(name, "hcal_activeL%d", (k + 1));
0985     tree_->Branch(name, &hcalDepthActiveLength_[k]);
0986     sprintf(name, "hcal_edepthHot%d", (k + 1));
0987     tree_->Branch(name, &hcalDepthEnergyHot_[k]);
0988     sprintf(name, "hcal_activeHotL%d", (k + 1));
0989     tree_->Branch(name, &hcalDepthActiveLengthHot_[k]);
0990     sprintf(name, "hcal_cdepthHot%d", (k + 1));
0991     tree_->Branch(name, &hcalDepthChargeHot_[k]);
0992     sprintf(name, "hcal_cdepthHotBG%d", (k + 1));
0993     tree_->Branch(name, &hcalDepthChargeHotBG_[k]);
0994     sprintf(name, "hcal_edepthCorrect%d", (k + 1));
0995     tree_->Branch(name, &hcalDepthEnergyCorr_[k]);
0996     sprintf(name, "hcal_edepthHotCorrect%d", (k + 1));
0997     tree_->Branch(name, &hcalDepthEnergyHotCorr_[k]);
0998     sprintf(name, "hcal_depthMatch%d", (k + 1));
0999     tree_->Branch(name, &hcalDepthMatch_[k]);
1000     sprintf(name, "hcal_depthMatchHot%d", (k + 1));
1001     tree_->Branch(name, &hcalDepthMatchHot_[k]);
1002   }
1003 
1004   tree_->Branch("activeLength", &hcalActiveLength_);
1005   tree_->Branch("activeLengthHot", &hcalActiveLengthHot_);
1006 
1007   tree_->Branch("hltresults", &hltresults_);
1008   tree_->Branch("all_triggers", &all_triggers_);
1009 }
1010 
1011 // ------------ method called when starting to processes a run  ------------
1012 void HcalHBHEMuonAnalyzer::beginRun(edm::Run const& iRun, edm::EventSetup const& iSetup) {
1013   hdc_ = &iSetup.getData(tok_ddrec_);
1014   actHB.clear();
1015   actHE.clear();
1016   actHB = hdc_->getThickActive(0);
1017   actHE = hdc_->getThickActive(1);
1018 #ifdef EDM_ML_DEBUG
1019   unsigned int k1(0), k2(0);
1020   edm::LogVerbatim("HBHEMuon") << actHB.size() << " Active Length for HB";
1021   for (const auto& act : actHB) {
1022     edm::LogVerbatim("HBHEMuon") << "[" << k1 << "] ieta " << act.ieta << " depth " << act.depth << " zside "
1023                                  << act.zside << " type " << act.stype << " phi " << act.iphis.size() << ":"
1024                                  << act.iphis[0] << " L " << act.thick;
1025     HcalDetId hcid1(HcalBarrel, (act.ieta) * (act.zside), act.iphis[0], act.depth);
1026     HcalDetId hcid2 = mergedDepth_ ? hdc_->mergedDepthDetId(hcid1) : hcid1;
1027     edm::LogVerbatim("HBHEMuon") << hcid1 << " | " << hcid2 << " L "
1028                                  << activeLength(DetId(hcid2), ((verbosity_ / 100000) % 10));
1029     ++k1;
1030   }
1031   edm::LogVerbatim("HBHEMuon") << actHE.size() << " Active Length for HE";
1032   for (const auto& act : actHE) {
1033     edm::LogVerbatim("HBHEMuon") << "[" << k2 << "] ieta " << act.ieta << " depth " << act.depth << " zside "
1034                                  << act.zside << " type " << act.stype << " phi " << act.iphis.size() << ":"
1035                                  << act.iphis[0] << " L " << act.thick;
1036     HcalDetId hcid1(HcalEndcap, (act.ieta) * (act.zside), act.iphis[0], act.depth);
1037     HcalDetId hcid2 = mergedDepth_ ? hdc_->mergedDepthDetId(hcid1) : hcid1;
1038     edm::LogVerbatim("HBHEMuon") << hcid1 << " | " << hcid2 << " L "
1039                                  << activeLength(DetId(hcid2), ((verbosity_ / 100000) % 10));
1040     ++k2;
1041   }
1042 #endif
1043 
1044   bool changed = true;
1045   all_triggers_.clear();
1046   if (hltConfig_.init(iRun, iSetup, "HLT", changed)) {
1047     // if init returns TRUE, initialisation has succeeded!
1048 #ifdef EDM_ML_DEBUG
1049     edm::LogVerbatim("HBHEMuon") << "HLT config with process name HLT successfully extracted";
1050 #endif
1051     unsigned int ntriggers = hltConfig_.size();
1052     for (unsigned int t = 0; t < ntriggers; ++t) {
1053       std::string hltname(hltConfig_.triggerName(t));
1054       for (unsigned int ik = 0; ik < triggers_.size(); ++ik) {
1055         if (hltname.find(triggers_[ik]) != std::string::npos) {
1056           all_triggers_.push_back(hltname);
1057           break;
1058         }
1059       }
1060     }  //loop over ntriggers
1061     edm::LogVerbatim("HBHEMuon") << "All triggers size in begin run " << all_triggers_.size();
1062   } else {
1063     edm::LogError("HBHEMuon") << "Error! HLT config extraction with process name HLT failed";
1064   }
1065 
1066   theHBHETopology_ = &iSetup.getData(tok_htopo_);
1067   const HcalRespCorrs* resp = &iSetup.getData(tok_respcorr_);
1068   respCorrs_ = new HcalRespCorrs(*resp);
1069   respCorrs_->setTopo(theHBHETopology_);
1070   geo_ = &iSetup.getData(tok_geom_);
1071 
1072   // Write correction factors for all HB/HE events
1073   if (writeRespCorr_) {
1074     const HcalGeometry* gHcal = static_cast<const HcalGeometry*>(geo_->getSubdetectorGeometry(DetId::Hcal, HcalBarrel));
1075     const std::vector<DetId>& ids = gHcal->getValidDetIds(DetId::Hcal, 0);
1076     edm::LogVerbatim("HBHEMuon") << "\nTable of Correction Factors for Run " << iRun.run() << "\n";
1077     for (auto const& id : ids) {
1078       if ((id.det() == DetId::Hcal) && ((id.subdetId() == HcalBarrel) || (id.subdetId() == HcalEndcap))) {
1079         edm::LogVerbatim("HBHEMuon") << HcalDetId(id) << " " << id.rawId() << " "
1080                                      << (respCorrs_->getValues(id))->getValue();
1081       }
1082     }
1083   }
1084 }
1085 
1086 // ------------ method fills 'descriptions' with the allowed parameters for the module  ------------
1087 void HcalHBHEMuonAnalyzer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
1088   edm::ParameterSetDescription desc;
1089   desc.add<edm::InputTag>("hlTriggerResults", edm::InputTag("TriggerResults", "", "HLT"));
1090   desc.add<edm::InputTag>("labelEBRecHit", edm::InputTag("ecalRecHit", "EcalRecHitsEB"));
1091   desc.add<edm::InputTag>("labelEERecHit", edm::InputTag("ecalRecHit", "EcalRecHitsEE"));
1092   desc.add<edm::InputTag>("labelHBHERecHit", edm::InputTag("hbhereco"));
1093   desc.add<std::string>("labelVertex", "offlinePrimaryVertices");
1094   desc.add<std::string>("labelMuon", "muons");
1095   std::vector<std::string> trig = {};
1096   desc.add<std::vector<std::string>>("triggers", trig);
1097   desc.add<double>("pMinMuon", 5.0);
1098   desc.addUntracked<int>("verbosity", 0);
1099   desc.add<int>("useRaw", 0);
1100   desc.add<bool>("unCorrect", true);
1101   desc.add<bool>("getCharge", true);
1102   desc.add<bool>("collapseDepth", false);
1103   desc.add<bool>("isItPlan1", false);
1104   desc.addUntracked<bool>("ignoreHECorr", false);
1105   desc.addUntracked<bool>("isItPreRecHit", false);
1106   desc.addUntracked<std::string>("moduleName", "");
1107   desc.addUntracked<std::string>("processName", "");
1108   desc.addUntracked<int>("maxDepth", 7);
1109   desc.addUntracked<std::string>("fileInCorr", "");
1110   desc.addUntracked<bool>("writeRespCorr", false);
1111   desc.add<bool>("usePFThreshold", true);
1112   descriptions.add("hcalHBHEMuon", desc);
1113 }
1114 
1115 void HcalHBHEMuonAnalyzer::clearVectors() {
1116   ///clearing vectots
1117   eventNumber_ = -99999;
1118   runNumber_ = -99999;
1119   lumiNumber_ = -99999;
1120   bxNumber_ = -99999;
1121   goodVertex_ = -99999;
1122 
1123   muon_is_good_ = false;
1124   muon_global_ = false;
1125   muon_tracker_ = false;
1126   ptGlob_ = 0;
1127   etaGlob_ = 0;
1128   phiGlob_ = 0;
1129   energyMuon_ = 0;
1130   pMuon_ = 0;
1131   muon_trkKink_ = 0;
1132   muon_chi2LocalPosition_ = 0;
1133   muon_segComp_ = 0;
1134   muon_is_tight_ = false;
1135   muon_is_medium_ = false;
1136 
1137   trackerLayer_ = 0;
1138   numPixelLayers_ = 0;
1139   tight_PixelHits_ = 0;
1140   innerTrack_ = false;
1141   chiTracker_ = 0;
1142   dxyTracker_ = 0;
1143   dzTracker_ = 0;
1144   innerTrackpt_ = 0;
1145   innerTracketa_ = 0;
1146   innerTrackphi_ = 0;
1147   tight_validFraction_ = 0;
1148 
1149   outerTrack_ = false;
1150   outerTrackPt_ = 0;
1151   outerTrackEta_ = 0;
1152   outerTrackPhi_ = 0;
1153   outerTrackHits_ = 0;
1154   outerTrackRHits_ = 0;
1155   outerTrackChi_ = 0;
1156 
1157   globalTrack_ = false;
1158   globalTrckPt_ = 0;
1159   globalTrckEta_ = 0;
1160   globalTrckPhi_ = 0;
1161   globalMuonHits_ = 0;
1162   matchedStat_ = 0;
1163   chiGlobal_ = 0;
1164   tight_LongPara_ = 0;
1165   tight_TransImpara_ = 0;
1166 
1167   isolationR04_ = 0;
1168   isolationR03_ = 0;
1169   ecalEnergy_ = 0;
1170   hcalEnergy_ = 0;
1171   hoEnergy_ = 0;
1172   matchedId_ = false;
1173   hcalHot_ = false;
1174   ecal3x3Energy_ = 0;
1175   hcal1x1Energy_ = 0;
1176   ecalDetId_ = 0;
1177   hcalDetId_ = 0;
1178   ehcalDetId_ = 0;
1179   hcal_ieta_ = 0;
1180   hcal_iphi_ = 0;
1181   for (int i = 0; i < maxDepth_; ++i) {
1182     hcalDepthEnergy_[i] = 0;
1183     hcalDepthActiveLength_[i] = 0;
1184     hcalDepthEnergyHot_[i] = 0;
1185     hcalDepthActiveLengthHot_[i] = 0;
1186     hcalDepthChargeHot_[i] = 0;
1187     hcalDepthChargeHotBG_[i] = 0;
1188     hcalDepthEnergyCorr_[i] = 0;
1189     hcalDepthEnergyHotCorr_[i] = 0;
1190     hcalDepthMatch_[i] = false;
1191     hcalDepthMatchHot_[i] = false;
1192   }
1193   hcalActiveLength_ = 0;
1194   hcalActiveLengthHot_ = 0;
1195   hltresults_.clear();
1196 }
1197 
1198 int HcalHBHEMuonAnalyzer::matchId(const HcalDetId& id1, const HcalDetId& id2) {
1199   HcalDetId kd1(id1.subdet(), id1.ieta(), id1.iphi(), 1);
1200   HcalDetId kd2(id1.subdet(), id2.ieta(), id2.iphi(), 1);
1201   int match = ((kd1 == kd2) ? 1 : 0);
1202   return match;
1203 }
1204 
1205 double HcalHBHEMuonAnalyzer::activeLength(const DetId& id_, bool debug) {
1206   HcalDetId id(id_);
1207   int ieta = id.ietaAbs();
1208   int zside = id.zside();
1209   int iphi = id.iphi();
1210   std::vector<int> dpths;
1211   if (mergedDepth_) {
1212     std::vector<HcalDetId> ids;
1213     hdc_->unmergeDepthDetId(id, ids);
1214     for (auto idh : ids)
1215       dpths.emplace_back(idh.depth());
1216   } else {
1217     dpths.emplace_back(id.depth());
1218   }
1219   double lx(0);
1220   if (debug) {
1221     std::ostringstream st1;
1222     st1 << "Subdet:zside:ieta:iphi:ndepth " << id.subdet() << " : " << zside << " : " << ieta << " : " << iphi << " : "
1223         << dpths.size() << " depths";
1224     for (unsigned int k = 0; k < dpths.size(); ++k)
1225       st1 << " : " << dpths[k];
1226     edm::LogVerbatim("HBHEMuon") << st1.str();
1227   }
1228   if (id.subdet() == HcalBarrel) {
1229     for (unsigned int i = 0; i < actHB.size(); ++i) {
1230       if ((ieta == actHB[i].ieta) && (zside == actHB[i].zside) &&
1231           (std::find(dpths.begin(), dpths.end(), actHB[i].depth) != dpths.end()) &&
1232           (std::find(actHB[i].iphis.begin(), actHB[i].iphis.end(), iphi) != actHB[i].iphis.end())) {
1233         lx += actHB[i].thick;
1234         if (debug) {
1235           edm::LogVerbatim("HBHEMuon") << "Eta: " << (ieta == actHB[i].ieta) << " zside " << (zside == actHB[i].zside)
1236                                        << " Depth " << actHB[i].depth << ":"
1237                                        << (std::find(dpths.begin(), dpths.end(), actHB[i].depth) != dpths.end())
1238                                        << " Phi: "
1239                                        << (std::find(actHB[i].iphis.begin(), actHB[i].iphis.end(), iphi) !=
1240                                            actHB[i].iphis.end())
1241                                        << " Length: " << actHB[i].thick << " : " << lx;
1242         }
1243       }
1244     }
1245   } else {
1246     for (unsigned int i = 0; i < actHE.size(); ++i) {
1247       if ((ieta == actHE[i].ieta) && (zside == actHE[i].zside) &&
1248           (std::find(dpths.begin(), dpths.end(), actHE[i].depth) != dpths.end()) &&
1249           (std::find(actHE[i].iphis.begin(), actHE[i].iphis.end(), iphi) != actHE[i].iphis.end())) {
1250         lx += actHE[i].thick;
1251         if (debug) {
1252           edm::LogVerbatim("HBHEMuon") << "Eta: " << (ieta == actHE[i].ieta) << " zside " << (zside == actHE[i].zside)
1253                                        << " Depth " << actHE[i].depth << ":"
1254                                        << (std::find(dpths.begin(), dpths.end(), actHE[i].depth) != dpths.end())
1255                                        << " Phi: "
1256                                        << (std::find(actHE[i].iphis.begin(), actHE[i].iphis.end(), iphi) !=
1257                                            actHE[i].iphis.end())
1258                                        << " Length: " << actHE[i].thick << " : " << lx;
1259         }
1260       }
1261     }
1262   }
1263   return lx;
1264 }
1265 
1266 bool HcalHBHEMuonAnalyzer::isGoodVertex(const reco::Vertex& vtx) {
1267   if (vtx.isFake())
1268     return false;
1269   if (vtx.ndof() < 4)
1270     return false;
1271   if (vtx.position().Rho() > 2.)
1272     return false;
1273   if (fabs(vtx.position().Z()) > 24.)
1274     return false;
1275   return true;
1276 }
1277 
1278 double HcalHBHEMuonAnalyzer::respCorr(const DetId& id) {
1279   double cfac(1.0);
1280   if (useMyCorr_) {
1281     auto itr = corrValue_.find(id);
1282     if (itr != corrValue_.end())
1283       cfac = itr->second;
1284   } else if (respCorrs_ != nullptr) {
1285     cfac = (respCorrs_->getValues(id))->getValue();
1286   }
1287   return cfac;
1288 }
1289 
1290 double HcalHBHEMuonAnalyzer::gainFactor(const HcalDbService* conditions, const HcalDetId& id) {
1291   double gain(0.0);
1292   const HcalCalibrations& calibs = conditions->getHcalCalibrations(id);
1293   for (int capid = 0; capid < 4; ++capid)
1294     gain += (0.25 * calibs.respcorrgain(capid));
1295   return gain;
1296 }
1297 
1298 int HcalHBHEMuonAnalyzer::depth16HE(int ieta, int iphi) {
1299   // Transition between HB/HE is special
1300   // For Run 1 or for Plan1 standard reconstruction it is 3
1301   // For runs beyond 2018 or in Plan1 for HEP17 it is 4
1302   int zside = (ieta > 0) ? 1 : -1;
1303   int depth = theHBHETopology_->dddConstants()->getMinDepth(1, 16, iphi, zside);
1304   if (isItPlan1_ && (!isItPreRecHit_))
1305     depth = 3;
1306 #ifdef EDM_ML_DEBUG
1307   edm::LogVerbatim("HBHEMuon") << "Plan1 " << isItPlan1_ << " PreRecHit " << isItPreRecHit_ << " phi " << iphi
1308                                << " depth " << depth;
1309 #endif
1310   return depth;
1311 }
1312 
1313 bool HcalHBHEMuonAnalyzer::goodCell(const HcalDetId& hcid,
1314                                     const reco::Track* pTrack,
1315                                     const CaloGeometry* geo,
1316                                     const MagneticField* bField) {
1317   std::pair<double, double> rz = hdc_->getRZ(hcid);
1318   bool typeRZ = (hcid.subdet() == HcalEndcap) ? false : true;
1319   bool match = spr::propagateHCAL(pTrack, geo, bField, typeRZ, rz, (((verbosity_ / 10000) % 10) > 0));
1320   return match;
1321 }
1322 
1323 //define this as a plug-in
1324 #include "FWCore/Framework/interface/MakerMacros.h"
1325 
1326 DEFINE_FWK_MODULE(HcalHBHEMuonAnalyzer);