Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 14:24:32

0001 #include "RecoLocalCalo/HGCalRecAlgos/interface/HGCalDepthPreClusterer.h"
0002 #include "DataFormats/Math/interface/deltaR.h"
0003 
0004 #include <list>
0005 
0006 namespace {
0007   std::vector<size_t> sorted_indices(const reco::HGCalMultiCluster::ClusterCollection &v) {
0008     // initialize original index locations
0009     std::vector<size_t> idx(v.size());
0010     for (size_t i = 0; i != idx.size(); ++i)
0011       idx[i] = i;
0012 
0013     // sort indices based on comparing values in v
0014     std::sort(idx.begin(), idx.end(), [&v](size_t i1, size_t i2) { return (*v[i1]) > (*v[i2]); });
0015 
0016     return idx;
0017   }
0018 
0019   float dist2(const edm::Ptr<reco::BasicCluster> &a, const edm::Ptr<reco::BasicCluster> &b) {
0020     return reco::deltaR2(*a, *b);
0021   }
0022 
0023   //get distance between cluster and multicluster axis (defined by remaning cluster with highest energy)
0024   // N.B. the order of the clusters matters
0025   float distAxisCluster2(const edm::Ptr<reco::BasicCluster> &a, const edm::Ptr<reco::BasicCluster> &b) {
0026     float tanTheta = tan(2 * atan(exp(-1 * a->eta())));
0027     float ax = b->z() * tanTheta * cos(a->phi());
0028     float ay = b->z() * tanTheta * sin(a->phi());
0029     return (ax - b->x()) * (ax - b->x()) + (ay - b->y()) * (ay - b->y());
0030   }
0031 }  // namespace
0032 
0033 std::vector<reco::HGCalMultiCluster> HGCalDepthPreClusterer::makePreClusters(
0034     const reco::HGCalMultiCluster::ClusterCollection &thecls) const {
0035   std::vector<reco::HGCalMultiCluster> thePreClusters;
0036   std::vector<size_t> es = sorted_indices(thecls);
0037   std::vector<int> vused(es.size(), 0);
0038   unsigned int used = 0;
0039 
0040   for (unsigned int i = 0; i < es.size(); ++i) {
0041     if (vused[i] == 0) {
0042       reco::HGCalMultiCluster temp;
0043       temp.push_back(thecls[es[i]]);
0044       vused[i] = (thecls[es[i]]->z() > 0) ? 1 : -1;
0045       ++used;
0046       for (unsigned int j = i + 1; j < es.size(); ++j) {
0047         if (vused[j] == 0) {
0048           float distanceCheck = 9999.;
0049           if (realSpaceCone)
0050             distanceCheck = distAxisCluster2(thecls[es[i]], thecls[es[j]]);
0051           else
0052             distanceCheck = dist2(thecls[es[i]], thecls[es[j]]);
0053           DetId detid = thecls[es[j]]->hitsAndFractions()[0].first();
0054           unsigned int layer = clusterTools->getLayer(detid);
0055           float radius = radii[2];
0056           if (layer <= rhtools_.lastLayerEE())
0057             radius = radii[0];
0058           else if (layer < rhtools_.firstLayerBH())
0059             radius = radii[1];
0060           float radius2 = radius * radius;
0061           if (distanceCheck<radius2 &&int(thecls[es[j]]->z() * vused[i])> 0) {
0062             temp.push_back(thecls[es[j]]);
0063             vused[j] = vused[i];
0064             ++used;
0065           }
0066         }
0067       }
0068       if (temp.size() > minClusters) {
0069         thePreClusters.push_back(temp);
0070         auto &back = thePreClusters.back();
0071         back.setPosition(clusterTools->getMultiClusterPosition(back));
0072         back.setEnergy(clusterTools->getMultiClusterEnergy(back));
0073       }
0074     }
0075   }
0076 
0077   return thePreClusters;
0078 }