Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-02-18 08:23:38

0001 #include "L1Trigger/L1THGCal/interface/backend/HGCalHistoClusteringImpl_SA.h"
0002 
0003 #include <map>
0004 #include <cmath>
0005 
0006 std::vector<l1thgcfirmware::HGCalMulticluster> HGCalHistoClusteringImplSA::clusterSeedMulticluster_SA(
0007     const std::vector<l1thgcfirmware::HGCalCluster>& clusters,
0008     const std::vector<l1thgcfirmware::HGCalSeed>& seeds,
0009     std::vector<l1thgcfirmware::HGCalCluster>& rejected_clusters,
0010     const l1thgcfirmware::ClusterAlgoConfig& configuration) const {
0011   std::map<int, l1thgcfirmware::HGCalMulticluster> mapSeedMulticluster;
0012   std::vector<l1thgcfirmware::HGCalMulticluster> multiclustersOut;
0013 
0014   for (const auto& clu : clusters) {
0015     int z_side = clu.zside();
0016 
0017     double radiusCoefficientA = configuration.dr_byLayer_coefficientA().empty()
0018                                     ? configuration.dr()
0019                                     : configuration.dr_byLayer_coefficientA()[clu.layer()];
0020     double radiusCoefficientB =
0021         configuration.dr_byLayer_coefficientB().empty() ? 0 : configuration.dr_byLayer_coefficientB()[clu.layer()];
0022 
0023     double minDistSqrd = radiusCoefficientA + radiusCoefficientB * (configuration.midRadius() - std::abs(clu.eta()));
0024     minDistSqrd *= minDistSqrd;
0025 
0026     std::vector<std::pair<int, double>> targetSeedsEnergy;
0027 
0028     unsigned int iseed = 0;
0029     for (const auto& seed : seeds) {
0030       if (z_side * seed.z() < 0) {
0031         ++iseed;
0032         continue;
0033       }
0034 
0035       double seedEnergy = seed.energy();
0036 
0037       double d = (clu.x() - seed.x()) * (clu.x() - seed.x()) + (clu.y() - seed.y()) * (clu.y() - seed.y());
0038 
0039       if (d < minDistSqrd) {
0040         // NearestNeighbour
0041         minDistSqrd = d;
0042 
0043         if (targetSeedsEnergy.empty()) {
0044           targetSeedsEnergy.emplace_back(iseed, seedEnergy);
0045         } else {
0046           targetSeedsEnergy.at(0).first = iseed;
0047           targetSeedsEnergy.at(0).second = seedEnergy;
0048         }
0049       }
0050       ++iseed;
0051     }
0052 
0053     if (targetSeedsEnergy.empty()) {
0054       rejected_clusters.emplace_back(clu);
0055       continue;
0056     }
0057 
0058     // N.B. as I have only implemented NearestNeighbour option
0059     // then targetSeedsEnergy has at most 1 seed for this cluster
0060     // Leaving in some redundant functionality in case we need
0061     // EnergySplit option
0062 
0063     //Loop over target seeds and divide up the clusters energy
0064     double totalTargetSeedEnergy = 0;
0065     for (const auto& energy : targetSeedsEnergy) {
0066       totalTargetSeedEnergy += energy.second;
0067     }
0068 
0069     for (const auto& energy : targetSeedsEnergy) {
0070       double seedWeight = 1;
0071       if (mapSeedMulticluster[energy.first].size() == 0) {
0072         mapSeedMulticluster[energy.first] = l1thgcfirmware::HGCalMulticluster(clu, 1);
0073       } else {
0074         mapSeedMulticluster[energy.first].addConstituent(clu, true, seedWeight);
0075       }
0076     }
0077   }
0078 
0079   multiclustersOut.reserve(mapSeedMulticluster.size());
0080   for (const auto& mclu : mapSeedMulticluster)
0081     multiclustersOut.emplace_back(mclu.second);
0082 
0083   return multiclustersOut;
0084 }
0085 
0086 void HGCalHistoClusteringImplSA::finalizeClusters_SA(
0087     const std::vector<l1thgcfirmware::HGCalMulticluster>& multiclusters_in,
0088     const std::vector<l1thgcfirmware::HGCalCluster>& rejected_clusters_in,
0089     std::vector<l1thgcfirmware::HGCalMulticluster>& multiclusters_out,
0090     std::vector<l1thgcfirmware::HGCalCluster>& rejected_clusters_out,
0091     const l1thgcfirmware::ClusterAlgoConfig& configuration) const {
0092   for (const auto& tc : rejected_clusters_in) {
0093     rejected_clusters_out.push_back(tc);
0094   }
0095 
0096   for (const auto& multicluster : multiclusters_in) {
0097     if (multicluster.sumPt() > configuration.ptC3dThreshold()) {
0098       multiclusters_out.push_back(multicluster);
0099     } else {
0100       for (const auto& tc : multicluster.constituents()) {
0101         rejected_clusters_out.push_back(tc);
0102       }
0103     }
0104   }
0105 }