EgammaHLTProducers/plugins/HLTRecHitInAllL1RegionsProducer.cc

0001 #include <memory>
0002
0003 #include "FWCore/Framework/interface/Frameworkfwd.h"
0004 #include "FWCore/Framework/interface/stream/EDProducer.h"
0005 #include "FWCore/Framework/interface/Event.h"
0006 #include "FWCore/Framework/interface/EventSetup.h"
0007 #include "FWCore/Framework/interface/MakerMacros.h"
0008 #include "FWCore/Framework/interface/ConsumesCollector.h"
0009 #include "FWCore/Framework/interface/ESHandle.h"
0010
0011 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0012
0013 // Reco candidates (might not need)
0014 #include "DataFormats/RecoCandidate/interface/RecoCandidate.h"
0015 #include "DataFormats/RecoCandidate/interface/RecoChargedCandidate.h"
0016
0017 // Geometry and topology
0018 #include "Geometry/Records/interface/CaloGeometryRecord.h"
0019 #include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h"
0020 #include "Geometry/CaloGeometry/interface/CaloCellGeometry.h"
0021 #include "Geometry/CaloGeometry/interface/CaloGeometry.h"
0022 #include "DataFormats/Math/interface/RectangularEtaPhiRegion.h"
0023
0024 // Level 1 Trigger
0025 #include "DataFormats/L1Trigger/interface/L1EmParticleFwd.h"
0026 #include "DataFormats/L1Trigger/interface/L1JetParticleFwd.h"
0027 #include "DataFormats/L1Trigger/interface/L1MuonParticleFwd.h"
0028 #include "DataFormats/L1Trigger/interface/L1EmParticle.h"
0029 #include "DataFormats/L1Trigger/interface/L1JetParticle.h"
0030 #include "DataFormats/L1Trigger/interface/L1MuonParticle.h"
0031
0032 #include "DataFormats/L1Trigger/interface/EGamma.h"
0033 #include "DataFormats/L1Trigger/interface/Jet.h"
0034 #include "DataFormats/L1Trigger/interface/Muon.h"
0035
0036 #include "CondFormats/L1TObjects/interface/L1CaloGeometry.h"
0037 #include "CondFormats/DataRecord/interface/L1CaloGeometryRecord.h"
0038
0039 #include "DataFormats/EcalRecHit/interface/EcalRecHit.h"
0040 #include "DataFormats/EcalRecHit/interface/EcalUncalibratedRecHit.h"
0041
0042 #include "L1Trigger/L1TCalorimeter/interface/CaloTools.h"
0043
0044 #include "HLTrigger/HLTcore/interface/defaultModuleLabel.h"
0045
0046 //this is a re-write of HLTRechitInRegionsProducer to be able to handle arbitary L1 collections as inputs
0047 //in the process, some of the cruft was cleaned up but it mantains almost all the old behaviour
0048 //think the only difference now is that it wont throw if its not ECALBarrel, ECALEndcap or ECAL PS rec-hit type
0049 class L1RegionDataBase {
0050 public:
0051   virtual ~L1RegionDataBase() {}
0052   virtual void getEtaPhiRegions(const edm::Event&,
0053                                 const edm::EventSetup&,
0054                                 std::vector<RectangularEtaPhiRegion>&) const = 0;
0055 };
0056
0057 template <typename T1>
0058 class L1RegionData : public L1RegionDataBase {
0059 private:
0060   double const minEt_;
0061   double const maxEt_;
0062   double const regionEtaMargin_;
0063   double const regionPhiMargin_;
0064   edm::EDGetTokenT<T1> const token_;
0065   edm::ESGetToken<L1CaloGeometry, L1CaloGeometryRecord> l1CaloGeometryToken_;
0066
0067   void eventSetupConsumes(edm::ConsumesCollector& consumesColl);
0068
0069 public:
0070   L1RegionData(const edm::ParameterSet& para, edm::ConsumesCollector& consumesColl)
0071       : minEt_(para.getParameter<double>("minEt")),
0072         maxEt_(para.getParameter<double>("maxEt")),
0073         regionEtaMargin_(para.getParameter<double>("regionEtaMargin")),
0074         regionPhiMargin_(para.getParameter<double>("regionPhiMargin")),
0075         token_(consumesColl.consumes<T1>(para.getParameter<edm::InputTag>("inputColl"))) {
0076     eventSetupConsumes(consumesColl);
0077   }
0078
0079   void getEtaPhiRegions(const edm::Event&,
0080                         const edm::EventSetup&,
0081                         std::vector<RectangularEtaPhiRegion>&) const override;
0082   template <typename T2>
0083   static typename T2::const_iterator beginIt(const T2& coll) {
0084     return coll.begin();
0085   }
0086   template <typename T2>
0087   static typename T2::const_iterator endIt(const T2& coll) {
0088     return coll.end();
0089   }
0090   template <typename T2>
0091   static typename BXVector<T2>::const_iterator beginIt(const BXVector<T2>& coll) {
0092     return coll.begin(0);
0093   }
0094   template <typename T2>
0095   static typename BXVector<T2>::const_iterator endIt(const BXVector<T2>& coll) {
0096     return coll.end(0);
0097   }
0098 };
0099
0100 template <typename RecHitType>
0101 class HLTRecHitInAllL1RegionsProducer : public edm::stream::EDProducer<> {
0102   using RecHitCollectionType = edm::SortedCollection<RecHitType>;
0103
0104 public:
0105   HLTRecHitInAllL1RegionsProducer(const edm::ParameterSet& ps);
0106   ~HLTRecHitInAllL1RegionsProducer() override {}
0107
0108   void produce(edm::Event&, const edm::EventSetup&) override;
0109   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0110
0111 private:
0112   L1RegionDataBase* createL1RegionData(const std::string&,
0113                                        const edm::ParameterSet&,
0114                                        edm::ConsumesCollector&&);  //calling function owns this
0115
0116   std::vector<std::unique_ptr<L1RegionDataBase>> l1RegionData_;
0117
0118   std::vector<edm::InputTag> recHitLabels_;
0119   std::vector<std::string> productLabels_;
0120
0121   std::vector<edm::EDGetTokenT<RecHitCollectionType>> recHitTokens_;
0122
0123   const edm::ESGetToken<CaloGeometry, CaloGeometryRecord> caloGeometryToken_;
0124 };
0125
0126 template <typename RecHitType>
0127 HLTRecHitInAllL1RegionsProducer<RecHitType>::HLTRecHitInAllL1RegionsProducer(const edm::ParameterSet& para)
0128     : caloGeometryToken_{esConsumes()} {
0129   const std::vector<edm::ParameterSet> l1InputRegions =
0130       para.getParameter<std::vector<edm::ParameterSet>>("l1InputRegions");
0131   for (auto& pset : l1InputRegions) {
0132     const std::string type = pset.getParameter<std::string>("type");
0133     // meh I was going to use a factory but it was going to be overly complex for my needs
0134     l1RegionData_.emplace_back(createL1RegionData(type, pset, consumesCollector()));
0135   }
0136   recHitLabels_ = para.getParameter<std::vector<edm::InputTag>>("recHitLabels");
0137   productLabels_ = para.getParameter<std::vector<std::string>>("productLabels");
0138
0139   for (unsigned int collNr = 0; collNr < recHitLabels_.size(); collNr++) {
0140     recHitTokens_.push_back(consumes<RecHitCollectionType>(recHitLabels_[collNr]));
0141     produces<RecHitCollectionType>(productLabels_[collNr]);
0142   }
0143 }
0144
0145 template <typename RecHitType>
0146 void HLTRecHitInAllL1RegionsProducer<RecHitType>::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0147   edm::ParameterSetDescription desc;
0148   std::vector<std::string> productTags;
0149   productTags.push_back("EcalRegionalRecHitsEB");
0150   productTags.push_back("EcalRegionalRecHitsEE");
0151   desc.add<std::vector<std::string>>("productLabels", productTags);
0152   std::vector<edm::InputTag> recHitLabels;
0153   recHitLabels.push_back(edm::InputTag("hltEcalRegionalEgammaRecHit:EcalRecHitsEB"));
0154   recHitLabels.push_back(edm::InputTag("hltEcalRegionalEgammaRecHit:EcalRecHitsEE"));
0155   recHitLabels.push_back(edm::InputTag("hltESRegionalEgammaRecHit:EcalRecHitsES"));
0156   desc.add<std::vector<edm::InputTag>>("recHitLabels", recHitLabels);
0157   std::vector<edm::ParameterSet> l1InputRegions;
0158
0159   edm::ParameterSet emIsoPSet;
0160   emIsoPSet.addParameter<std::string>("type", "L1EmParticle");
0161   emIsoPSet.addParameter<double>("minEt", 5);
0162   emIsoPSet.addParameter<double>("maxEt", 999);
0163   emIsoPSet.addParameter<double>("regionEtaMargin", 0.14);
0164   emIsoPSet.addParameter<double>("regionPhiMargin", 0.4);
0165   emIsoPSet.addParameter<edm::InputTag>("inputColl", edm::InputTag("hltL1extraParticles:NonIsolated"));
0166   l1InputRegions.push_back(emIsoPSet);
0167   edm::ParameterSet emNonIsoPSet;
0168   emNonIsoPSet.addParameter<std::string>("type", "L1EmParticle");
0169   emNonIsoPSet.addParameter<double>("minEt", 5);
0170   emNonIsoPSet.addParameter<double>("maxEt", 999);
0171   emNonIsoPSet.addParameter<double>("regionEtaMargin", 0.14);
0172   emNonIsoPSet.addParameter<double>("regionPhiMargin", 0.4);
0173   emNonIsoPSet.addParameter<edm::InputTag>("inputColl", edm::InputTag("hltL1extraParticles:Isolated"));
0174   l1InputRegions.push_back(emNonIsoPSet);
0175
0176   // Why no Central Jets here? They are present in the python config, e.g. OnLine_HLT_GRun.py
0177   // SHarper: because these are the default parameters designed to reproduce the original (no jets) behaviour
0178   //
0179   edm::ParameterSet egPSet;
0180   egPSet.addParameter<std::string>("type", "EGamma");
0181   egPSet.addParameter<double>("minEt", 5);
0182   egPSet.addParameter<double>("maxEt", 999);
0183   egPSet.addParameter<double>("regionEtaMargin", 0.4);
0184   egPSet.addParameter<double>("regionPhiMargin", 0.5);
0185   egPSet.addParameter<edm::InputTag>("inputColl", edm::InputTag("hltCaloStage2Digis"));
0186   l1InputRegions.push_back(egPSet);
0187
0188   edm::ParameterSet jetPSet;
0189   jetPSet.addParameter<std::string>("type", "EGamma");
0190   jetPSet.addParameter<double>("minEt", 200);
0191   jetPSet.addParameter<double>("maxEt", 999);
0192   jetPSet.addParameter<double>("regionEtaMargin", 0.4);
0193   jetPSet.addParameter<double>("regionPhiMargin", 0.5);
0194   jetPSet.addParameter<edm::InputTag>("inputColl", edm::InputTag("hltCaloStage2Digis"));
0195   l1InputRegions.push_back(jetPSet);
0196
0197   edm::ParameterSetDescription l1InputRegionDesc;
0198   l1InputRegionDesc.add<std::string>("type");
0199   l1InputRegionDesc.add<double>("minEt");
0200   l1InputRegionDesc.add<double>("maxEt");
0201   l1InputRegionDesc.add<double>("regionEtaMargin");
0202   l1InputRegionDesc.add<double>("regionPhiMargin");
0203   l1InputRegionDesc.add<edm::InputTag>("inputColl");
0204   desc.addVPSet("l1InputRegions", l1InputRegionDesc, l1InputRegions);
0205
0206   descriptions.add(defaultModuleLabel<HLTRecHitInAllL1RegionsProducer<RecHitType>>(), desc);
0207 }
0208
0209 template <typename RecHitType>
0210 void HLTRecHitInAllL1RegionsProducer<RecHitType>::produce(edm::Event& event, const edm::EventSetup& setup) {
0211   // get the collection geometry:
0212   auto const& caloGeom = setup.getData(caloGeometryToken_);
0213
0214   std::vector<RectangularEtaPhiRegion> regions;
0215   std::for_each(l1RegionData_.begin(),
0216                 l1RegionData_.end(),
0217                 [&event, &setup, &regions](const std::unique_ptr<L1RegionDataBase>& input) {
0218                   input->getEtaPhiRegions(event, setup, regions);
0219                 });
0220
0221   for (size_t recHitCollNr = 0; recHitCollNr < recHitTokens_.size(); recHitCollNr++) {
0222     edm::Handle<RecHitCollectionType> recHits;
0223     event.getByToken(recHitTokens_[recHitCollNr], recHits);
0224
0225     if (!(recHits.isValid())) {
0226       edm::LogError("ProductNotFound") << "could not get a handle on the " << typeid(RecHitCollectionType).name()
0227                                        << " named " << recHitLabels_[recHitCollNr].encode() << std::endl;
0228       continue;
0229     }
0230
0231     auto filteredRecHits = std::make_unique<RecHitCollectionType>();
0232
0233     if (!recHits->empty()) {
0234       const CaloSubdetectorGeometry* subDetGeom = caloGeom.getSubdetectorGeometry(recHits->front().id());
0235       if (!regions.empty()) {
0236         for (const RecHitType& recHit : *recHits) {
0237           auto this_cell = subDetGeom->getGeometry(recHit.id());
0238           for (const auto& region : regions) {
0239             if (region.inRegion(this_cell->etaPos(), this_cell->phiPos())) {
0240               filteredRecHits->push_back(recHit);
0241               break;
0242             }
0243           }
0244         }
0245       }  //end check of empty regions
0246     }    //end check of empty rec-hits
0247     //   std::cout <<"putting fileter coll in "<<filteredRecHits->size()<<std::endl;
0248     event.put(std::move(filteredRecHits), productLabels_[recHitCollNr]);
0249   }  //end loop over all rec hit collections
0250 }
0251
0252 template <typename RecHitType>
0253 L1RegionDataBase* HLTRecHitInAllL1RegionsProducer<RecHitType>::createL1RegionData(
0254     const std::string& type, const edm::ParameterSet& para, edm::ConsumesCollector&& consumesColl) {
0255   if (type == "L1EmParticle") {
0256     return new L1RegionData<l1extra::L1EmParticleCollection>(para, consumesColl);
0257   } else if (type == "L1JetParticle") {
0258     return new L1RegionData<l1extra::L1JetParticleCollection>(para, consumesColl);
0259   } else if (type == "L1MuonParticle") {
0260     return new L1RegionData<l1extra::L1MuonParticleCollection>(para, consumesColl);
0261   } else if (type == "EGamma") {
0262     return new L1RegionData<l1t::EGammaBxCollection>(para, consumesColl);
0263   } else if (type == "Jet") {
0264     return new L1RegionData<l1t::JetBxCollection>(para, consumesColl);
0265   } else if (type == "Muon") {
0266     return new L1RegionData<l1t::MuonBxCollection>(para, consumesColl);
0267   } else if (type == "Tau") {
0268     return new L1RegionData<l1t::TauBxCollection>(para, consumesColl);
0269   } else {
0270     //this is a major issue and could lead to rather subtle efficiency losses, so if its incorrectly configured, we're aborting the job!
0271     throw cms::Exception("InvalidConfig")
0272         << " type " << type
0273         << " is not recognised, this means the rec-hit you think you are keeping may not be and you should fix this "
0274            "error as it can lead to hard to find efficiency loses"
0275         << std::endl;
0276   }
0277 }
0278
0279 template <typename L1CollType>
0280 void L1RegionData<L1CollType>::eventSetupConsumes(edm::ConsumesCollector&) {}
0281
0282 template <typename L1CollType>
0283 void L1RegionData<L1CollType>::getEtaPhiRegions(const edm::Event& event,
0284                                                 const edm::EventSetup&,
0285                                                 std::vector<RectangularEtaPhiRegion>& regions) const {
0286   edm::Handle<L1CollType> l1Cands;
0287   event.getByToken(token_, l1Cands);
0288
0289   for (auto l1CandIt = beginIt(*l1Cands); l1CandIt != endIt(*l1Cands); ++l1CandIt) {
0290     if (l1CandIt->et() >= minEt_ && l1CandIt->et() < maxEt_) {
0291       double etaLow = l1CandIt->eta() - regionEtaMargin_;
0292       double etaHigh = l1CandIt->eta() + regionEtaMargin_;
0293       double phiLow = l1CandIt->phi() - regionPhiMargin_;
0294       double phiHigh = l1CandIt->phi() + regionPhiMargin_;
0295
0296       regions.push_back(RectangularEtaPhiRegion(etaLow, etaHigh, phiLow, phiHigh));
0297     }
0298   }
0299 }
0300
0301 template <>
0302 void L1RegionData<l1extra::L1JetParticleCollection>::eventSetupConsumes(edm::ConsumesCollector& consumesColl) {
0303   l1CaloGeometryToken_ = consumesColl.esConsumes();
0304 }
0305
0306 template <>
0307 void L1RegionData<l1extra::L1JetParticleCollection>::getEtaPhiRegions(
0308     const edm::Event& event, const edm::EventSetup& setup, std::vector<RectangularEtaPhiRegion>& regions) const {
0309   edm::Handle<l1extra::L1JetParticleCollection> l1Cands;
0310   event.getByToken(token_, l1Cands);
0311
0312   auto const& l1CaloGeom = setup.getData(l1CaloGeometryToken_);
0313
0314   for (const auto& l1Cand : *l1Cands) {
0315     if (l1Cand.et() >= minEt_ && l1Cand.et() < maxEt_) {
0316       // Access the GCT hardware object corresponding to the L1Extra EM object.
0317       int etaIndex = l1Cand.gctJetCand()->etaIndex();
0318       int phiIndex = l1Cand.gctJetCand()->phiIndex();
0319
0320       // Use the L1CaloGeometry to find the eta, phi bin boundaries.
0321       double etaLow = l1CaloGeom.etaBinLowEdge(etaIndex);
0322       double etaHigh = l1CaloGeom.etaBinHighEdge(etaIndex);
0323       double phiLow = l1CaloGeom.emJetPhiBinLowEdge(phiIndex);
0324       double phiHigh = l1CaloGeom.emJetPhiBinHighEdge(phiIndex);
0325
0326       etaLow -= regionEtaMargin_;
0327       etaHigh += regionEtaMargin_;
0328       phiLow -= regionPhiMargin_;
0329       phiHigh += regionPhiMargin_;
0330
0331       regions.push_back(RectangularEtaPhiRegion(etaLow, etaHigh, phiLow, phiHigh));
0332     }
0333   }
0334 }
0335
0336 template <>
0337 void L1RegionData<l1extra::L1EmParticleCollection>::eventSetupConsumes(edm::ConsumesCollector& consumesColl) {
0338   l1CaloGeometryToken_ = consumesColl.esConsumes();
0339 }
0340
0341 template <>
0342 void L1RegionData<l1extra::L1EmParticleCollection>::getEtaPhiRegions(
0343     const edm::Event& event, const edm::EventSetup& setup, std::vector<RectangularEtaPhiRegion>& regions) const {
0344   edm::Handle<l1extra::L1EmParticleCollection> l1Cands;
0345   event.getByToken(token_, l1Cands);
0346
0347   auto const& l1CaloGeom = setup.getData(l1CaloGeometryToken_);
0348
0349   for (const auto& l1Cand : *l1Cands) {
0350     if (l1Cand.et() >= minEt_ && l1Cand.et() < maxEt_) {
0351       // Access the GCT hardware object corresponding to the L1Extra EM object.
0352       int etaIndex = l1Cand.gctEmCand()->etaIndex();
0353       int phiIndex = l1Cand.gctEmCand()->phiIndex();
0354
0355       // Use the L1CaloGeometry to find the eta, phi bin boundaries.
0356       double etaLow = l1CaloGeom.etaBinLowEdge(etaIndex);
0357       double etaHigh = l1CaloGeom.etaBinHighEdge(etaIndex);
0358       double phiLow = l1CaloGeom.emJetPhiBinLowEdge(phiIndex);
0359       double phiHigh = l1CaloGeom.emJetPhiBinHighEdge(phiIndex);
0360
0361       etaLow -= regionEtaMargin_;
0362       etaHigh += regionEtaMargin_;
0363       phiLow -= regionPhiMargin_;
0364       phiHigh += regionPhiMargin_;
0365
0366       regions.push_back(RectangularEtaPhiRegion(etaLow, etaHigh, phiLow, phiHigh));
0367     }
0368   }
0369 }
0370
0371 typedef HLTRecHitInAllL1RegionsProducer<EcalRecHit> HLTEcalRecHitInAllL1RegionsProducer;
0372 DEFINE_FWK_MODULE(HLTEcalRecHitInAllL1RegionsProducer);
0373
0374 typedef HLTRecHitInAllL1RegionsProducer<EcalUncalibratedRecHit> HLTEcalUncalibratedRecHitInAllL1RegionsProducer;
0375 DEFINE_FWK_MODULE(HLTEcalUncalibratedRecHitInAllL1RegionsProducer);