File indexing completed on 2024-10-16 05:06:33
0001
0002 #include <unordered_map>
0003
0004 #include "FWCore/Framework/interface/Frameworkfwd.h"
0005 #include "FWCore/Framework/interface/global/EDProducer.h"
0006 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0007
0008 #include "FWCore/Framework/interface/Event.h"
0009 #include "FWCore/Framework/interface/ESHandle.h"
0010 #include "FWCore/Framework/interface/MakerMacros.h"
0011 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0012 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0013
0014 #include "DataFormats/HGCRecHit/interface/HGCRecHitCollections.h"
0015 #include "DataFormats/ParticleFlowReco/interface/PFRecHit.h"
0016 #include "CommonTools/RecoAlgos/interface/MultiVectorManager.h"
0017
0018 class RecHitMapProducer : public edm::global::EDProducer<> {
0019 public:
0020 RecHitMapProducer(const edm::ParameterSet&);
0021 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0022
0023 void produce(edm::StreamID, edm::Event&, const edm::EventSetup&) const override;
0024
0025 private:
0026 const edm::EDGetTokenT<HGCRecHitCollection> hits_ee_token_;
0027 const edm::EDGetTokenT<HGCRecHitCollection> hits_fh_token_;
0028 const edm::EDGetTokenT<HGCRecHitCollection> hits_bh_token_;
0029 const edm::EDGetTokenT<reco::PFRecHitCollection> hits_eb_token_;
0030 const edm::EDGetTokenT<reco::PFRecHitCollection> hits_hb_token_;
0031 const edm::EDGetTokenT<reco::PFRecHitCollection> hits_ho_token_;
0032 bool hgcalOnly_;
0033 };
0034
0035 DEFINE_FWK_MODULE(RecHitMapProducer);
0036
0037 using DetIdRecHitMap = std::unordered_map<DetId, const unsigned int>;
0038
0039 RecHitMapProducer::RecHitMapProducer(const edm::ParameterSet& ps)
0040 : hits_ee_token_(consumes<HGCRecHitCollection>(ps.getParameter<edm::InputTag>("EEInput"))),
0041 hits_fh_token_(consumes<HGCRecHitCollection>(ps.getParameter<edm::InputTag>("FHInput"))),
0042 hits_bh_token_(consumes<HGCRecHitCollection>(ps.getParameter<edm::InputTag>("BHInput"))),
0043 hits_eb_token_(consumes<reco::PFRecHitCollection>(ps.getParameter<edm::InputTag>("EBInput"))),
0044 hits_hb_token_(consumes<reco::PFRecHitCollection>(ps.getParameter<edm::InputTag>("HBInput"))),
0045 hits_ho_token_(consumes<reco::PFRecHitCollection>(ps.getParameter<edm::InputTag>("HOInput"))),
0046 hgcalOnly_(ps.getParameter<bool>("hgcalOnly")) {
0047 produces<DetIdRecHitMap>("hgcalRecHitMap");
0048 if (!hgcalOnly_)
0049 produces<DetIdRecHitMap>("barrelRecHitMap");
0050 }
0051
0052 void RecHitMapProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0053 edm::ParameterSetDescription desc;
0054 desc.add<edm::InputTag>("EEInput", {"HGCalRecHit", "HGCEERecHits"});
0055 desc.add<edm::InputTag>("FHInput", {"HGCalRecHit", "HGCHEFRecHits"});
0056 desc.add<edm::InputTag>("BHInput", {"HGCalRecHit", "HGCHEBRecHits"});
0057 desc.add<edm::InputTag>("EBInput", {"particleFlowRecHitECAL", ""});
0058 desc.add<edm::InputTag>("HBInput", {"particleFlowRecHitHBHE", ""});
0059 desc.add<edm::InputTag>("HOInput", {"particleFlowRecHitHO", ""});
0060 desc.add<bool>("hgcalOnly", true);
0061 descriptions.add("recHitMapProducer", desc);
0062 }
0063
0064 void RecHitMapProducer::produce(edm::StreamID, edm::Event& evt, const edm::EventSetup& es) const {
0065 auto hitMapHGCal = std::make_unique<DetIdRecHitMap>();
0066 const auto& ee_hits = evt.get(hits_ee_token_);
0067 const auto& fh_hits = evt.get(hits_fh_token_);
0068 const auto& bh_hits = evt.get(hits_bh_token_);
0069
0070 MultiVectorManager<HGCRecHit> rechitManager;
0071 rechitManager.addVector(ee_hits);
0072 rechitManager.addVector(fh_hits);
0073 rechitManager.addVector(bh_hits);
0074
0075 for (unsigned int i = 0; i < rechitManager.size(); ++i) {
0076 const auto recHitDetId = rechitManager[i].detid();
0077 hitMapHGCal->emplace(recHitDetId, i);
0078 }
0079
0080 evt.put(std::move(hitMapHGCal), "hgcalRecHitMap");
0081
0082 if (!hgcalOnly_) {
0083 auto hitMapBarrel = std::make_unique<DetIdRecHitMap>();
0084 MultiVectorManager<reco::PFRecHit> barrelRechitManager;
0085 barrelRechitManager.addVector(evt.get(hits_eb_token_));
0086 barrelRechitManager.addVector(evt.get(hits_hb_token_));
0087 barrelRechitManager.addVector(evt.get(hits_ho_token_));
0088 for (unsigned int i = 0; i < barrelRechitManager.size(); ++i) {
0089 const auto recHitDetId = barrelRechitManager[i].detId();
0090 hitMapBarrel->emplace(recHitDetId, i);
0091 }
0092 evt.put(std::move(hitMapBarrel), "barrelRecHitMap");
0093 }
0094 }