EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
SeedFilter.ipp
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file SeedFilter.ipp
1 // This file is part of the Acts project.
2 //
3 // Copyright (C) 2018 CERN for the benefit of the Acts project
4 //
5 // This Source Code Form is subject to the terms of the Mozilla Public
6 // License, v. 2.0. If a copy of the MPL was not distributed with this
7 // file, You can obtain one at http://mozilla.org/MPL/2.0/.
8 
9 #include <utility>
10 
11 namespace Acts {
12 // constructor
13 template <typename external_spacepoint_t>
17  : m_cfg(config), m_experimentCuts(expCuts) {}
18 
19 // function to filter seeds based on all seeds with same bottom- and
20 // middle-spacepoint.
21 // return vector must contain weight of each seed
22 template <typename external_spacepoint_t>
23 std::vector<std::pair<
24  float, std::unique_ptr<const InternalSeed<external_spacepoint_t>>>>
29  std::vector<float>& invHelixDiameterVec,
30  std::vector<float>& impactParametersVec, float zOrigin) const {
31  std::vector<std::pair<
32  float, std::unique_ptr<const InternalSeed<external_spacepoint_t>>>>
33  selectedSeeds;
34 
35  for (size_t i = 0; i < topSpVec.size(); i++) {
36  // if two compatible seeds with high distance in r are found, compatible
37  // seeds span 5 layers
38  // -> very good seed
39  std::vector<float> compatibleSeedR;
40 
41  float invHelixDiameter = invHelixDiameterVec[i];
42  float lowerLimitCurv = invHelixDiameter - m_cfg.deltaInvHelixDiameter;
43  float upperLimitCurv = invHelixDiameter + m_cfg.deltaInvHelixDiameter;
44  float currentTop_r = topSpVec[i]->radius();
45  float impact = impactParametersVec[i];
46 
47  float weight = -(impact * m_cfg.impactWeightFactor);
48  for (size_t j = 0; j < topSpVec.size(); j++) {
49  if (i == j) {
50  continue;
51  }
52  // compared top SP should have at least deltaRMin distance
53  float otherTop_r = topSpVec[j]->radius();
54  float deltaR = currentTop_r - otherTop_r;
55  if (std::abs(deltaR) < m_cfg.deltaRMin) {
56  continue;
57  }
58  // curvature difference within limits?
59  // TODO: how much slower than sorting all vectors by curvature
60  // and breaking out of loop? i.e. is vector size large (e.g. in jets?)
61  if (invHelixDiameterVec[j] < lowerLimitCurv) {
62  continue;
63  }
64  if (invHelixDiameterVec[j] > upperLimitCurv) {
65  continue;
66  }
67  bool newCompSeed = true;
68  for (float previousDiameter : compatibleSeedR) {
69  // original ATLAS code uses higher min distance for 2nd found compatible
70  // seed (20mm instead of 5mm)
71  // add new compatible seed only if distance larger than rmin to all
72  // other compatible seeds
73  if (std::abs(previousDiameter - otherTop_r) < m_cfg.deltaRMin) {
74  newCompSeed = false;
75  break;
76  }
77  }
78  if (newCompSeed) {
79  compatibleSeedR.push_back(otherTop_r);
80  weight += m_cfg.compatSeedWeight;
81  }
82  if (compatibleSeedR.size() >= m_cfg.compatSeedLimit) {
83  break;
84  }
85  }
86 
87  if (m_experimentCuts != nullptr) {
88  // add detector specific considerations on the seed weight
89  weight += m_experimentCuts->seedWeight(bottomSP, middleSP, *topSpVec[i]);
90  // discard seeds according to detector specific cuts (e.g.: weight)
91  if (!m_experimentCuts->singleSeedCut(weight, bottomSP, middleSP,
92  *topSpVec[i])) {
93  continue;
94  }
95  }
96  selectedSeeds.push_back(std::make_pair(
97  weight, std::make_unique<const InternalSeed<external_spacepoint_t>>(
98  bottomSP, middleSP, *topSpVec[i], zOrigin)));
99  }
100  return selectedSeeds;
101 }
102 
103 // after creating all seeds with a common middle space point, filter again
104 template <typename external_spacepoint_t>
106  std::vector<std::pair<
107  float, std::unique_ptr<const InternalSeed<external_spacepoint_t>>>>&
108  seedsPerSpM,
109  std::vector<Seed<external_spacepoint_t>>& outVec) const {
110  // sort by weight and iterate only up to configured max number of seeds per
111  // middle SP
112  std::sort((seedsPerSpM.begin()), (seedsPerSpM.end()),
113  [](const std::pair<float, std::unique_ptr<const Acts::InternalSeed<
114  external_spacepoint_t>>>& i1,
115  const std::pair<float, std::unique_ptr<const Acts::InternalSeed<
116  external_spacepoint_t>>>& i2) {
117  if (i1.first != i2.first) {
118  return i1.first > i2.first;
119  } else {
120  // This is for the case when the weights from different seeds
121  // are same. This makes cpu & cuda results same
122  float seed1_sum = 0;
123  float seed2_sum = 0;
124  for (int i = 0; i < 3; i++) {
125  seed1_sum += pow(i1.second->sp[i]->sp().y(), 2) +
126  pow(i1.second->sp[i]->sp().z(), 2);
127  seed2_sum += pow(i2.second->sp[i]->sp().y(), 2) +
128  pow(i2.second->sp[i]->sp().z(), 2);
129  }
130  return seed1_sum > seed2_sum;
131  }
132  });
133  if (m_experimentCuts != nullptr) {
134  seedsPerSpM = m_experimentCuts->cutPerMiddleSP(std::move(seedsPerSpM));
135  }
136  unsigned int maxSeeds = seedsPerSpM.size();
137 
138  if (maxSeeds > m_cfg.maxSeedsPerSpM) {
139  maxSeeds = m_cfg.maxSeedsPerSpM + 1;
140  }
141  auto itBegin = seedsPerSpM.begin();
142  auto it = seedsPerSpM.begin();
143  // default filter removes the last seeds if maximum amount exceeded
144  // ordering by weight by filterSeeds_2SpFixed means these are the lowest
145  // weight seeds
146  for (; it < itBegin + maxSeeds; ++it) {
147  outVec.push_back(Seed<external_spacepoint_t>(
148  (*it).second->sp[0]->sp(), (*it).second->sp[1]->sp(),
149  (*it).second->sp[2]->sp(), (*it).second->z()));
150  }
151 }
152 
153 } // namespace Acts