Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-10-01 22:40:16

0001 // -*- C++ -*-
0002 //#define EDM_ML_DEBUG
0003 
0004 // system include files
0005 #include <atomic>
0006 #include <memory>
0007 #include <string>
0008 #include <cmath>
0009 #include <iostream>
0010 #include <sstream>
0011 #include <fstream>
0012 #include <vector>
0013 #include <boost/regex.hpp>
0014 
0015 // user include files
0016 #include "FWCore/Framework/interface/Frameworkfwd.h"
0017 #include "FWCore/Framework/interface/stream/EDProducer.h"
0018 #include "FWCore/Framework/interface/Event.h"
0019 #include "FWCore/Framework/interface/EventSetup.h"
0020 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0021 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0022 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0023 #include "FWCore/Common/interface/TriggerNames.h"
0024 #include "FWCore/Utilities/interface/Exception.h"
0025 
0026 #include "DataFormats/Common/interface/Handle.h"
0027 #include "DataFormats/Common/interface/Ref.h"
0028 #include "DataFormats/DetId/interface/DetId.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/EcalRecHit/interface/EcalRecHitCollections.h"
0033 #include "DataFormats/HcalDetId/interface/HcalSubdetector.h"
0034 #include "DataFormats/HcalDetId/interface/HcalDetId.h"
0035 #include "DataFormats/HcalRecHit/interface/HcalRecHitCollections.h"
0036 #include "DataFormats/HcalIsolatedTrack/interface/HcalIsolatedTrackCandidate.h"
0037 #include "DataFormats/HLTReco/interface/TriggerEvent.h"
0038 #include "DataFormats/HLTReco/interface/TriggerFilterObjectWithRefs.h"
0039 #include "DataFormats/L1Trigger/interface/L1JetParticleFwd.h"
0040 #include "DataFormats/L1Trigger/interface/L1JetParticle.h"
0041 #include "DataFormats/Provenance/interface/ProductID.h"
0042 #include "DataFormats/TrackReco/interface/Track.h"
0043 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0044 #include "DataFormats/TrackReco/interface/TrackExtra.h"
0045 #include "DataFormats/TrackReco/interface/TrackExtraFwd.h"
0046 #include "DataFormats/TrackReco/interface/HitPattern.h"
0047 #include "DataFormats/TrackReco/interface/TrackBase.h"
0048 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0049 #include "DataFormats/VertexReco/interface/Vertex.h"
0050 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0051 #include "DataFormats/Common/interface/TriggerResults.h"
0052 
0053 #include "HLTrigger/HLTcore/interface/HLTConfigProvider.h"
0054 
0055 #include "Calibration/IsolatedParticles/interface/CaloPropagateTrack.h"
0056 #include "Calibration/IsolatedParticles/interface/ChargeIsolation.h"
0057 #include "Calibration/IsolatedParticles/interface/eCone.h"
0058 #include "Calibration/IsolatedParticles/interface/TrackSelection.h"
0059 
0060 #include "MagneticField/Engine/interface/MagneticField.h"
0061 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0062 #include "Geometry/CaloGeometry/interface/CaloGeometry.h"
0063 #include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h"
0064 #include "Geometry/Records/interface/CaloTopologyRecord.h"
0065 #include "Geometry/CaloTopology/interface/CaloSubdetectorTopology.h"
0066 #include "Geometry/CaloTopology/interface/HcalTopology.h"
0067 #include "Geometry/CaloTopology/interface/CaloTopology.h"
0068 #include "Geometry/CaloTopology/interface/EcalTrigTowerConstituentsMap.h"
0069 #include "Geometry/Records/interface/CaloGeometryRecord.h"
0070 #include "CondFormats/DataRecord/interface/EcalChannelStatusRcd.h"
0071 #include "RecoLocalCalo/EcalRecAlgos/interface/EcalSeverityLevelAlgo.h"
0072 #include "RecoLocalCalo/EcalRecAlgos/interface/EcalSeverityLevelAlgoRcd.h"
0073 
0074 //#define EDM_ML_DEBUG
0075 //
0076 // class declaration
0077 //
0078 
0079 namespace AlCaIsoTracks {
0080   struct Counters {
0081     Counters() : nAll_(0), nGood_(0), nRange_(0) {}
0082     mutable std::atomic<unsigned int> nAll_, nGood_, nRange_;
0083   };
0084 }  // namespace AlCaIsoTracks
0085 
0086 class AlCaIsoTracksProducer : public edm::stream::EDProducer<edm::GlobalCache<AlCaIsoTracks::Counters> > {
0087 public:
0088   explicit AlCaIsoTracksProducer(edm::ParameterSet const&, const AlCaIsoTracks::Counters* count);
0089   ~AlCaIsoTracksProducer() override;
0090 
0091   static std::unique_ptr<AlCaIsoTracks::Counters> initializeGlobalCache(edm::ParameterSet const&) {
0092     return std::make_unique<AlCaIsoTracks::Counters>();
0093   }
0094 
0095   void produce(edm::Event&, edm::EventSetup const&) override;
0096   void endStream() override;
0097   static void globalEndJob(const AlCaIsoTracks::Counters* counters);
0098   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0099 
0100 private:
0101   void beginRun(edm::Run const&, edm::EventSetup const&) override;
0102   void endRun(edm::Run const&, edm::EventSetup const&) override;
0103   reco::HcalIsolatedTrackCandidateCollection* select(edm::Handle<edm::TriggerResults> const& triggerResults,
0104                                                      const std::vector<std::string>& triggerNames_,
0105                                                      edm::Handle<reco::TrackCollection>& trkCollection,
0106                                                      math::XYZPoint& leadPV,
0107                                                      edm::Handle<EcalRecHitCollection>& barrelRecHitsHandle,
0108                                                      edm::Handle<EcalRecHitCollection>& endcapRecHitsHandle,
0109                                                      edm::Handle<HBHERecHitCollection>& hbhe,
0110                                                      double ptL1,
0111                                                      double etaL1,
0112                                                      double phiL1);
0113   void setPtEtaPhi(std::vector<edm::Ref<l1extra::L1JetParticleCollection> >& objref,
0114                    double& ptL1,
0115                    double& etaL1,
0116                    double& phiL1);
0117 
0118   // ----------member data ---------------------------
0119   HLTConfigProvider hltConfig_;
0120   unsigned int nRun_, nAll_, nGood_, nRange_;
0121   spr::trackSelectionParameters selectionParameter_;
0122   const std::vector<std::string> trigNames_;
0123   const std::string theTrackQuality_, processName_;
0124   const double a_coneR_, a_mipR_;
0125   const double maxRestrictionP_, slopeRestrictionP_;
0126   const double pTrackMin_, eEcalMax_, eIsolate_;
0127   const double pTrackLow_, pTrackHigh_;
0128   const int preScale_;
0129   const edm::InputTag labelGenTrack_, labelRecVtx_, labelBS_;
0130   const edm::InputTag labelEB_, labelEE_, labelHBHE_, labelHltGT_;
0131   const edm::InputTag labelTriggerEvent_, labelTriggerResults_;
0132   const std::string labelIsoTk_;
0133   double a_charIsoR_;
0134   const MagneticField* bField;
0135   const CaloGeometry* geo;
0136 
0137   edm::EDGetTokenT<trigger::TriggerFilterObjectWithRefs> tok_hltGT_;
0138   edm::EDGetTokenT<trigger::TriggerEvent> tok_trigEvt_;
0139   edm::EDGetTokenT<edm::TriggerResults> tok_trigRes_;
0140   edm::EDGetTokenT<reco::TrackCollection> tok_genTrack_;
0141   edm::EDGetTokenT<reco::VertexCollection> tok_recVtx_;
0142   edm::EDGetTokenT<reco::BeamSpot> tok_bs_;
0143   edm::EDGetTokenT<EcalRecHitCollection> tok_EB_;
0144   edm::EDGetTokenT<EcalRecHitCollection> tok_EE_;
0145   edm::EDGetTokenT<HBHERecHitCollection> tok_hbhe_;
0146   edm::ESGetToken<CaloGeometry, CaloGeometryRecord> tok_geom_;
0147   edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> tok_magField_;
0148 };
0149 
0150 AlCaIsoTracksProducer::AlCaIsoTracksProducer(edm::ParameterSet const& iConfig, const AlCaIsoTracks::Counters* counters)
0151     : nRun_(0),
0152       nAll_(0),
0153       nGood_(0),
0154       nRange_(0),
0155       trigNames_(iConfig.getParameter<std::vector<std::string> >("triggers")),
0156       theTrackQuality_(iConfig.getParameter<std::string>("trackQuality")),
0157       processName_(iConfig.getParameter<std::string>("processName")),
0158       a_coneR_(iConfig.getParameter<double>("coneRadius")),
0159       a_mipR_(iConfig.getParameter<double>("coneRadiusMIP")),
0160       maxRestrictionP_(iConfig.getParameter<double>("maxTrackP")),
0161       slopeRestrictionP_(iConfig.getParameter<double>("slopeTrackP")),
0162       pTrackMin_(iConfig.getParameter<double>("minimumTrackP")),
0163       eEcalMax_(iConfig.getParameter<double>("maximumEcalEnergy")),
0164       eIsolate_(iConfig.getParameter<double>("isolationEnergy")),
0165       pTrackLow_(iConfig.getParameter<double>("momentumRangeLow")),
0166       pTrackHigh_(iConfig.getParameter<double>("momentumRangeHigh")),
0167       preScale_(iConfig.getParameter<int>("preScaleFactor")),
0168       labelGenTrack_(iConfig.getParameter<edm::InputTag>("TrackLabel")),
0169       labelRecVtx_(iConfig.getParameter<edm::InputTag>("VertexLabel")),
0170       labelBS_(iConfig.getParameter<edm::InputTag>("BeamSpotLabel")),
0171       labelEB_(iConfig.getParameter<edm::InputTag>("EBRecHitLabel")),
0172       labelEE_(iConfig.getParameter<edm::InputTag>("EERecHitLabel")),
0173       labelHBHE_(iConfig.getParameter<edm::InputTag>("HBHERecHitLabel")),
0174       labelHltGT_(iConfig.getParameter<edm::InputTag>("L1GTSeedLabel")),
0175       labelTriggerEvent_(iConfig.getParameter<edm::InputTag>("TriggerEventLabel")),
0176       labelTriggerResults_(iConfig.getParameter<edm::InputTag>("TriggerResultLabel")),
0177       labelIsoTk_(iConfig.getParameter<std::string>("IsoTrackLabel")) {
0178   // Get the run parameters
0179   // Different isolation cuts are described in DN-2016/029
0180   // Tight cut uses 2 GeV; Loose cut uses 10 GeV
0181   // Eta dependent cut uses (maxRestrictionP_ * exp(|ieta|*log(2.5)/18))
0182   // with the factor for exponential slopeRestrictionP_ = log(2.5)/18
0183   // maxRestrictionP_ = 8 GeV as came from a study
0184   const double isolationRadius(28.9);
0185   selectionParameter_.minPt = iConfig.getParameter<double>("minTrackPt");
0186   ;
0187   selectionParameter_.minQuality = reco::TrackBase::qualityByName(theTrackQuality_);
0188   selectionParameter_.maxDxyPV = iConfig.getParameter<double>("maxDxyPV");
0189   selectionParameter_.maxDzPV = iConfig.getParameter<double>("maxDzPV");
0190   selectionParameter_.maxChi2 = iConfig.getParameter<double>("maxChi2");
0191   selectionParameter_.maxDpOverP = iConfig.getParameter<double>("maxDpOverP");
0192   selectionParameter_.minOuterHit = iConfig.getParameter<int>("minOuterHit");
0193   selectionParameter_.minLayerCrossed = iConfig.getParameter<int>("minLayerCrossed");
0194   selectionParameter_.maxInMiss = iConfig.getParameter<int>("maxInMiss");
0195   selectionParameter_.maxOutMiss = iConfig.getParameter<int>("maxOutMiss");
0196   a_charIsoR_ = a_coneR_ + isolationRadius;
0197 
0198   // define tokens for access
0199   tok_hltGT_ = consumes<trigger::TriggerFilterObjectWithRefs>(labelHltGT_);
0200   tok_trigEvt_ = consumes<trigger::TriggerEvent>(labelTriggerEvent_);
0201   tok_trigRes_ = consumes<edm::TriggerResults>(labelTriggerResults_);
0202   tok_genTrack_ = consumes<reco::TrackCollection>(labelGenTrack_);
0203   tok_recVtx_ = consumes<reco::VertexCollection>(labelRecVtx_);
0204   tok_bs_ = consumes<reco::BeamSpot>(labelBS_);
0205   tok_EB_ = consumes<EcalRecHitCollection>(labelEB_);
0206   tok_EE_ = consumes<EcalRecHitCollection>(labelEE_);
0207   tok_hbhe_ = consumes<HBHERecHitCollection>(labelHBHE_);
0208   // for event setup
0209   tok_geom_ = esConsumes<CaloGeometry, CaloGeometryRecord, edm::Transition::BeginRun>();
0210   tok_magField_ = esConsumes<MagneticField, IdealMagneticFieldRecord, edm::Transition::BeginRun>();
0211 
0212   edm::LogVerbatim("HcalIsoTrack") << "Parameters read from config file \n"
0213                                    << "\t minPt " << selectionParameter_.minPt << "\t theTrackQuality "
0214                                    << theTrackQuality_ << "\t minQuality " << selectionParameter_.minQuality
0215                                    << "\t maxDxyPV " << selectionParameter_.maxDxyPV << "\t maxDzPV "
0216                                    << selectionParameter_.maxDzPV << "\t maxChi2 " << selectionParameter_.maxChi2
0217                                    << "\t maxDpOverP " << selectionParameter_.maxDpOverP << "\t minOuterHit "
0218                                    << selectionParameter_.minOuterHit << "\t minLayerCrossed "
0219                                    << selectionParameter_.minLayerCrossed << "\t maxInMiss "
0220                                    << selectionParameter_.maxInMiss << "\t maxOutMiss "
0221                                    << selectionParameter_.maxOutMiss << "\n"
0222                                    << "\t a_coneR " << a_coneR_ << "\t a_charIsoR " << a_charIsoR_ << "\t a_mipR "
0223                                    << a_mipR_ << "\t pTrackMin " << pTrackMin_ << "\t eEcalMax " << eEcalMax_
0224                                    << "\t maxRestrictionP_ " << maxRestrictionP_ << "\t slopeRestrictionP_ "
0225                                    << slopeRestrictionP_ << "\t eIsolate_ " << eIsolate_ << "\t Process "
0226                                    << processName_ << "\n"
0227                                    << "\t Precale factor " << preScale_ << "\t in momentum range " << pTrackLow_ << ":"
0228                                    << pTrackHigh_;
0229   for (unsigned int k = 0; k < trigNames_.size(); ++k)
0230     edm::LogVerbatim("HcalIsoTrack") << "Trigger[" << k << "] " << trigNames_[k];
0231 
0232   //create also IsolatedPixelTrackCandidateCollection which contains isolation info and reference to primary track
0233   produces<reco::HcalIsolatedTrackCandidateCollection>(labelIsoTk_);
0234   produces<reco::VertexCollection>(labelRecVtx_.label());
0235   produces<EcalRecHitCollection>(labelEB_.instance());
0236   produces<EcalRecHitCollection>(labelEE_.instance());
0237   produces<HBHERecHitCollection>(labelHBHE_.label());
0238 
0239   edm::LogVerbatim("HcalIsoTrack") << " Expected to produce the collections:\n"
0240                                    << "reco::HcalIsolatedTrackCandidateCollection "
0241                                    << " with label HcalIsolatedTrackCollection\n"
0242                                    << "reco::VertexCollection with label " << labelRecVtx_.label() << "\n"
0243                                    << "EcalRecHitCollection with label EcalRecHitsEB\n"
0244                                    << "EcalRecHitCollection with label EcalRecHitsEE\n"
0245                                    << "HBHERecHitCollection with label " << labelHBHE_.label();
0246 }
0247 
0248 AlCaIsoTracksProducer::~AlCaIsoTracksProducer() {}
0249 
0250 void AlCaIsoTracksProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0251   edm::ParameterSetDescription desc;
0252   // producer for  (HCAL isolated tracks)
0253   desc.add<edm::InputTag>("TrackLabel", edm::InputTag("generalTracks"));
0254   desc.add<edm::InputTag>("VertexLabel", edm::InputTag("offlinePrimaryVertices"));
0255   desc.add<edm::InputTag>("BeamSpotLabel", edm::InputTag("offlineBeamSpot"));
0256   desc.add<edm::InputTag>("EBRecHitLabel", edm::InputTag("ecalRecHit", "EcalRecHitsEB"));
0257   desc.add<edm::InputTag>("EERecHitLabel", edm::InputTag("ecalRecHit", "EcalRecHitsEE"));
0258   desc.add<edm::InputTag>("HBHERecHitLabel", edm::InputTag("hbhereco"));
0259   desc.add<edm::InputTag>("L1GTSeedLabel", edm::InputTag("hltL1sV0SingleJet60"));
0260   desc.add<edm::InputTag>("TriggerEventLabel", edm::InputTag("hltTriggerSummaryAOD", "", "HLT"));
0261   desc.add<edm::InputTag>("TriggerResultLabel", edm::InputTag("TriggerResults", "", "HLT"));
0262   desc.add<std::string>("IsoTrackLabel", "HcalIsolatedTrackCollection");
0263   std::vector<std::string> triggers = {"HLT_IsoTrackHB", "HLT_IsoTrackHE"};
0264   desc.add<std::vector<std::string> >("triggers", triggers);
0265   desc.add<std::string>("processName", "HLT");
0266   // following 10 parameters are parameters to select good tracks
0267   desc.add<std::string>("trackQuality", "highPurity");
0268   desc.add<double>("minTrackPt", 1.0);
0269   desc.add<double>("maxDxyPV", 10.0);
0270   desc.add<double>("maxDzPV", 100.0);
0271   desc.add<double>("maxChi2", 5.0);
0272   desc.add<double>("maxDpOverP", 0.1);
0273   desc.add<int>("minOuterHit", 4);
0274   desc.add<int>("minLayerCrossed", 8);
0275   desc.add<int>("maxInMiss", 2);
0276   desc.add<int>("maxOutMiss", 2);
0277   // Minimum momentum of selected isolated track and signal zone
0278   desc.add<double>("coneRadius", 34.98);
0279   desc.add<double>("minimumTrackP", 20.0);
0280   // signal zone in ECAL and MIP energy cutoff
0281   desc.add<double>("coneRadiusMIP", 14.0);
0282   desc.add<double>("maximumEcalEnergy", 2.0);
0283   // following 3 parameters are for isolation cuts and described in the code
0284   desc.add<double>("maxTrackP", 8.0);
0285   desc.add<double>("slopeTrackP", 0.05090504066);
0286   desc.add<double>("isolationEnergy", 10.0);
0287   // Prescale events only containing isolated tracks in the range
0288   desc.add<double>("momentumRangeLow", 20.0);
0289   desc.add<double>("momentumRangeHigh", 40.0);
0290   desc.add<int>("preScaleFactor", 10);
0291   descriptions.add("alcaisotrk", desc);
0292 }
0293 
0294 void AlCaIsoTracksProducer::produce(edm::Event& iEvent, edm::EventSetup const& iSetup) {
0295   nAll_++;
0296 #ifdef EDM_ML_DEBUG
0297   edm::LogVerbatim("HcalIsoTrack") << "Run " << iEvent.id().run() << " Event " << iEvent.id().event() << " Luminosity "
0298                                    << iEvent.luminosityBlock() << " Bunch " << iEvent.bunchCrossing();
0299 #endif
0300   bool valid(true);
0301   //Step1: Get all the relevant containers
0302   trigger::TriggerEvent triggerEvent;
0303   auto const& triggerEventHandle = iEvent.getHandle(tok_trigEvt_);
0304   if (!triggerEventHandle.isValid()) {
0305     edm::LogWarning("HcalIsoTrack") << "Cannot access the collection " << labelTriggerEvent_;
0306     valid = false;
0307   }
0308   auto const& triggerResults = iEvent.getHandle(tok_trigRes_);
0309   if (!triggerResults.isValid()) {
0310     edm::LogWarning("HcalIsoTrack") << "Cannot access the collection " << labelTriggerResults_;
0311     valid = false;
0312   }
0313 
0314   auto trkCollection = iEvent.getHandle(tok_genTrack_);
0315   if (!trkCollection.isValid()) {
0316     edm::LogWarning("HcalIsoTrack") << "Cannot access the collection " << labelGenTrack_;
0317     valid = false;
0318   }
0319 
0320   auto const& recVtxs = iEvent.getHandle(tok_recVtx_);
0321   if (!recVtxs.isValid()) {
0322     edm::LogWarning("HcalIsoTrack") << "Cannot access the collection " << labelGenTrack_;
0323     valid = false;
0324   }
0325 
0326   auto const& beamSpotH = iEvent.getHandle(tok_bs_);
0327   math::XYZPoint leadPV(0, 0, 0);
0328   if (valid) {
0329     if (!recVtxs->empty() && !((*recVtxs)[0].isFake())) {
0330       leadPV = math::XYZPoint((*recVtxs)[0].x(), (*recVtxs)[0].y(), (*recVtxs)[0].z());
0331     } else if (beamSpotH.isValid()) {
0332       leadPV = beamSpotH->position();
0333     }
0334   }
0335 #ifdef EDM_ML_DEBUG
0336   edm::LogVerbatim("HcalIsoTrack") << "Primary Vertex " << leadPV;
0337 #endif
0338 
0339   auto barrelRecHitsHandle = iEvent.getHandle(tok_EB_);
0340   if (!barrelRecHitsHandle.isValid()) {
0341     edm::LogWarning("HcalIsoTrack") << "Cannot access the collection " << labelEB_;
0342     valid = false;
0343   }
0344   auto endcapRecHitsHandle = iEvent.getHandle(tok_EE_);
0345   if (!endcapRecHitsHandle.isValid()) {
0346     edm::LogWarning("HcalIsoTrack") << "Cannot access the collection " << labelEE_;
0347     valid = false;
0348   }
0349   auto hbhe = iEvent.getHandle(tok_hbhe_);
0350   if (!hbhe.isValid()) {
0351     edm::LogWarning("HcalIsoTrack") << "Cannot access the collection " << labelHBHE_;
0352     valid = false;
0353   }
0354 
0355   //Get L1 trigger object
0356   double ptL1(0), etaL1(0), phiL1(0);
0357   auto const& l1trigobj = iEvent.getHandle(tok_hltGT_);
0358 
0359   if (l1trigobj.isValid()) {
0360     std::vector<edm::Ref<l1extra::L1JetParticleCollection> > l1tauobjref;
0361     l1trigobj->getObjects(trigger::TriggerL1TauJet, l1tauobjref);
0362     setPtEtaPhi(l1tauobjref, ptL1, etaL1, phiL1);
0363 
0364     std::vector<edm::Ref<l1extra::L1JetParticleCollection> > l1jetobjref;
0365     l1trigobj->getObjects(trigger::TriggerL1CenJet, l1jetobjref);
0366     setPtEtaPhi(l1jetobjref, ptL1, etaL1, phiL1);
0367 
0368     std::vector<edm::Ref<l1extra::L1JetParticleCollection> > l1forjetobjref;
0369     l1trigobj->getObjects(trigger::TriggerL1ForJet, l1forjetobjref);
0370     setPtEtaPhi(l1forjetobjref, ptL1, etaL1, phiL1);
0371   } else {
0372     valid = false;
0373   }
0374 
0375   auto outputHcalIsoTrackColl = std::make_unique<reco::HcalIsolatedTrackCandidateCollection>();
0376   auto outputVColl = std::make_unique<reco::VertexCollection>();
0377   auto outputEBColl = std::make_unique<EBRecHitCollection>();
0378   auto outputEEColl = std::make_unique<EERecHitCollection>();
0379   auto outputHBHEColl = std::make_unique<HBHERecHitCollection>();
0380 
0381   //For valid HLT record
0382   if (!valid) {
0383     edm::LogWarning("HcalIsoTrack") << "Error! Can't get some of the products";
0384   } else {
0385     trigger::TriggerEvent triggerEvent = *(triggerEventHandle.product());
0386     if (triggerResults.isValid()) {
0387       const edm::TriggerNames& triggerNames = iEvent.triggerNames(*triggerResults);
0388       const std::vector<std::string>& triggerNames_ = triggerNames.triggerNames();
0389       reco::HcalIsolatedTrackCandidateCollection* isotk = select(triggerResults,
0390                                                                  triggerNames_,
0391                                                                  trkCollection,
0392                                                                  leadPV,
0393                                                                  barrelRecHitsHandle,
0394                                                                  endcapRecHitsHandle,
0395                                                                  hbhe,
0396                                                                  ptL1,
0397                                                                  etaL1,
0398                                                                  phiL1);
0399 #ifdef EDM_ML_DEBUG
0400       edm::LogVerbatim("HcalIsoTrack") << "AlCaIsoTracksProducer::select returns " << isotk->size()
0401                                        << " isolated tracks";
0402 #endif
0403 
0404       if (!isotk->empty()) {
0405         int ntrin(0), ntrout(0);
0406         for (reco::HcalIsolatedTrackCandidateCollection::const_iterator itr = isotk->begin(); itr != isotk->end();
0407              ++itr) {
0408           if (itr->p() > pTrackLow_ && itr->p() < pTrackHigh_)
0409             ntrin++;
0410           else
0411             ntrout++;
0412         }
0413         bool selectEvent = ntrout > 0;
0414         if (!selectEvent && ntrin > 0) {
0415           ++nRange_;
0416           if (preScale_ <= 1)
0417             selectEvent = true;
0418           else if (nRange_ % preScale_ == 1)
0419             selectEvent = true;
0420         }
0421         if (selectEvent) {
0422           for (reco::HcalIsolatedTrackCandidateCollection::const_iterator itr = isotk->begin(); itr != isotk->end();
0423                ++itr)
0424             outputHcalIsoTrackColl->push_back(*itr);
0425 
0426           for (reco::VertexCollection::const_iterator vtx = recVtxs->begin(); vtx != recVtxs->end(); ++vtx)
0427             outputVColl->push_back(*vtx);
0428 
0429           for (edm::SortedCollection<EcalRecHit>::const_iterator ehit = barrelRecHitsHandle->begin();
0430                ehit != barrelRecHitsHandle->end();
0431                ++ehit)
0432             outputEBColl->push_back(*ehit);
0433 
0434           for (edm::SortedCollection<EcalRecHit>::const_iterator ehit = endcapRecHitsHandle->begin();
0435                ehit != endcapRecHitsHandle->end();
0436                ++ehit)
0437             outputEEColl->push_back(*ehit);
0438 
0439           for (std::vector<HBHERecHit>::const_iterator hhit = hbhe->begin(); hhit != hbhe->end(); ++hhit)
0440             outputHBHEColl->push_back(*hhit);
0441           ++nGood_;
0442         }
0443       }
0444     }
0445   }
0446   iEvent.put(std::move(outputHcalIsoTrackColl), labelIsoTk_);
0447   iEvent.put(std::move(outputVColl), labelRecVtx_.label());
0448   iEvent.put(std::move(outputEBColl), labelEB_.instance());
0449   iEvent.put(std::move(outputEEColl), labelEE_.instance());
0450   iEvent.put(std::move(outputHBHEColl), labelHBHE_.label());
0451 }
0452 
0453 void AlCaIsoTracksProducer::endStream() {
0454   globalCache()->nAll_ += nAll_;
0455   globalCache()->nGood_ += nGood_;
0456   globalCache()->nRange_ += nRange_;
0457 }
0458 
0459 void AlCaIsoTracksProducer::globalEndJob(const AlCaIsoTracks::Counters* count) {
0460   edm::LogVerbatim("HcalIsoTrack") << "Finds " << count->nGood_ << " good tracks in " << count->nAll_ << " events and "
0461                                    << count->nRange_ << " events in the momentum raange";
0462 }
0463 
0464 void AlCaIsoTracksProducer::beginRun(edm::Run const& iRun, edm::EventSetup const& iSetup) {
0465   bool changed(false);
0466   edm::LogVerbatim("HcalIsoTrack") << "Run[" << nRun_ << "] " << iRun.run() << " hltconfig.init "
0467                                    << hltConfig_.init(iRun, iSetup, processName_, changed);
0468 
0469   bField = &(iSetup.getData(tok_magField_));
0470   geo = &(iSetup.getData(tok_geom_));
0471 }
0472 
0473 void AlCaIsoTracksProducer::endRun(edm::Run const& iRun, edm::EventSetup const&) {
0474   edm::LogVerbatim("HcalIsoTrack") << "endRun [" << nRun_ << "] " << iRun.run();
0475   ++nRun_;
0476 }
0477 
0478 reco::HcalIsolatedTrackCandidateCollection* AlCaIsoTracksProducer::select(
0479     edm::Handle<edm::TriggerResults> const& triggerResults,
0480     const std::vector<std::string>& triggerNames_,
0481     edm::Handle<reco::TrackCollection>& trkCollection,
0482     math::XYZPoint& leadPV,
0483     edm::Handle<EcalRecHitCollection>& barrelRecHitsHandle,
0484     edm::Handle<EcalRecHitCollection>& endcapRecHitsHandle,
0485     edm::Handle<HBHERecHitCollection>& hbhe,
0486     double ptL1,
0487     double etaL1,
0488     double phiL1) {
0489   reco::HcalIsolatedTrackCandidateCollection* trackCollection = new reco::HcalIsolatedTrackCandidateCollection;
0490   bool ok(false);
0491 
0492   // Find a good HLT trigger
0493   for (unsigned int iHLT = 0; iHLT < triggerResults->size(); iHLT++) {
0494     int hlt = triggerResults->accept(iHLT);
0495     for (unsigned int i = 0; i < trigNames_.size(); ++i) {
0496       if (triggerNames_[iHLT].find(trigNames_[i]) != std::string::npos) {
0497         if (hlt > 0) {
0498           ok = true;
0499         }
0500         edm::LogVerbatim("HcalIsoTrack") << "The trigger we are looking for " << triggerNames_[iHLT] << " Flag " << hlt
0501                                          << ":" << ok;
0502       }
0503     }
0504   }
0505 
0506   //Propagate tracks to calorimeter surface)
0507   std::vector<spr::propagatedTrackDirection> trkCaloDirections;
0508   spr::propagateCALO(trkCollection, geo, bField, theTrackQuality_, trkCaloDirections, false);
0509 
0510   std::vector<spr::propagatedTrackDirection>::const_iterator trkDetItr;
0511   unsigned int nTracks(0), nselTracks(0);
0512   for (trkDetItr = trkCaloDirections.begin(), nTracks = 0; trkDetItr != trkCaloDirections.end();
0513        trkDetItr++, nTracks++) {
0514     const reco::Track* pTrack = &(*(trkDetItr->trkItr));
0515     math::XYZTLorentzVector v4(pTrack->px(), pTrack->py(), pTrack->pz(), pTrack->p());
0516 #ifdef EDM_ML_DEBUG
0517     edm::LogVerbatim("HcalIsoTrack") << "This track : " << nTracks << " (pt|eta|phi|p) :" << pTrack->pt() << "|"
0518                                      << pTrack->eta() << "|" << pTrack->phi() << "|" << pTrack->p();
0519 #endif
0520     //Selection of good track
0521     bool qltyFlag = spr::goodTrack(pTrack, leadPV, selectionParameter_, false);
0522 #ifdef EDM_ML_DEBUG
0523     edm::LogVerbatim("HcalIsoTrack") << "qltyFlag|okECAL|okHCAL : " << qltyFlag << "|" << trkDetItr->okECAL << "|"
0524                                      << trkDetItr->okHCAL;
0525 #endif
0526     if (qltyFlag && trkDetItr->okECAL && trkDetItr->okHCAL) {
0527       double t_p = pTrack->p();
0528       nselTracks++;
0529       int nRH_eMipDR(0), nNearTRKs(0);
0530       double eMipDR = spr::eCone_ecal(geo,
0531                                       barrelRecHitsHandle,
0532                                       endcapRecHitsHandle,
0533                                       trkDetItr->pointHCAL,
0534                                       trkDetItr->pointECAL,
0535                                       a_mipR_,
0536                                       trkDetItr->directionECAL,
0537                                       nRH_eMipDR);
0538       double hmaxNearP = spr::chargeIsolationCone(nTracks, trkCaloDirections, a_charIsoR_, nNearTRKs, false);
0539       HcalDetId detId = (HcalDetId)(trkDetItr->detIdHCAL);
0540       int ieta = detId.ietaAbs();
0541       double eIsolation = (maxRestrictionP_ * exp(slopeRestrictionP_ * ((double)(ieta))));
0542       if (eIsolation < eIsolate_)
0543         eIsolation = eIsolate_;
0544 #ifdef EDM_ML_DEBUG
0545       edm::LogVerbatim("HcalIsoTrack") << "This track : " << nTracks << " (pt|eta|phi|p) :" << pTrack->pt() << "|"
0546                                        << pTrack->eta() << "|" << pTrack->phi() << "|" << t_p << " e_MIP " << eMipDR
0547                                        << " Chg Isolation " << hmaxNearP << ":" << eIsolation;
0548 #endif
0549       if (t_p > pTrackMin_ && eMipDR < eEcalMax_ && hmaxNearP < eIsolation) {
0550         reco::HcalIsolatedTrackCandidate newCandidate(v4);
0551         newCandidate.SetMaxP(hmaxNearP);
0552         newCandidate.SetEnergyEcal(eMipDR);
0553         newCandidate.setL1(ptL1, etaL1, phiL1);
0554         newCandidate.SetEtaPhiEcal((trkDetItr->pointECAL).eta(), (trkDetItr->pointECAL).phi());
0555         HcalDetId detId = HcalDetId(trkDetItr->detIdHCAL);
0556         newCandidate.SetEtaPhiHcal(
0557             (trkDetItr->pointHCAL).eta(), (trkDetItr->pointHCAL).phi(), detId.ieta(), detId.iphi());
0558         int indx(0);
0559         for (reco::TrackCollection::const_iterator trkItr1 = trkCollection->begin(); trkItr1 != trkCollection->end();
0560              ++trkItr1, ++indx) {
0561           const reco::Track* pTrack1 = &(*trkItr1);
0562           if (pTrack1 == pTrack) {
0563             reco::TrackRef tRef = reco::TrackRef(trkCollection, indx);
0564             newCandidate.setTrack(tRef);
0565             break;
0566           }
0567         }
0568         trackCollection->push_back(newCandidate);
0569       }
0570     }
0571   }
0572   return trackCollection;
0573 }
0574 
0575 void AlCaIsoTracksProducer::setPtEtaPhi(std::vector<edm::Ref<l1extra::L1JetParticleCollection> >& objref,
0576                                         double& ptL1,
0577                                         double& etaL1,
0578                                         double& phiL1) {
0579   for (unsigned int p = 0; p < objref.size(); p++) {
0580     if (objref[p]->pt() > ptL1) {
0581       ptL1 = objref[p]->pt();
0582       phiL1 = objref[p]->phi();
0583       etaL1 = objref[p]->eta();
0584     }
0585   }
0586 }
0587 
0588 #include "FWCore/Framework/interface/MakerMacros.h"
0589 
0590 DEFINE_FWK_MODULE(AlCaIsoTracksProducer);