ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
HelixHough_findPairs.cpp
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file HelixHough_findPairs.cpp
1 #include "HelixHough.h"
2 #include "HelixRange.h"
3 
4 #include <cstddef>
5 #include <map>
6 #include <memory>
7 #include <set>
8 #include <sys/time.h>
9 #include <utility>
10 #include <vector>
11 
12 class SimpleTrack3D;
13 
14 using namespace std;
15 
16 static inline void setClusterRange(HelixRange& range1, HelixRange& range2,
17  ParRange& prange, unsigned int n_phi,
18  unsigned int n_d, unsigned int n_k,
19  unsigned int n_dzdl, unsigned int n_z0) {
20  float dzdl_size = (range1.max_dzdl - range1.min_dzdl) / ((float)n_dzdl);
21  float z0_size = (range1.max_z0 - range1.min_z0) / ((float)n_z0);
22  float phi_size = (range1.max_phi - range1.min_phi) / ((float)n_phi);
23  float d_size = (range1.max_d - range1.min_d) / ((float)n_d);
24  float k_size = (range1.max_k - range1.min_k) / ((float)n_k);
25 
26  range2.min_phi = range1.min_phi + phi_size * ((float)(prange.min_phi));
27  range2.max_phi = range1.min_phi + phi_size * ((float)(prange.max_phi + 1));
28  range2.min_d = range1.min_d + d_size * ((float)(prange.min_d));
29  range2.max_d = range1.min_d + d_size * ((float)(prange.max_d + 1));
30  range2.min_k = range1.min_k + k_size * ((float)(prange.min_k));
31  range2.max_k = range1.min_k + k_size * ((float)(prange.max_k + 1));
32  range2.min_dzdl = range1.min_dzdl + dzdl_size * ((float)(prange.min_dzdl));
33  range2.max_dzdl =
34  range1.min_dzdl + dzdl_size * ((float)(prange.max_dzdl + 1));
35  range2.min_z0 = range1.min_z0 + z0_size * ((float)(prange.min_z0));
36  range2.max_z0 = range1.min_z0 + z0_size * ((float)(prange.max_z0 + 1));
37 }
38 
39 // static inline void setRange(const BinEntryPair5D& bp, HelixRange& range1,
40 // HelixRange& range2, unsigned int n_phi,
41 // unsigned int n_d, unsigned int n_k,
42 // unsigned int n_dzdl, unsigned int n_z0) {
43 // float dzdl_size = (range1.max_dzdl - range1.min_dzdl) / ((float)n_dzdl);
44 // float z0_size = (range1.max_z0 - range1.min_z0) / ((float)n_z0);
45 // float phi_size = (range1.max_phi - range1.min_phi) / ((float)n_phi);
46 // float d_size = (range1.max_d - range1.min_d) / ((float)n_d);
47 // float k_size = (range1.max_k - range1.min_k) / ((float)n_k);
48 
49 // unsigned int z0_bin = 0;
50 // unsigned int dzdl_bin = 0;
51 // unsigned int k_bin = 0;
52 // unsigned int d_bin = 0;
53 // unsigned int phi_bin = 0;
54 
55 // bp.bin5D(n_d, n_k, n_dzdl, n_z0, phi_bin, d_bin, k_bin, dzdl_bin, z0_bin);
56 // range2.min_phi = range1.min_phi + phi_size * ((float)(phi_bin));
57 // range2.max_phi = range2.min_phi + phi_size;
58 // range2.min_d = range1.min_d + d_size * ((float)(d_bin));
59 // range2.max_d = range2.min_d + d_size;
60 // range2.min_k = range1.min_k + k_size * ((float)(k_bin));
61 // range2.max_k = range2.min_k + k_size;
62 // range2.min_dzdl = range1.min_dzdl + dzdl_size * ((float)(dzdl_bin));
63 // range2.max_dzdl = range2.min_dzdl + dzdl_size;
64 // range2.min_z0 = range1.min_z0 + z0_size * ((float)(z0_bin));
65 // range2.max_z0 = range2.min_z0 + z0_size;
66 // }
67 
68 void HelixHough::findHelicesByPairsBegin(unsigned int min_hits,
69  unsigned int max_hits,
70  vector<SimpleTrack3D>& tracks,
71  unsigned int maxtracks,
72  unsigned int zoomlevel) {
73  pairs_vec[zoomlevel]->clear();
74  pair<unsigned int, unsigned int> onepair;
75  // make list of all pairs from the hits in this zoomlevel
76  for (unsigned int i = 0; i < (hits_vec[zoomlevel]->size()); ++i) {
77  // if( (*(hits_vec[zoomlevel]))[i].get_layer() != 0 ){continue;}
78  onepair.first = i;
79  for (unsigned int j = 0; j < (hits_vec[zoomlevel]->size()); ++j) {
80  if (((*(hits_vec[zoomlevel]))[j].get_layer() <=
81  (*(hits_vec[zoomlevel]))[i].get_layer())) {
82  continue;
83  }
84 
85  // if( (*(hits_vec[zoomlevel]))[i].get_layer() != 4 ){continue;}
86  onepair.second = j;
87  pairs_vec[zoomlevel]->push_back(onepair);
88  }
89  }
90  findHelicesByPairs(min_hits, max_hits, tracks, maxtracks, zoomlevel);
91 }
92 
93 void HelixHough::findHelicesByPairs(unsigned int min_hits,
94  unsigned int max_hits,
95  vector<SimpleTrack3D>& tracks,
96  unsigned int maxtracks,
97  unsigned int zoomlevel) {
98  if ((maxtracks != 0) && (tracks.size() >= max_tracks)) {
99  return;
100  }
101 // unsigned int tracks_at_start = tracks.size();
102 
103  timeval t1, t2;
104  double time1 = 0.;
105  double time2 = 0.;
106  if (print_timings == true) {
107  gettimeofday(&t1, NULL);
108  }
109  vote_pairs(zoomlevel);
110  if (print_timings == true) {
111  gettimeofday(&t2, NULL);
112  time1 = ((double)(t1.tv_sec) + (double)(t1.tv_usec) / 1000000.);
113  time2 = ((double)(t2.tv_sec) + (double)(t2.tv_usec) / 1000000.);
114  vote_time += (time2 - time1);
115  }
116 
117  unsigned int n_phi = n_phi_bins[zoomlevel];
118  unsigned int n_d = n_d_bins[zoomlevel];
119  unsigned int n_k = n_k_bins[zoomlevel];
120  unsigned int n_dzdl = n_dzdl_bins[zoomlevel];
121  unsigned int n_z0 = n_z0_bins[zoomlevel];
122 
123  unsigned int n_entries = bins_vec[zoomlevel]->size();
124  if (n_entries == 0) {
125  return;
126  }
127 
128  if (print_timings == true) {
129  gettimeofday(&t1, NULL);
130  }
131  num_clusters[zoomlevel] = 0;
132  bool use_clusters = true;
133  bool is_super_bin = false;
134  makeClusters(zoomlevel, pairs_vec[zoomlevel]->size(), n_phi, n_d, n_k, n_dzdl,
135  n_z0, min_hits, *(clusters_vec[zoomlevel]), use_clusters,
136  is_super_bin);
137  if (print_timings) {
138  gettimeofday(&t2, NULL);
139  time1 = ((double)(t1.tv_sec) + (double)(t1.tv_usec) / 1000000.);
140  time2 = ((double)(t2.tv_sec) + (double)(t2.tv_usec) / 1000000.);
141  cluster_time += (time2 - time1);
142  }
143 
144  pair<unsigned int, unsigned int> onepair;
145  for (unsigned int i = 0, size = num_clusters[zoomlevel]; i < size; ++i) {
146  temp_pairs.clear();
147  new_hits.clear();
148  old_to_new.clear();
149  hits_vec[zoomlevel + 1]->clear();
150  pairs_vec[zoomlevel + 1]->clear();
151  vector<unsigned int>::iterator index_iter;
152  for (index_iter = (*(clusters_vec[zoomlevel]))[i].hit_indexes.begin();
153  index_iter != (*(clusters_vec[zoomlevel]))[i].hit_indexes.end();
154  ++index_iter) {
155  if (remove_hits == true) {
156  if (((*hit_used)[(*(hits_vec[zoomlevel]))
157  [(*(pairs_vec[zoomlevel]))[*index_iter].first]
158  .get_id()] == false) &&
159  ((*hit_used)[(*(hits_vec[zoomlevel]))
160  [(*(pairs_vec[zoomlevel]))[*index_iter].second]
161  .get_id()] == false)) {
162  temp_pairs.push_back((*(pairs_vec[zoomlevel]))[*index_iter]);
163  new_hits.insert((temp_pairs.back()).first);
164  new_hits.insert((temp_pairs.back()).second);
165  }
166  } else {
167  temp_pairs.push_back((*(pairs_vec[zoomlevel]))[*index_iter]);
168  new_hits.insert((temp_pairs.back()).first);
169  new_hits.insert((temp_pairs.back()).second);
170  }
171  }
172 
173  unsigned int temp_index = 0;
174  set<unsigned int>::iterator it;
175  for (it = new_hits.begin(); it != new_hits.end(); ++it) {
176  hits_vec[zoomlevel + 1]->push_back((*(hits_vec[zoomlevel]))[*it]);
177  old_to_new[*it] = temp_index;
178  temp_index += 1;
179  }
180  for (unsigned int j = 0, size = temp_pairs.size(); j != size; ++j) {
181  onepair.first = old_to_new[temp_pairs[j].first];
182  onepair.second = old_to_new[temp_pairs[j].second];
183  pairs_vec[zoomlevel + 1]->push_back(onepair);
184  }
185  setClusterRange(zoomranges[zoomlevel], zoomranges[zoomlevel + 1],
186  (*(clusters_vec[zoomlevel]))[i].range, n_phi, n_d, n_k,
187  n_dzdl, n_z0);
188  if ((breakRecursion(*(hits_vec[zoomlevel + 1]),
189  zoomranges[zoomlevel + 1]) == true)) {
190  } else if (((zoomlevel + 1) == max_zoom) ||
191  (hits_vec[zoomlevel + 1]->size() <= max_hits)) {
192  // findTracksByPairs(*(hits_vec[zoomlevel+1]),
193  // *(pairs_vec[zoomlevel+1]), tracks, zoomranges[zoomlevel+1]);
194  findTracks(*(hits_vec[zoomlevel + 1]), tracks, zoomranges[zoomlevel + 1]);
195  } else {
196  findHelicesByPairs(min_hits, max_hits, tracks, maxtracks, zoomlevel + 1);
197  }
198  }
199 }