Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-09-07 22:45:14

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   // Get collections
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   //There shouldn't be any SimTrackster without vertices, but maybe they will be added later.
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   // Initialize tssInSimTrackster. It contains the simTracksterOnLayer structure for all CaloParticles in each layer and
0038   // among other the information to compute the SimTrackster-To-Trackster score. It is one of the two objects that
0039   // build the output of the makeConnections function.
0040   // tssInSimTrackster[stId]
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   // Fill lcToSimTracksterId_Map and update tssInSimTrackster
0050   // The lcToSimTracksterId_Map is used to connect a LayerCluster, via its id (key), with all the SimTracksters that
0051   // contributed to that LayerCluster by storing the SimTrackster id and the fraction of the LayerCluster's energy
0052   // in which the SimTrackster contributed.
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   }  // end of loop over SimTracksters
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   // Fill stsInTrackster and update tssInSimTrackster
0108   // this contains the ids of the simTracksters contributing with at least one
0109   // hit to the Trackster. To be returned since this contains the information
0110   // to compute the Trackster-To-SimTrackster score.
0111   hgcal::tracksterToSimTrackster stsInTrackster;  // tsId->(stId,score)
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         //Loops through all the simTracksters that contain the layer cluster under study
0123         //Here is time to update the tssInSimTrackster and connect the SimTrackster with all
0124         //the Tracksters that have the current layer cluster matched.
0125         for (const auto& st : lc_find_in_ST->second) {
0126           //tssInSimTrackster[simTracksterId][layerclusterId]-> (energy,score)
0127           //ST_i - > TS_j, TS_k, ...
0128           tssInSimTrackster[st.clusterId].tracksterIdToEnergyAndScore[tsId].first +=
0129               lcFractionInTs * st.fraction * layerClusters[lcId].energy();
0130           //TS_i -> ST_j, ST_k, ...
0131           stsInTrackster[tsId].emplace_back(st.clusterId, 0.f);
0132         }
0133       }
0134     }  // End loop over LayerClusters in Trackster
0135   }    // End of loop over Tracksters
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       // This vector will store, for each hit in the Layercluster, the index of
0144       // the SimTrackster that contributed the most, in terms of energy, to it.
0145       // Special values are:
0146       //
0147       // -2  --> the reconstruction fraction of the RecHit is 0 (used in the past to monitor Halo Hits)
0148       // -3  --> same as before with the added condition that no SimTrackster has been linked to this RecHit
0149       // -1  --> the reco fraction is >0, but no SimTrackster has been linked to it
0150       // >=0 --> index of the linked SimTrackster
0151       std::vector<int> hitsToSimTracksterId(numberOfHitsInLC);
0152       // This will store the index of the SimTrackster linked to the LayerCluster that has the largest number of hits in common.
0153       int maxSTId_byNumberOfHits = -1;
0154       // This will store the maximum number of shared hits between a LayerCluster and a SimTrackster
0155       unsigned int maxSTNumberOfHitsInLC = 0;
0156       // This will store the index of the SimTrackster linked to the LayerCluster that has the largest energy in common.
0157       int maxSTId_byEnergy = -1;
0158       // This will store the maximum number of shared energy between a LayerCluster and a SimTrackster
0159       float maxEnergySharedLCandST = 0.f;
0160       // This will store the fraction of the LayerCluster energy shared with the best(energy) SimTrackster: e_shared/lc_energy
0161       float energyFractionOfLCinST = 0.f;
0162       // This will store the fraction of the SimTrackster energy shared with the Trackster: e_shared/sc_energy
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         // if the fraction is zero or the hit does not belong to any SimTrackster,
0172         // set the SimTrackster Id for the hit to -1; this will
0173         // contribute to the number of noise hits
0174 
0175         // MR Remove the case in which the fraction is 0, since this could be a
0176         // real hit that has been marked as halo.
0177         if (rhFraction == 0.) {
0178           hitsToSimTracksterId[hitId] = -2;
0179         }
0180         //Now check if there are SimTracksters linked to this rechit of the layercluster
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           //Loop through all the linked SimTracksters
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             // Keep track of which SimTrackster contributed the most, in terms
0191             // of energy, to this specific Layer Cluster.
0192             if (STEnergyInLC[stId] > maxSTEnergyInLC) {
0193               maxSTEnergyInLC = STEnergyInLC[stId];
0194               maxSTId = stId;
0195             }
0196           }
0197           hitsToSimTracksterId[hitId] = maxSTId;
0198         }
0199       }  // End loop over hits on a LayerCluster
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     }  // End of loop over LayerClusters in Trackster
0247   }    // End of loop over Tracksters
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   // Update stsInTrackster; compute the score Trackster-to-SimTrackster,
0285   // together with the returned AssociationMap
0286   for (unsigned int tsId = 0; tsId < nTracksters; ++tsId) {
0287     // The SimTracksters contributing to the Trackster's LayerClusters should already be unique.
0288     // find the unique SimTracksters id contributing to the Trackster's LayerClusters
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     // If a reconstructed Trackster has energy 0 but is linked to a
0294     // SimTrackster, assigned score 1
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       // Compute the correct normalization
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     }  // End of loop over LayerClusters in Trackster
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   }  // End of loop over Tracksters
0352 
0353   // Compute the SimTrackster-To-Trackster score
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     //energy of the most energetic TS from all that were linked to ST
0364     float maxEnergyTSinST = 0.f;
0365     float STenergy = tssInSimTrackster[stId].energy;
0366     //most energetic TS from all TSs linked to ST over ST energy.
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     // Compute the correct normalization
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;  // hopefully this should never happen
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         }  // End of loop over Trackster's LayerClusters
0418       }    // End of loop over Tracksters linked to hits of this SimTrackster
0419     }      // End of loop over hits of SimTrackster on a Layer
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   }  // End loop over SimTrackster indices
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       // Fill AssociationMap
0449       returnValue.insert(edm::Ref<ticl::TracksterCollection>(tCH, tsId),  // Ref to TS
0450                          std::make_pair(edm::Ref<ticl::TracksterCollection>(sTCH, stPair.first),
0451                                         stPair.second)  // Pair <Ref to ST, score>
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),                           // Ref to ST
0469           std::make_pair(edm::Ref<ticl::TracksterCollection>(tCH, tsPair.first),     // Pair <Ref to TS,
0470                          std::make_pair(tsPair.second.first, tsPair.second.second))  // pair <energy, score> >
0471       );
0472     }
0473   }
0474   return returnValue;
0475 }