File indexing completed on 2023-03-17 11:23:54
0001 #include "TSToSCAssociatorByEnergyScoreImpl.h"
0002
0003 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0004 #include "SimDataFormats/CaloAnalysis/interface/SimCluster.h"
0005
0006 TSToSCAssociatorByEnergyScoreImpl::TSToSCAssociatorByEnergyScoreImpl(
0007 edm::EDProductGetter const& productGetter,
0008 bool hardScatterOnly,
0009 std::shared_ptr<hgcal::RecHitTools> recHitTools,
0010 const std::unordered_map<DetId, const HGCRecHit*>* hitMap)
0011 : hardScatterOnly_(hardScatterOnly), recHitTools_(recHitTools), hitMap_(hitMap), productGetter_(&productGetter) {
0012 layers_ = recHitTools_->lastLayerBH();
0013 }
0014
0015 hgcal::association TSToSCAssociatorByEnergyScoreImpl::makeConnections(
0016 const edm::Handle<ticl::TracksterCollection>& tCH,
0017 const edm::Handle<reco::CaloClusterCollection>& lCCH,
0018 const edm::Handle<SimClusterCollection>& sCCH) const {
0019
0020 const auto& tracksters = *tCH.product();
0021 const auto& layerClusters = *lCCH.product();
0022 const auto& simClusters = *sCCH.product();
0023 auto nTracksters = tracksters.size();
0024
0025
0026
0027 auto nSimClusters = simClusters.size();
0028 std::vector<size_t> sCIndices;
0029 for (unsigned int scId = 0; scId < nSimClusters; ++scId) {
0030 if (hardScatterOnly_ && (simClusters[scId].g4Tracks()[0].eventId().event() != 0 or
0031 simClusters[scId].g4Tracks()[0].eventId().bunchCrossing() != 0)) {
0032 LogDebug("TSToSCAssociatorByEnergyScoreImpl")
0033 << "Excluding SimCluster from event: " << simClusters[scId].g4Tracks()[0].eventId().event()
0034 << " with BX: " << simClusters[scId].g4Tracks()[0].eventId().bunchCrossing() << std::endl;
0035 continue;
0036 }
0037 sCIndices.emplace_back(scId);
0038 }
0039 nSimClusters = sCIndices.size();
0040
0041
0042
0043
0044 hgcal::simClusterToTrackster tssInSimCluster;
0045 tssInSimCluster.resize(nSimClusters);
0046 for (unsigned int i = 0; i < nSimClusters; ++i) {
0047 tssInSimCluster[i].simClusterId = i;
0048 tssInSimCluster[i].energy = 0.f;
0049 tssInSimCluster[i].hits_and_fractions.clear();
0050 }
0051
0052
0053 std::unordered_map<DetId, std::vector<hgcal::detIdInfoInCluster>> detIdToSimClusterId_Map;
0054 for (const auto& scId : sCIndices) {
0055 const auto& hits_and_fractions = simClusters[scId].hits_and_fractions();
0056 for (const auto& it_haf : hits_and_fractions) {
0057 const auto hitid = it_haf.first;
0058 const auto itcheck = hitMap_->find(hitid);
0059 if (itcheck != hitMap_->end()) {
0060 const auto hit_find_it = detIdToSimClusterId_Map.find(hitid);
0061 if (hit_find_it == detIdToSimClusterId_Map.end()) {
0062 detIdToSimClusterId_Map[hitid] = std::vector<hgcal::detIdInfoInCluster>();
0063 }
0064 detIdToSimClusterId_Map[hitid].emplace_back(scId, it_haf.second);
0065
0066 const HGCRecHit* hit = itcheck->second;
0067 tssInSimCluster[scId].energy += it_haf.second * hit->energy();
0068 tssInSimCluster[scId].hits_and_fractions.emplace_back(hitid, it_haf.second);
0069 }
0070 }
0071 }
0072
0073 #ifdef EDM_ML_DEBUG
0074 LogDebug("TSToSCAssociatorByEnergyScoreImpl")
0075 << "tssInSimCluster INFO (Only SimCluster filled at the moment)" << std::endl;
0076 for (size_t sc = 0; sc < tssInSimCluster.size(); ++sc) {
0077 LogDebug("TSToSCAssociatorByEnergyScoreImpl") << "For SimCluster Idx: " << sc << " we have: " << std::endl;
0078 LogDebug("TSToSCAssociatorByEnergyScoreImpl")
0079 << "\tSimClusterIdx:\t" << tssInSimCluster[sc].simClusterId << std::endl;
0080 LogDebug("TSToSCAssociatorByEnergyScoreImpl") << "\tEnergy:\t" << tssInSimCluster[sc].energy << std::endl;
0081 LogDebug("TSToSCAssociatorByEnergyScoreImpl") << "\t# of clusters:\t" << layerClusters.size() << std::endl;
0082 double tot_energy = 0.;
0083 for (auto const& haf : tssInSimCluster[sc].hits_and_fractions) {
0084 LogDebug("TSToSCAssociatorByEnergyScoreImpl")
0085 << "\tHits/fraction/energy: " << (uint32_t)haf.first << "/" << haf.second << "/"
0086 << haf.second * hitMap_->at(haf.first)->energy() << std::endl;
0087 tot_energy += haf.second * hitMap_->at(haf.first)->energy();
0088 }
0089 LogDebug("TSToSCAssociatorByEnergyScoreImpl") << "\tTot Sum haf: " << tot_energy << std::endl;
0090 for (auto const& ts : tssInSimCluster[sc].tracksterIdToEnergyAndScore) {
0091 LogDebug("TSToSCAssociatorByEnergyScoreImpl")
0092 << "\ttsIdx/energy/score: " << ts.first << "/" << ts.second.first << "/" << ts.second.second << std::endl;
0093 }
0094 }
0095
0096 LogDebug("TSToSCAssociatorByEnergyScoreImpl") << "detIdToSimClusterId_Map INFO" << std::endl;
0097 for (auto const& detId : detIdToSimClusterId_Map) {
0098 LogDebug("TSToSCAssociatorByEnergyScoreImpl")
0099 << "For detId: " << (uint32_t)detId.first
0100 << " we have found the following connections with SimClusters:" << std::endl;
0101 for (auto const& sc : detId.second) {
0102 LogDebug("TSToSCAssociatorByEnergyScoreImpl")
0103 << "\tSimCluster Id: " << sc.clusterId << " with fraction: " << sc.fraction
0104 << " and energy: " << sc.fraction * hitMap_->at(detId.first)->energy() << std::endl;
0105 }
0106 }
0107 #endif
0108
0109
0110 std::unordered_map<DetId, std::vector<hgcal::detIdInfoInCluster>> detIdToLayerClusterId_Map;
0111
0112
0113
0114 hgcal::tracksterToSimCluster scsInTrackster;
0115 scsInTrackster.resize(nTracksters);
0116
0117 for (unsigned int tsId = 0; tsId < nTracksters; ++tsId) {
0118 for (unsigned int i = 0; i < tracksters[tsId].vertices().size(); ++i) {
0119 const auto lcId = tracksters[tsId].vertices(i);
0120 const auto lcFractionInTs = 1.f / tracksters[tsId].vertex_multiplicity(i);
0121
0122 const std::vector<std::pair<DetId, float>>& hits_and_fractions = layerClusters[lcId].hitsAndFractions();
0123 unsigned int numberOfHitsInLC = hits_and_fractions.size();
0124
0125 for (unsigned int hitId = 0; hitId < numberOfHitsInLC; hitId++) {
0126 const auto rh_detid = hits_and_fractions[hitId].first;
0127 const auto rhFraction = hits_and_fractions[hitId].second;
0128
0129 const auto hit_find_in_LC = detIdToLayerClusterId_Map.find(rh_detid);
0130 if (hit_find_in_LC == detIdToLayerClusterId_Map.end()) {
0131 detIdToLayerClusterId_Map[rh_detid] = std::vector<hgcal::detIdInfoInCluster>();
0132 }
0133 detIdToLayerClusterId_Map[rh_detid].emplace_back(lcId, rhFraction);
0134
0135 const auto hit_find_in_SC = detIdToSimClusterId_Map.find(rh_detid);
0136
0137 if (hit_find_in_SC != detIdToSimClusterId_Map.end()) {
0138 const auto itcheck = hitMap_->find(rh_detid);
0139 const HGCRecHit* hit = itcheck->second;
0140
0141
0142
0143 for (const auto& h : hit_find_in_SC->second) {
0144
0145
0146 tssInSimCluster[h.clusterId].tracksterIdToEnergyAndScore[tsId].first +=
0147 lcFractionInTs * h.fraction * hit->energy();
0148
0149 scsInTrackster[tsId].emplace_back(h.clusterId, 0.f);
0150 }
0151 }
0152 }
0153 }
0154 }
0155
0156 #ifdef EDM_ML_DEBUG
0157 for (unsigned int tsId = 0; tsId < nTracksters; ++tsId) {
0158 for (const auto& lcId : tracksters[tsId].vertices()) {
0159 const auto& hits_and_fractions = layerClusters[lcId].hitsAndFractions();
0160 unsigned int numberOfHitsInLC = hits_and_fractions.size();
0161
0162
0163
0164
0165
0166
0167
0168
0169
0170 std::vector<int> hitsToSimClusterId(numberOfHitsInLC);
0171
0172 int maxSCId_byNumberOfHits = -1;
0173
0174 unsigned int maxSCNumberOfHitsInLC = 0;
0175
0176 int maxSCId_byEnergy = -1;
0177
0178 float maxEnergySharedLCandSC = 0.f;
0179
0180 float energyFractionOfLCinSC = 0.f;
0181
0182 float energyFractionOfSCinLC = 0.f;
0183 std::unordered_map<unsigned, unsigned> occurrencesSCinLC;
0184 unsigned int numberOfNoiseHitsInLC = 0;
0185 std::unordered_map<unsigned, float> SCEnergyInLC;
0186
0187 for (unsigned int hitId = 0; hitId < numberOfHitsInLC; hitId++) {
0188 const auto rh_detid = hits_and_fractions[hitId].first;
0189 const auto rhFraction = hits_and_fractions[hitId].second;
0190
0191 const auto hit_find_in_SC = detIdToSimClusterId_Map.find(rh_detid);
0192
0193
0194
0195
0196
0197
0198
0199 if (rhFraction == 0.) {
0200 hitsToSimClusterId[hitId] = -2;
0201 }
0202
0203 if (hit_find_in_SC == detIdToSimClusterId_Map.end()) {
0204 hitsToSimClusterId[hitId] -= 1;
0205 } else {
0206 const auto itcheck = hitMap_->find(rh_detid);
0207 const HGCRecHit* hit = itcheck->second;
0208 auto maxSCEnergyInLC = 0.f;
0209 auto maxSCId = -1;
0210
0211 for (const auto& h : hit_find_in_SC->second) {
0212 SCEnergyInLC[h.clusterId] += h.fraction * hit->energy();
0213
0214
0215 if (SCEnergyInLC[h.clusterId] > maxSCEnergyInLC) {
0216 maxSCEnergyInLC = SCEnergyInLC[h.clusterId];
0217 maxSCId = h.clusterId;
0218 }
0219 }
0220 hitsToSimClusterId[hitId] = maxSCId;
0221 }
0222 }
0223
0224 for (const auto& c : hitsToSimClusterId) {
0225 if (c < 0) {
0226 numberOfNoiseHitsInLC++;
0227 } else {
0228 occurrencesSCinLC[c]++;
0229 }
0230 }
0231
0232 for (const auto& c : occurrencesSCinLC) {
0233 if (c.second > maxSCNumberOfHitsInLC) {
0234 maxSCId_byNumberOfHits = c.first;
0235 maxSCNumberOfHitsInLC = c.second;
0236 }
0237 }
0238
0239 for (const auto& c : SCEnergyInLC) {
0240 if (c.second > maxEnergySharedLCandSC) {
0241 maxSCId_byEnergy = c.first;
0242 maxEnergySharedLCandSC = c.second;
0243 }
0244 }
0245
0246 float totalSCEnergyOnLayer = 0.f;
0247 if (maxSCId_byEnergy >= 0) {
0248 totalSCEnergyOnLayer = tssInSimCluster[maxSCId_byEnergy].energy;
0249 energyFractionOfSCinLC = maxEnergySharedLCandSC / totalSCEnergyOnLayer;
0250 if (tracksters[tsId].raw_energy() > 0.f) {
0251 energyFractionOfLCinSC = maxEnergySharedLCandSC / tracksters[tsId].raw_energy();
0252 }
0253 }
0254
0255 LogDebug("TSToSCAssociatorByEnergyScoreImpl")
0256 << std::setw(12) << "TracksterID:\t" << std::setw(12) << "layerCluster\t" << std::setw(10) << "lc energy\t"
0257 << std::setw(5) << "nhits\t" << std::setw(12) << "noise hits\t" << std::setw(22) << "maxSCId_byNumberOfHits\t"
0258 << std::setw(8) << "nhitsSC\t" << std::setw(13) << "maxSCId_byEnergy\t" << std::setw(20)
0259 << "maxEnergySharedLCandSC\t" << std::setw(22) << "totalSCEnergyOnLayer\t" << std::setw(22)
0260 << "energyFractionOfLCinSC\t" << std::setw(25) << "energyFractionOfSCinLC\t"
0261 << "\n";
0262 LogDebug("TSToSCAssociatorByEnergyScoreImpl")
0263 << std::setw(12) << tsId << "\t" << std::setw(12) << lcId << "\t" << std::setw(10)
0264 << tracksters[tsId].raw_energy() << "\t" << std::setw(5) << numberOfHitsInLC << "\t" << std::setw(12)
0265 << numberOfNoiseHitsInLC << "\t" << std::setw(22) << maxSCId_byNumberOfHits << "\t" << std::setw(8)
0266 << maxSCNumberOfHitsInLC << "\t" << std::setw(13) << maxSCId_byEnergy << "\t" << std::setw(20)
0267 << maxEnergySharedLCandSC << "\t" << std::setw(22) << totalSCEnergyOnLayer << "\t" << std::setw(22)
0268 << energyFractionOfLCinSC << "\t" << std::setw(25) << energyFractionOfSCinLC << "\n";
0269 }
0270 }
0271
0272 LogDebug("TSToSCAssociatorByEnergyScoreImpl")
0273 << "Improved tssInSimCluster INFO (Now containing the linked tracksters id and energy - score still empty)"
0274 << std::endl;
0275 for (size_t sc = 0; sc < tssInSimCluster.size(); ++sc) {
0276 LogDebug("TSToSCAssociatorByEnergyScoreImpl") << "For SimCluster Idx: " << sc << " we have: " << std::endl;
0277 LogDebug("TSToSCAssociatorByEnergyScoreImpl")
0278 << " SimClusterIdx: " << tssInSimCluster[sc].simClusterId << std::endl;
0279 LogDebug("TSToSCAssociatorByEnergyScoreImpl") << "\tEnergy:\t" << tssInSimCluster[sc].energy << std::endl;
0280 double tot_energy = 0.;
0281 for (auto const& haf : tssInSimCluster[sc].hits_and_fractions) {
0282 LogDebug("TSToSCAssociatorByEnergyScoreImpl")
0283 << "\tHits/fraction/energy: " << (uint32_t)haf.first << "/" << haf.second << "/"
0284 << haf.second * hitMap_->at(haf.first)->energy() << std::endl;
0285 tot_energy += haf.second * hitMap_->at(haf.first)->energy();
0286 }
0287 LogDebug("TSToSCAssociatorByEnergyScoreImpl") << "\tTot Sum haf: " << tot_energy << std::endl;
0288 for (auto const& ts : tssInSimCluster[sc].tracksterIdToEnergyAndScore) {
0289 LogDebug("TSToSCAssociatorByEnergyScoreImpl")
0290 << "\ttsIdx/energy/score: " << ts.first << "/" << ts.second.first << "/" << ts.second.second << std::endl;
0291 }
0292 }
0293
0294 LogDebug("TSToSCAssociatorByEnergyScoreImpl") << "Improved detIdToSimClusterId_Map INFO" << std::endl;
0295 for (auto const& sc : detIdToSimClusterId_Map) {
0296 LogDebug("TSToSCAssociatorByEnergyScoreImpl")
0297 << "For detId: " << (uint32_t)sc.first
0298 << " we have found the following connections with SimClusters:" << std::endl;
0299 for (auto const& sclu : sc.second) {
0300 LogDebug("TSToSCAssociatorByEnergyScoreImpl")
0301 << " SimCluster Id: " << sclu.clusterId << " with fraction: " << sclu.fraction
0302 << " and energy: " << sclu.fraction * hitMap_->at(sc.first)->energy() << std::endl;
0303 }
0304 }
0305 #endif
0306
0307
0308
0309 for (unsigned int tsId = 0; tsId < nTracksters; ++tsId) {
0310
0311
0312 std::sort(scsInTrackster[tsId].begin(), scsInTrackster[tsId].end());
0313 auto last = std::unique(scsInTrackster[tsId].begin(), scsInTrackster[tsId].end());
0314 scsInTrackster[tsId].erase(last, scsInTrackster[tsId].end());
0315
0316
0317
0318 if (tracksters[tsId].raw_energy() == 0. && !scsInTrackster[tsId].empty()) {
0319 for (auto& scPair : scsInTrackster[tsId]) {
0320 scPair.second = 1.;
0321 LogDebug("TSToSCAssociatorByEnergyScoreImpl")
0322 << "TracksterId:\t " << tsId << "\tSC id:\t" << scPair.first << "\tscore\t " << scPair.second << "\n";
0323 }
0324 continue;
0325 }
0326
0327 float invTracksterEnergyWeight = 0.f;
0328 for (unsigned int i = 0; i < tracksters[tsId].vertices().size(); ++i) {
0329 const auto lcId = tracksters[tsId].vertices(i);
0330 const auto lcFractionInTs = 1.f / tracksters[tsId].vertex_multiplicity(i);
0331
0332 const auto& hits_and_fractions = layerClusters[lcId].hitsAndFractions();
0333
0334 for (auto const& haf : hits_and_fractions) {
0335 invTracksterEnergyWeight += (lcFractionInTs * haf.second * hitMap_->at(haf.first)->energy()) *
0336 (lcFractionInTs * haf.second * hitMap_->at(haf.first)->energy());
0337 }
0338 }
0339 invTracksterEnergyWeight = 1.f / invTracksterEnergyWeight;
0340
0341 for (unsigned int i = 0; i < tracksters[tsId].vertices().size(); ++i) {
0342 const auto lcId = tracksters[tsId].vertices(i);
0343 const auto lcFractionInTs = 1.f / tracksters[tsId].vertex_multiplicity(i);
0344
0345 const auto& hits_and_fractions = layerClusters[lcId].hitsAndFractions();
0346 unsigned int numberOfHitsInLC = hits_and_fractions.size();
0347 for (unsigned int i = 0; i < numberOfHitsInLC; ++i) {
0348 DetId rh_detid = hits_and_fractions[i].first;
0349 float rhFraction = hits_and_fractions[i].second * lcFractionInTs;
0350
0351 const bool hitWithSC = (detIdToSimClusterId_Map.find(rh_detid) != detIdToSimClusterId_Map.end());
0352
0353 const auto itcheck = hitMap_->find(rh_detid);
0354 const HGCRecHit* hit = itcheck->second;
0355 float hitEnergyWeight = hit->energy() * hit->energy();
0356
0357 for (auto& scPair : scsInTrackster[tsId]) {
0358 float scFraction = 0.f;
0359 if (hitWithSC) {
0360 const auto findHitIt = std::find(detIdToSimClusterId_Map[rh_detid].begin(),
0361 detIdToSimClusterId_Map[rh_detid].end(),
0362 hgcal::detIdInfoInCluster{scPair.first, 0.f});
0363 if (findHitIt != detIdToSimClusterId_Map[rh_detid].end())
0364 scFraction = findHitIt->fraction;
0365 }
0366 scPair.second +=
0367 (rhFraction - scFraction) * (rhFraction - scFraction) * hitEnergyWeight * invTracksterEnergyWeight;
0368 #ifdef EDM_ML_DEBUG
0369 LogDebug("TSToSCAssociatorByEnergyScoreImpl")
0370 << "rh_detid:\t" << (uint32_t)rh_detid << "\ttracksterId:\t" << tsId << "\t"
0371 << "rhfraction,scFraction:\t" << rhFraction << ", " << scFraction << "\t"
0372 << "hitEnergyWeight:\t" << hitEnergyWeight << "\t"
0373 << "current score:\t" << scPair.second << "\t"
0374 << "invTracksterEnergyWeight:\t" << invTracksterEnergyWeight << "\n";
0375 #endif
0376 }
0377 }
0378 }
0379
0380 #ifdef EDM_ML_DEBUG
0381 if (scsInTrackster[tsId].empty())
0382 LogDebug("TSToSCAssociatorByEnergyScoreImpl") << "trackster Id:\t" << tsId << "\tSC id:\t-1"
0383 << "\tscore\t-1\n";
0384 #endif
0385 }
0386
0387
0388 for (const auto& scId : sCIndices) {
0389 float invSCEnergyWeight = 0.f;
0390
0391 const unsigned int SCNumberOfHits = tssInSimCluster[scId].hits_and_fractions.size();
0392 if (SCNumberOfHits == 0)
0393 continue;
0394 #ifdef EDM_ML_DEBUG
0395 int tsWithMaxEnergyInSC = -1;
0396
0397 float maxEnergyTSinSC = 0.f;
0398 float SCenergy = tssInSimCluster[scId].energy;
0399
0400 float SCEnergyFractionInTS = 0.f;
0401 for (const auto& ts : tssInSimCluster[scId].tracksterIdToEnergyAndScore) {
0402 if (ts.second.first > maxEnergyTSinSC) {
0403 maxEnergyTSinSC = ts.second.first;
0404 tsWithMaxEnergyInSC = ts.first;
0405 }
0406 }
0407 if (SCenergy > 0.f)
0408 SCEnergyFractionInTS = maxEnergyTSinSC / SCenergy;
0409
0410 LogDebug("TSToSCAssociatorByEnergyScoreImpl")
0411 << std::setw(12) << "simcluster\t" << std::setw(15) << "sc total energy\t" << std::setw(15)
0412 << "scEnergyOnLayer\t" << std::setw(14) << "SCNhitsOnLayer\t" << std::setw(18) << "tsWithMaxEnergyInSC\t"
0413 << std::setw(15) << "maxEnergyTSinSC\t" << std::setw(20) << "SCEnergyFractionInTS"
0414 << "\n";
0415 LogDebug("TSToSCAssociatorByEnergyScoreImpl")
0416 << std::setw(12) << scId << "\t" << std::setw(15) << simClusters[scId].energy() << "\t" << std::setw(15)
0417 << SCenergy << "\t" << std::setw(14) << SCNumberOfHits << "\t" << std::setw(18) << tsWithMaxEnergyInSC << "\t"
0418 << std::setw(15) << maxEnergyTSinSC << "\t" << std::setw(20) << SCEnergyFractionInTS << "\n";
0419 #endif
0420
0421 for (auto const& haf : tssInSimCluster[scId].hits_and_fractions) {
0422 invSCEnergyWeight += std::pow(haf.second * hitMap_->at(haf.first)->energy(), 2);
0423 }
0424 invSCEnergyWeight = 1.f / invSCEnergyWeight;
0425
0426 for (unsigned int i = 0; i < SCNumberOfHits; ++i) {
0427 auto& sc_hitDetId = tssInSimCluster[scId].hits_and_fractions[i].first;
0428 auto& scFraction = tssInSimCluster[scId].hits_and_fractions[i].second;
0429
0430 bool hitWithLC = false;
0431 if (scFraction == 0.f)
0432 continue;
0433 const auto hit_find_in_LC = detIdToLayerClusterId_Map.find(sc_hitDetId);
0434 if (hit_find_in_LC != detIdToLayerClusterId_Map.end())
0435 hitWithLC = true;
0436 const auto itcheck = hitMap_->find(sc_hitDetId);
0437 const HGCRecHit* hit = itcheck->second;
0438 float hitEnergyWeight = hit->energy() * hit->energy();
0439 for (auto& tsPair : tssInSimCluster[scId].tracksterIdToEnergyAndScore) {
0440 unsigned int tsId = tsPair.first;
0441 float tsFraction = 0.f;
0442
0443 for (unsigned int i = 0; i < tracksters[tsId].vertices().size(); ++i) {
0444 const auto lcId = tracksters[tsId].vertices(i);
0445 const auto lcFractionInTs = 1.f / tracksters[tsId].vertex_multiplicity(i);
0446
0447 if (hitWithLC) {
0448 const auto findHitIt = std::find(detIdToLayerClusterId_Map[sc_hitDetId].begin(),
0449 detIdToLayerClusterId_Map[sc_hitDetId].end(),
0450 hgcal::detIdInfoInCluster{lcId, 0.f});
0451 if (findHitIt != detIdToLayerClusterId_Map[sc_hitDetId].end())
0452 tsFraction = findHitIt->fraction * lcFractionInTs;
0453 }
0454 tsPair.second.second +=
0455 (tsFraction - scFraction) * (tsFraction - scFraction) * hitEnergyWeight * invSCEnergyWeight;
0456 #ifdef EDM_ML_DEBUG
0457 LogDebug("TSToSCAssociatorByEnergyScoreImpl")
0458 << "SCDetId:\t" << (uint32_t)sc_hitDetId << "\tTracksterId:\t" << tsId << "\t"
0459 << "tsFraction, scFraction:\t" << tsFraction << ", " << scFraction << "\t"
0460 << "hitEnergyWeight:\t" << hitEnergyWeight << "\t"
0461 << "current score:\t" << tsPair.second.second << "\t"
0462 << "invSCEnergyWeight:\t" << invSCEnergyWeight << "\n";
0463 #endif
0464 }
0465 }
0466 }
0467 #ifdef EDM_ML_DEBUG
0468 if (tssInSimCluster[scId].tracksterIdToEnergyAndScore.empty())
0469 LogDebug("TSToSCAssociatorByEnergyScoreImpl") << "SC Id:\t" << scId << "\tTS id:\t-1 "
0470 << "\tscore\t-1\n";
0471
0472 for (const auto& tsPair : tssInSimCluster[scId].tracksterIdToEnergyAndScore) {
0473 LogDebug("TSToSCAssociatorByEnergyScoreImpl")
0474 << "SC Id: \t" << scId << "\t TS id: \t" << tsPair.first << "\t score \t" << tsPair.second.second
0475 << "\t shared energy:\t" << tsPair.second.first << "\t shared energy fraction:\t"
0476 << (tsPair.second.first / SCenergy) << "\n";
0477 }
0478 #endif
0479 }
0480 return {scsInTrackster, tssInSimCluster};
0481 }
0482
0483 hgcal::RecoToSimCollectionTracksters TSToSCAssociatorByEnergyScoreImpl::associateRecoToSim(
0484 const edm::Handle<ticl::TracksterCollection>& tCH,
0485 const edm::Handle<reco::CaloClusterCollection>& lCCH,
0486 const edm::Handle<SimClusterCollection>& sCCH) const {
0487 hgcal::RecoToSimCollectionTracksters returnValue(productGetter_);
0488 const auto& links = makeConnections(tCH, lCCH, sCCH);
0489
0490 const auto& scsInTrackster = std::get<0>(links);
0491 for (size_t tsId = 0; tsId < scsInTrackster.size(); ++tsId) {
0492 for (auto& scPair : scsInTrackster[tsId]) {
0493 LogDebug("TSToSCAssociatorByEnergyScoreImpl")
0494 << "Trackster Id:\t" << tsId << "\tSimCluster id:\t" << scPair.first << "\tscore:\t" << scPair.second << "\n";
0495
0496 returnValue.insert(edm::Ref<ticl::TracksterCollection>(tCH, tsId),
0497 std::make_pair(edm::Ref<SimClusterCollection>(sCCH, scPair.first),
0498 scPair.second)
0499 );
0500 }
0501 }
0502 return returnValue;
0503 }
0504
0505 hgcal::SimToRecoCollectionTracksters TSToSCAssociatorByEnergyScoreImpl::associateSimToReco(
0506 const edm::Handle<ticl::TracksterCollection>& tCH,
0507 const edm::Handle<reco::CaloClusterCollection>& lCCH,
0508 const edm::Handle<SimClusterCollection>& sCCH) const {
0509 hgcal::SimToRecoCollectionTracksters returnValue(productGetter_);
0510 const auto& links = makeConnections(tCH, lCCH, sCCH);
0511 const auto& tssInSimCluster = std::get<1>(links);
0512 for (size_t scId = 0; scId < tssInSimCluster.size(); ++scId) {
0513 for (auto& tsPair : tssInSimCluster[scId].tracksterIdToEnergyAndScore) {
0514 returnValue.insert(
0515 edm::Ref<SimClusterCollection>(sCCH, scId),
0516 std::make_pair(edm::Ref<ticl::TracksterCollection>(tCH, tsPair.first),
0517 std::make_pair(tsPair.second.first, tsPair.second.second))
0518 );
0519 }
0520 }
0521 return returnValue;
0522 }