@@ -23,35 +23,37 @@ class RecHitMapProducer : public edm::global::EDProducer<> {
2323 void produce (edm::StreamID, edm::Event&, const edm::EventSetup&) const override ;
2424
2525private:
26- const edm::EDGetTokenT<HGCRecHitCollection> hits_ee_token_;
27- const edm::EDGetTokenT<HGCRecHitCollection> hits_fh_token_;
28- const edm::EDGetTokenT<HGCRecHitCollection> hits_bh_token_;
29- const edm::EDGetTokenT<reco::PFRecHitCollection> hits_eb_token_;
30- const edm::EDGetTokenT<reco::PFRecHitCollection> hits_hb_token_;
31- const edm::EDGetTokenT<reco::PFRecHitCollection> hits_ho_token_;
26+ std::vector<edm::EDGetTokenT<HGCRecHitCollection>> hgcal_hits_token_;
27+ std::vector<edm::EDGetTokenT<reco::PFRecHitCollection>> barrel_hits_token_;
28+
29+ bool hgcalOnly_;
3230};
3331
3432DEFINE_FWK_MODULE (RecHitMapProducer);
3533
3634using DetIdRecHitMap = std::unordered_map<DetId, const unsigned int >;
3735
38- RecHitMapProducer::RecHitMapProducer (const edm::ParameterSet& ps)
39- : hits_ee_token_(consumes<HGCRecHitCollection>(ps.getParameter<edm::InputTag>(" EEInput" ))),
40- hits_fh_token_(consumes<HGCRecHitCollection>(ps.getParameter<edm::InputTag>(" FHInput" ))),
41- hits_bh_token_(consumes<HGCRecHitCollection>(ps.getParameter<edm::InputTag>(" BHInput" ))),
42- hits_eb_token_(consumes<reco::PFRecHitCollection>(ps.getParameter<edm::InputTag>(" EBInput" ))),
43- hits_hb_token_(consumes<reco::PFRecHitCollection>(ps.getParameter<edm::InputTag>(" HBInput" ))) {
36+ RecHitMapProducer::RecHitMapProducer (const edm::ParameterSet& ps) : hgcalOnly_(ps.getParameter<bool >(" hgcalOnly" )) {
37+ std::vector<edm::InputTag> tags = ps.getParameter <std::vector<edm::InputTag>>(" hits" );
38+ for (auto & tag : tags) {
39+ if (tag.label ().find (" HGCalRecHit" ) != std::string::npos) {
40+ hgcal_hits_token_.push_back (consumes<HGCRecHitCollection>(tag));
41+ } else {
42+ barrel_hits_token_.push_back (consumes<reco::PFRecHitCollection>(tag));
43+ }
44+ }
45+
4446 produces<DetIdRecHitMap>(" hgcalRecHitMap" );
45- produces<DetIdRecHitMap>(" barrelRecHitMap" );
47+ if (!hgcalOnly_)
48+ produces<DetIdRecHitMap>(" barrelRecHitMap" );
4649}
4750
4851void RecHitMapProducer::fillDescriptions (edm::ConfigurationDescriptions& descriptions) {
4952 edm::ParameterSetDescription desc;
50- desc.add <edm::InputTag>(" EEInput" , {" HGCalRecHit" , " HGCEERecHits" });
51- desc.add <edm::InputTag>(" FHInput" , {" HGCalRecHit" , " HGCHEFRecHits" });
52- desc.add <edm::InputTag>(" BHInput" , {" HGCalRecHit" , " HGCHEBRecHits" });
53- desc.add <edm::InputTag>(" EBInput" , {" particleFlowRecHitECAL" , " " });
54- desc.add <edm::InputTag>(" HBInput" , {" particleFlowRecHitHBHE" , " " });
53+ desc.add <std::vector<edm::InputTag>>(" hits" ,
54+ {edm::InputTag (" HGCalRecHit" , " HGCEERecHits" ),
55+ edm::InputTag (" HGCalRecHit" , " HGCHEFRecHits" ),
56+ edm::InputTag (" HGCalRecHit" , " HGCHEBRecHits" )});
5557 desc.add <bool >(" hgcalOnly" , true );
5658 descriptions.add (" recHitMapProducer" , desc);
5759}
@@ -60,9 +62,9 @@ void RecHitMapProducer::produce(edm::StreamID, edm::Event& evt, const edm::Event
6062 auto hitMapHGCal = std::make_unique<DetIdRecHitMap>();
6163
6264 // Retrieve collections
63- const auto & ee_hits = evt.getHandle (hits_ee_token_ );
64- const auto & fh_hits = evt.getHandle (hits_fh_token_ );
65- const auto & bh_hits = evt.getHandle (hits_bh_token_ );
65+ const auto & ee_hits = evt.getHandle (hgcal_hits_token_[ 0 ] );
66+ const auto & fh_hits = evt.getHandle (hgcal_hits_token_[ 1 ] );
67+ const auto & bh_hits = evt.getHandle (hgcal_hits_token_[ 2 ] );
6668
6769 // Check validity of all handles
6870 if (!ee_hits.isValid () || !fh_hits.isValid () || !bh_hits.isValid ()) {
@@ -71,6 +73,8 @@ void RecHitMapProducer::produce(edm::StreamID, edm::Event& evt, const edm::Event
7173 return ;
7274 }
7375
76+ // TODO may be worth to avoid dependency on the order
77+ // of the collections, maybe using a map
7478 MultiVectorManager<HGCRecHit> rechitManager;
7579 rechitManager.addVector (*ee_hits);
7680 rechitManager.addVector (*fh_hits);
@@ -83,13 +87,15 @@ void RecHitMapProducer::produce(edm::StreamID, edm::Event& evt, const edm::Event
8387
8488 evt.put (std::move (hitMapHGCal), " hgcalRecHitMap" );
8589
86- auto hitMapBarrel = std::make_unique<DetIdRecHitMap>();
87- MultiVectorManager<reco::PFRecHit> barrelRechitManager;
88- barrelRechitManager.addVector (evt.get (hits_eb_token_));
89- barrelRechitManager.addVector (evt.get (hits_hb_token_));
90- for (unsigned int i = 0 ; i < barrelRechitManager.size (); ++i) {
91- const auto recHitDetId = barrelRechitManager[i].detId ();
92- hitMapBarrel->emplace (recHitDetId, i);
90+ if (!hgcalOnly_) {
91+ auto hitMapBarrel = std::make_unique<DetIdRecHitMap>();
92+ MultiVectorManager<reco::PFRecHit> barrelRechitManager;
93+ barrelRechitManager.addVector (evt.get (barrel_hits_token_[0 ]));
94+ barrelRechitManager.addVector (evt.get (barrel_hits_token_[1 ]));
95+ for (unsigned int i = 0 ; i < barrelRechitManager.size (); ++i) {
96+ const auto recHitDetId = barrelRechitManager[i].detId ();
97+ hitMapBarrel->emplace (recHitDetId, i);
98+ }
99+ evt.put (std::move (hitMapBarrel), " barrelRecHitMap" );
93100 }
94- evt.put (std::move (hitMapBarrel), " barrelRecHitMap" );
95101}
0 commit comments