CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
HitToTracksterAssociatorProducer.cc
Go to the documentation of this file.
1 // Author: Felice Pantaleo, felice.pantaleo@cern.ch 06/2024
2 
3 // user include files
17 
19  : LCCollectionToken_(consumes<std::vector<reco::CaloCluster>>(pset.getParameter<edm::InputTag>("layer_clusters"))),
20  tracksterCollectionToken_(consumes<std::vector<ticl::Trackster>>(pset.getParameter<edm::InputTag>("tracksters"))),
21  hitMapToken_(
22  consumes<std::unordered_map<DetId, const unsigned int>>(pset.getParameter<edm::InputTag>("hitMapTag"))) {
23  auto hitsTags = pset.getParameter<std::vector<edm::InputTag>>("hits");
24  for (const auto &tag : hitsTags) {
25  hitsTokens_.push_back(consumes<HGCRecHitCollection>(tag));
26  }
27  produces<ticl::AssociationMap<ticl::mapWithFraction>>("hitToTracksterMap");
28  produces<ticl::AssociationMap<ticl::mapWithFraction>>("tracksterToHitMap");
29 }
30 
32 
34  using namespace edm;
35 
38 
41 
43  iEvent.getByToken(hitMapToken_, hitMap);
44 
45  MultiVectorManager<HGCRecHit> rechitManager;
46  for (const auto &token : hitsTokens_) {
47  Handle<HGCRecHitCollection> hitsHandle;
48  iEvent.getByToken(token, hitsHandle);
49  rechitManager.addVector(*hitsHandle);
50  }
51 
52  // Create association map
53  auto hitToTracksterMap = std::make_unique<ticl::AssociationMap<ticl::mapWithFraction>>(rechitManager.size());
54  auto tracksterToHitMap = std::make_unique<ticl::AssociationMap<ticl::mapWithFraction>>(tracksters->size());
55 
56  // Loop over tracksters
57  for (unsigned int tracksterId = 0; tracksterId < tracksters->size(); ++tracksterId) {
58  const auto &trackster = (*tracksters)[tracksterId];
59  // Loop over vertices in trackster
60  for (unsigned int i = 0; i < trackster.vertices().size(); ++i) {
61  // Get layerCluster
62  const auto &lc = (*layer_clusters)[trackster.vertices()[i]];
63  float invMultiplicity = 1.0f / trackster.vertex_multiplicity()[i];
64 
65  for (const auto &hitAndFraction : lc.hitsAndFractions()) {
66  auto hitMapIter = hitMap->find(hitAndFraction.first);
67  if (hitMapIter != hitMap->end()) {
68  unsigned int rechitIndex = hitMapIter->second;
69  float fraction = hitAndFraction.second * invMultiplicity;
70  hitToTracksterMap->insert(rechitIndex, tracksterId, fraction);
71  tracksterToHitMap->insert(tracksterId, rechitIndex, fraction);
72  }
73  }
74  }
75  }
76  iEvent.put(std::move(hitToTracksterMap), "hitToTracksterMap");
77  iEvent.put(std::move(tracksterToHitMap), "tracksterToHitMap");
78 }
79 
82  desc.add<edm::InputTag>("layer_clusters", edm::InputTag("hgcalMergeLayerClusters"));
83  desc.add<edm::InputTag>("tracksters", edm::InputTag("ticlTracksters"));
84  desc.add<edm::InputTag>("hitMapTag", edm::InputTag("recHitMapProducer", "hgcalRecHitMap"));
85  desc.add<std::vector<edm::InputTag>>("hits",
86  {edm::InputTag("HGCalRecHit", "HGCEERecHits"),
87  edm::InputTag("HGCalRecHit", "HGCHEFRecHits"),
88  edm::InputTag("HGCalRecHit", "HGCHEBRecHits")});
89  descriptions.add("hitToTracksterAssociator", desc);
90 }
91 
92 // Define this as a plug-in
edm::EDGetTokenT< std::unordered_map< DetId, const unsigned int > > hitMapToken_
edm::EDGetTokenT< std::vector< reco::CaloCluster > > LCCollectionToken_
std::vector< edm::EDGetTokenT< HGCRecHitCollection > > hitsTokens_
int iEvent
Definition: GenABIO.cc:224
HitToTracksterAssociatorProducer(const edm::ParameterSet &)
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
void produce(edm::StreamID, edm::Event &, const edm::EventSetup &) const override
Definition: DetId.h:17
void add(std::string const &label, ParameterSetDescription const &psetDescription)
void addVector(std::span< const T > vec)
edm::EDGetTokenT< std::vector< ticl::Trackster > > tracksterCollectionToken_
fixed size matrix
HLT enums.
Definition: Common.h:10
def move(src, dest)
Definition: eostools.py:511