File indexing completed on 2023-03-17 11:23:55
0001 #include "TSToSimTSAssociatorByEnergyScoreImpl.h"
0002
0003 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0004
0005 TSToSimTSAssociatorByEnergyScoreImpl::TSToSimTSAssociatorByEnergyScoreImpl(
0006 edm::EDProductGetter const& productGetter,
0007 bool hardScatterOnly,
0008 std::shared_ptr<hgcal::RecHitTools> recHitTools,
0009 const std::unordered_map<DetId, const HGCRecHit*>* hitMap)
0010 : hardScatterOnly_(hardScatterOnly), recHitTools_(recHitTools), hitMap_(hitMap), productGetter_(&productGetter) {
0011 layers_ = recHitTools_->lastLayerBH();
0012 }
0013
0014 hgcal::association TSToSimTSAssociatorByEnergyScoreImpl::makeConnections(
0015 const edm::Handle<ticl::TracksterCollection>& tCH,
0016 const edm::Handle<reco::CaloClusterCollection>& lCCH,
0017 const edm::Handle<ticl::TracksterCollection>& sTCH) const {
0018
0019 const auto& tracksters = *tCH.product();
0020 const auto& layerClusters = *lCCH.product();
0021 const auto& simTracksters = *sTCH.product();
0022 auto nTracksters = tracksters.size();
0023
0024
0025 auto nSimTracksters = simTracksters.size();
0026 std::vector<size_t> sTIndices;
0027 for (unsigned int stId = 0; stId < nSimTracksters; ++stId) {
0028 if (simTracksters[stId].vertices().empty()) {
0029 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl")
0030 << "Excluding SimTrackster " << stId << " witH no vertices!" << std::endl;
0031 continue;
0032 }
0033 sTIndices.emplace_back(stId);
0034 }
0035 nSimTracksters = sTIndices.size();
0036
0037
0038
0039
0040
0041 hgcal::simTracksterToTrackster tssInSimTrackster;
0042 tssInSimTrackster.resize(nSimTracksters);
0043 for (unsigned int i = 0; i < nSimTracksters; ++i) {
0044 tssInSimTrackster[i].simTracksterId = i;
0045 tssInSimTrackster[i].energy = 0.f;
0046 tssInSimTrackster[i].lcs_and_fractions.clear();
0047 }
0048
0049
0050
0051
0052
0053 std::unordered_map<int, std::vector<hgcal::lcInfoInTrackster>> lcToSimTracksterId_Map;
0054 for (const auto& stId : sTIndices) {
0055 const auto& lcs = simTracksters[stId].vertices();
0056 int lcInSimTst = 0;
0057 for (const auto& lcId : lcs) {
0058 const auto fraction = 1.f / simTracksters[stId].vertex_multiplicity(lcInSimTst++);
0059
0060 const auto lc_find_it = lcToSimTracksterId_Map.find(lcId);
0061 if (lc_find_it == lcToSimTracksterId_Map.end()) {
0062 lcToSimTracksterId_Map[lcId] = std::vector<hgcal::lcInfoInTrackster>();
0063 }
0064 lcToSimTracksterId_Map[lcId].emplace_back(stId, fraction);
0065
0066 tssInSimTrackster[stId].energy += fraction * layerClusters[lcId].energy();
0067 tssInSimTrackster[stId].lcs_and_fractions.emplace_back(lcId, fraction);
0068 }
0069 }
0070
0071 #ifdef EDM_ML_DEBUG
0072 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl")
0073 << "tssInSimTrackster INFO (Only SimTrackster filled at the moment)" << std::endl;
0074 for (size_t st = 0; st < tssInSimTrackster.size(); ++st) {
0075 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl") << "For SimTrackster Idx: " << st << " we have: " << std::endl;
0076 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl")
0077 << "\tSimTracksterIdx:\t" << tssInSimTrackster[st].simTracksterId << std::endl;
0078 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl") << "\tEnergy:\t" << tssInSimTrackster[st].energy << std::endl;
0079 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl") << "\t# of clusters:\t" << layerClusters.size() << std::endl;
0080 double tot_energy = 0.;
0081 for (auto const& lcaf : tssInSimTrackster[st].lcs_and_fractions) {
0082 const auto lcEn = lcaf.second * layerClusters[lcaf.first].energy();
0083 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl")
0084 << "\tLC/fraction/energy: " << (uint32_t)lcaf.first << "/" << lcaf.second << "/" << lcEn << std::endl;
0085 tot_energy += lcEn;
0086 }
0087 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl") << "\tTot Sum lcaf: " << tot_energy << std::endl;
0088 for (auto const& ts : tssInSimTrackster[st].tracksterIdToEnergyAndScore) {
0089 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl")
0090 << "\ttsIdx/energy/score: " << ts.first << "/" << ts.second.first << "/" << ts.second.second << std::endl;
0091 }
0092 }
0093
0094 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl") << "lcToSimTracksterId_Map INFO" << std::endl;
0095 for (auto const& lc : lcToSimTracksterId_Map) {
0096 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl")
0097 << "For lcId: " << (uint32_t)lc.first
0098 << " we have found the following connections with SimTracksters:" << std::endl;
0099 for (auto const& st : lc.second) {
0100 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl")
0101 << "\tSimTrackster Id: " << st.clusterId << " with fraction: " << st.fraction
0102 << " and energy: " << st.fraction * layerClusters[lc.first].energy() << std::endl;
0103 }
0104 }
0105 #endif
0106
0107
0108
0109
0110
0111 hgcal::tracksterToSimTrackster stsInTrackster;
0112 stsInTrackster.resize(nTracksters);
0113
0114 for (unsigned int tsId = 0; tsId < nTracksters; ++tsId) {
0115 for (unsigned int i = 0; i < tracksters[tsId].vertices().size(); ++i) {
0116 const auto lcId = tracksters[tsId].vertices(i);
0117 const auto lcFractionInTs = 1.f / tracksters[tsId].vertex_multiplicity(i);
0118
0119 const auto lc_find_in_ST = lcToSimTracksterId_Map.find(lcId);
0120
0121 if (lc_find_in_ST != lcToSimTracksterId_Map.end()) {
0122
0123
0124
0125 for (const auto& st : lc_find_in_ST->second) {
0126
0127
0128 tssInSimTrackster[st.clusterId].tracksterIdToEnergyAndScore[tsId].first +=
0129 lcFractionInTs * st.fraction * layerClusters[lcId].energy();
0130
0131 stsInTrackster[tsId].emplace_back(st.clusterId, 0.f);
0132 }
0133 }
0134 }
0135 }
0136
0137 #ifdef EDM_ML_DEBUG
0138 for (unsigned int tsId = 0; tsId < nTracksters; ++tsId) {
0139 for (const auto& lcId : tracksters[tsId].vertices()) {
0140 const auto& hits_and_fractions = layerClusters[lcId].hitsAndFractions();
0141 unsigned int numberOfHitsInLC = hits_and_fractions.size();
0142
0143
0144
0145
0146
0147
0148
0149
0150
0151 std::vector<int> hitsToSimTracksterId(numberOfHitsInLC);
0152
0153 int maxSTId_byNumberOfHits = -1;
0154
0155 unsigned int maxSTNumberOfHitsInLC = 0;
0156
0157 int maxSTId_byEnergy = -1;
0158
0159 float maxEnergySharedLCandST = 0.f;
0160
0161 float energyFractionOfLCinST = 0.f;
0162
0163 float energyFractionOfSTinLC = 0.f;
0164 std::unordered_map<unsigned, unsigned> occurrencesSTinLC;
0165 unsigned int numberOfNoiseHitsInLC = 0;
0166 std::unordered_map<unsigned, float> STEnergyInLC;
0167
0168 const auto lc_find_in_ST = lcToSimTracksterId_Map.find(lcId);
0169 for (unsigned int hitId = 0; hitId < numberOfHitsInLC; hitId++) {
0170 const auto rhFraction = hits_and_fractions[hitId].second;
0171
0172
0173
0174
0175
0176
0177 if (rhFraction == 0.) {
0178 hitsToSimTracksterId[hitId] = -2;
0179 }
0180
0181 if (lc_find_in_ST == lcToSimTracksterId_Map.end()) {
0182 hitsToSimTracksterId[hitId] -= 1;
0183 } else {
0184 auto maxSTEnergyInLC = 0.f;
0185 auto maxSTId = -1;
0186
0187 for (const auto& st : lc_find_in_ST->second) {
0188 const auto stId = st.clusterId;
0189 STEnergyInLC[stId] += st.fraction * layerClusters[lcId].energy();
0190
0191
0192 if (STEnergyInLC[stId] > maxSTEnergyInLC) {
0193 maxSTEnergyInLC = STEnergyInLC[stId];
0194 maxSTId = stId;
0195 }
0196 }
0197 hitsToSimTracksterId[hitId] = maxSTId;
0198 }
0199 }
0200
0201 for (const auto& c : hitsToSimTracksterId) {
0202 if (c < 0) {
0203 numberOfNoiseHitsInLC++;
0204 } else {
0205 occurrencesSTinLC[c]++;
0206 }
0207 }
0208
0209 for (const auto& c : occurrencesSTinLC) {
0210 if (c.second > maxSTNumberOfHitsInLC) {
0211 maxSTId_byNumberOfHits = c.first;
0212 maxSTNumberOfHitsInLC = c.second;
0213 }
0214 }
0215
0216 for (const auto& c : STEnergyInLC) {
0217 if (c.second > maxEnergySharedLCandST) {
0218 maxSTId_byEnergy = c.first;
0219 maxEnergySharedLCandST = c.second;
0220 }
0221 }
0222
0223 float totalSTEnergyOnLayer = 0.f;
0224 if (maxSTId_byEnergy >= 0) {
0225 totalSTEnergyOnLayer = tssInSimTrackster[maxSTId_byEnergy].energy;
0226 energyFractionOfSTinLC = maxEnergySharedLCandST / totalSTEnergyOnLayer;
0227 if (tracksters[tsId].raw_energy() > 0.f) {
0228 energyFractionOfLCinST = maxEnergySharedLCandST / tracksters[tsId].raw_energy();
0229 }
0230 }
0231
0232 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl")
0233 << std::setw(12) << "TracksterID:\t" << std::setw(12) << "layerCluster\t" << std::setw(10) << "lc energy\t"
0234 << std::setw(5) << "nhits\t" << std::setw(12) << "noise hits\t" << std::setw(22) << "maxSTId_byNumberOfHits\t"
0235 << std::setw(8) << "nhitsST\t" << std::setw(13) << "maxSTId_byEnergy\t" << std::setw(20)
0236 << "maxEnergySharedLCandST\t" << std::setw(22) << "totalSTEnergyOnLayer\t" << std::setw(22)
0237 << "energyFractionOfLCinST\t" << std::setw(25) << "energyFractionOfSTinLC\t"
0238 << "\n";
0239 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl")
0240 << std::setw(12) << tsId << "\t" << std::setw(12) << lcId << "\t" << std::setw(10)
0241 << tracksters[tsId].raw_energy() << "\t" << std::setw(5) << numberOfHitsInLC << "\t" << std::setw(12)
0242 << numberOfNoiseHitsInLC << "\t" << std::setw(22) << maxSTId_byNumberOfHits << "\t" << std::setw(8)
0243 << maxSTNumberOfHitsInLC << "\t" << std::setw(13) << maxSTId_byEnergy << "\t" << std::setw(20)
0244 << maxEnergySharedLCandST << "\t" << std::setw(22) << totalSTEnergyOnLayer << "\t" << std::setw(22)
0245 << energyFractionOfLCinST << "\t" << std::setw(25) << energyFractionOfSTinLC << "\n";
0246 }
0247 }
0248
0249 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl")
0250 << "Improved tssInSimTrackster INFO (Now containing the linked tracksters id and energy - score still empty)"
0251 << std::endl;
0252 for (size_t st = 0; st < tssInSimTrackster.size(); ++st) {
0253 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl") << "For SimTrackster Idx: " << st << " we have: " << std::endl;
0254 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl")
0255 << " SimTracksterIdx: " << tssInSimTrackster[st].simTracksterId << std::endl;
0256 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl") << "\tEnergy:\t" << tssInSimTrackster[st].energy << std::endl;
0257 double tot_energy = 0.;
0258 for (auto const& lcaf : tssInSimTrackster[st].lcs_and_fractions) {
0259 const auto lcEn = lcaf.second * layerClusters[lcaf.first].energy();
0260 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl")
0261 << "\tLC/fraction/energy: " << (uint32_t)lcaf.first << "/" << lcaf.second << "/" << lcEn << std::endl;
0262 tot_energy += lcEn;
0263 }
0264 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl") << "\tTot Sum lcaf: " << tot_energy << std::endl;
0265 for (auto const& ts : tssInSimTrackster[st].tracksterIdToEnergyAndScore) {
0266 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl")
0267 << "\ttsIdx/energy/score: " << ts.first << "/" << ts.second.first << "/" << ts.second.second << std::endl;
0268 }
0269 }
0270
0271 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl") << "Improved lcToSimTracksterId_Map INFO" << std::endl;
0272 for (auto const& lc : lcToSimTracksterId_Map) {
0273 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl")
0274 << "For lcId: " << (uint32_t)lc.first
0275 << " we have found the following connections with SimTracksters:" << std::endl;
0276 for (auto const& st : lc.second) {
0277 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl")
0278 << " SimTrackster Id: " << st.clusterId << " with fraction: " << st.fraction
0279 << " and energy: " << st.fraction * layerClusters[lc.first].energy() << std::endl;
0280 }
0281 }
0282 #endif
0283
0284
0285
0286 for (unsigned int tsId = 0; tsId < nTracksters; ++tsId) {
0287
0288
0289 std::sort(stsInTrackster[tsId].begin(), stsInTrackster[tsId].end());
0290 auto last = std::unique(stsInTrackster[tsId].begin(), stsInTrackster[tsId].end());
0291 stsInTrackster[tsId].erase(last, stsInTrackster[tsId].end());
0292
0293
0294
0295 if (tracksters[tsId].raw_energy() == 0. && !stsInTrackster[tsId].empty()) {
0296 for (auto& stPair : stsInTrackster[tsId]) {
0297 stPair.second = 1.;
0298 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl")
0299 << "TracksterId:\t " << tsId << "\tST id:\t" << stPair.first << "\tscore\t " << stPair.second << "\n";
0300 }
0301 continue;
0302 }
0303
0304 float invTracksterEnergyWeight = 0.f;
0305 for (unsigned int i = 0; i < tracksters[tsId].vertices().size(); ++i) {
0306 const auto lcId = tracksters[tsId].vertices(i);
0307 const auto lcFractionInTs = 1.f / tracksters[tsId].vertex_multiplicity(i);
0308
0309 const auto& hits_and_fractions = layerClusters[lcId].hitsAndFractions();
0310
0311 for (auto const& haf : hits_and_fractions) {
0312 invTracksterEnergyWeight += std::pow(lcFractionInTs * haf.second * hitMap_->at(haf.first)->energy(), 2);
0313 }
0314 }
0315 invTracksterEnergyWeight = 1.f / invTracksterEnergyWeight;
0316
0317 for (unsigned int i = 0; i < tracksters[tsId].vertices().size(); ++i) {
0318 const auto lcId = tracksters[tsId].vertices(i);
0319 const auto lcFractionInTs = 1.f / tracksters[tsId].vertex_multiplicity(i);
0320
0321 const bool lcWithST = (lcToSimTracksterId_Map.find(lcId) != lcToSimTracksterId_Map.end());
0322
0323 float lcEnergyWeight = pow(layerClusters[lcId].energy(), 2);
0324
0325 for (auto& stPair : stsInTrackster[tsId]) {
0326 float stFraction = 0.f;
0327 if (lcWithST) {
0328 const auto findLCIt = std::find(lcToSimTracksterId_Map[lcId].begin(),
0329 lcToSimTracksterId_Map[lcId].end(),
0330 hgcal::lcInfoInTrackster{stPair.first, 0.f});
0331 if (findLCIt != lcToSimTracksterId_Map[lcId].end()) {
0332 stFraction = findLCIt->fraction;
0333 }
0334 }
0335 stPair.second +=
0336 (lcFractionInTs - stFraction) * (lcFractionInTs - stFraction) * lcEnergyWeight * invTracksterEnergyWeight;
0337 #ifdef EDM_ML_DEBUG
0338 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl")
0339 << "lcId:\t" << (uint32_t)lcId << "\ttracksterId:\t" << tsId << "\ttsFraction,stFraction:\t"
0340 << lcFractionInTs << ", " << stFraction << "\tlcEnergyWeight:\t" << lcEnergyWeight << "\tcurrent score:\t"
0341 << stPair.second << "\tinvTracksterEnergyWeight:\t" << invTracksterEnergyWeight << "\n";
0342 #endif
0343 }
0344 }
0345
0346 #ifdef EDM_ML_DEBUG
0347 if (stsInTrackster[tsId].empty())
0348 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl") << "trackster Id:\t" << tsId << "\tST id:\t-1"
0349 << "\tscore\t-1\n";
0350 #endif
0351 }
0352
0353
0354 for (const auto& stId : sTIndices) {
0355 float invSTEnergyWeight = 0.f;
0356
0357 const unsigned int STNumberOfLCs = tssInSimTrackster[stId].lcs_and_fractions.size();
0358 if (STNumberOfLCs == 0)
0359 continue;
0360
0361 #ifdef EDM_ML_DEBUG
0362 int tsWithMaxEnergyInST = -1;
0363
0364 float maxEnergyTSinST = 0.f;
0365 float STenergy = tssInSimTrackster[stId].energy;
0366
0367 float STEnergyFractionInTS = 0.f;
0368 for (const auto& ts : tssInSimTrackster[stId].tracksterIdToEnergyAndScore) {
0369 if (ts.second.first > maxEnergyTSinST) {
0370 maxEnergyTSinST = ts.second.first;
0371 tsWithMaxEnergyInST = ts.first;
0372 }
0373 }
0374 if (STenergy > 0.f)
0375 STEnergyFractionInTS = maxEnergyTSinST / STenergy;
0376
0377 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl")
0378 << std::setw(12) << "simTrackster\t" << std::setw(15) << "st total energy\t" << std::setw(15)
0379 << "stEnergyOnLayer\t" << std::setw(14) << "STNhitsOnLayer\t" << std::setw(18) << "tsWithMaxEnergyInST\t"
0380 << std::setw(15) << "maxEnergyTSinST\t" << std::setw(20) << "STEnergyFractionInTS"
0381 << "\n";
0382 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl")
0383 << std::setw(12) << stId << "\t" << std::setw(15) << simTracksters[stId].raw_energy() << "\t" << std::setw(15)
0384 << STenergy << "\t" << std::setw(14) << STNumberOfLCs << "\t" << std::setw(18) << tsWithMaxEnergyInST << "\t"
0385 << std::setw(15) << maxEnergyTSinST << "\t" << std::setw(20) << STEnergyFractionInTS << "\n";
0386 #endif
0387
0388
0389 for (auto const& lcaf : tssInSimTrackster[stId].lcs_and_fractions) {
0390 invSTEnergyWeight += std::pow(lcaf.second * layerClusters[lcaf.first].energy(), 2);
0391 }
0392 invSTEnergyWeight = 1.f / invSTEnergyWeight;
0393
0394 for (unsigned int i = 0; i < STNumberOfLCs; ++i) {
0395 auto& st_lcId = tssInSimTrackster[stId].lcs_and_fractions[i].first;
0396 auto& st_lcFraction = tssInSimTrackster[stId].lcs_and_fractions[i].second;
0397
0398 if (st_lcFraction == 0.f)
0399 continue;
0400 float lcEnergyWeight = pow(layerClusters[st_lcId].energy(), 2);
0401 for (auto& tsPair : tssInSimTrackster[stId].tracksterIdToEnergyAndScore) {
0402 unsigned int tsId = tsPair.first;
0403
0404 for (unsigned int i = 0; i < tracksters[tsId].vertices().size(); ++i) {
0405 const auto tsFraction = 1.f / tracksters[tsId].vertex_multiplicity(i);
0406
0407 tsPair.second.second +=
0408 (tsFraction - st_lcFraction) * (tsFraction - st_lcFraction) * lcEnergyWeight * invSTEnergyWeight;
0409 #ifdef EDM_ML_DEBUG
0410 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl")
0411 << "STLCId:\t" << (uint32_t)st_lcId << "\tTracksterId:\t" << tsId << "\t"
0412 << "tsFraction, stFraction:\t" << tsFraction << ", " << st_lcFraction << "\t"
0413 << "lcEnergyWeight:\t" << lcEnergyWeight << "\t"
0414 << "current score:\t" << tsPair.second.second << "\t"
0415 << "invSTEnergyWeight:\t" << invSTEnergyWeight << "\n";
0416 #endif
0417 }
0418 }
0419 }
0420 #ifdef EDM_ML_DEBUG
0421 if (tssInSimTrackster[stId].tracksterIdToEnergyAndScore.empty())
0422 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl") << "ST Id:\t" << stId << "\tTS id:\t-1 "
0423 << "\tscore\t-1\n";
0424
0425 for (const auto& tsPair : tssInSimTrackster[stId].tracksterIdToEnergyAndScore) {
0426 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl")
0427 << "ST Id: \t" << stId << "\t TS id: \t" << tsPair.first << "\t score \t" << tsPair.second.second
0428 << "\t shared energy:\t" << tsPair.second.first << "\t shared energy fraction:\t"
0429 << (tsPair.second.first / STenergy) << "\n";
0430 }
0431 #endif
0432 }
0433 return {stsInTrackster, tssInSimTrackster};
0434 }
0435
0436 hgcal::RecoToSimCollectionSimTracksters TSToSimTSAssociatorByEnergyScoreImpl::associateRecoToSim(
0437 const edm::Handle<ticl::TracksterCollection>& tCH,
0438 const edm::Handle<reco::CaloClusterCollection>& lCCH,
0439 const edm::Handle<ticl::TracksterCollection>& sTCH) const {
0440 hgcal::RecoToSimCollectionSimTracksters returnValue(productGetter_);
0441 const auto& links = makeConnections(tCH, lCCH, sTCH);
0442
0443 const auto& stsInTrackster = std::get<0>(links);
0444 for (size_t tsId = 0; tsId < stsInTrackster.size(); ++tsId) {
0445 for (auto& stPair : stsInTrackster[tsId]) {
0446 LogDebug("TSToSimTSAssociatorByEnergyScoreImpl") << "Trackster Id:\t" << tsId << "\tSimTrackster id:\t"
0447 << stPair.first << "\tscore:\t" << stPair.second << "\n";
0448
0449 returnValue.insert(edm::Ref<ticl::TracksterCollection>(tCH, tsId),
0450 std::make_pair(edm::Ref<ticl::TracksterCollection>(sTCH, stPair.first),
0451 stPair.second)
0452 );
0453 }
0454 }
0455 return returnValue;
0456 }
0457
0458 hgcal::SimToRecoCollectionSimTracksters TSToSimTSAssociatorByEnergyScoreImpl::associateSimToReco(
0459 const edm::Handle<ticl::TracksterCollection>& tCH,
0460 const edm::Handle<reco::CaloClusterCollection>& lCCH,
0461 const edm::Handle<ticl::TracksterCollection>& sTCH) const {
0462 hgcal::SimToRecoCollectionSimTracksters returnValue(productGetter_);
0463 const auto& links = makeConnections(tCH, lCCH, sTCH);
0464 const auto& tssInSimTrackster = std::get<1>(links);
0465 for (size_t stId = 0; stId < tssInSimTrackster.size(); ++stId) {
0466 for (auto& tsPair : tssInSimTrackster[stId].tracksterIdToEnergyAndScore) {
0467 returnValue.insert(
0468 edm::Ref<ticl::TracksterCollection>(sTCH, stId),
0469 std::make_pair(edm::Ref<ticl::TracksterCollection>(tCH, tsPair.first),
0470 std::make_pair(tsPair.second.first, tsPair.second.second))
0471 );
0472 }
0473 }
0474 return returnValue;
0475 }