33 vector<unsigned int>&
C) {
34 unsigned int SIZE = A.size();
35 unsigned int current = 0;
36 for (
unsigned int i = 0; i < SIZE; ++i) {
48 unsigned int n_d,
unsigned int n_k,
49 unsigned int n_dzdl,
unsigned int n_z0) {
51 float z0_size = (range1.
max_z0 - range1.
min_z0) / ((
float)n_z0);
52 float phi_size = (range1.
max_phi - range1.
min_phi) / ((
float)n_phi);
53 float d_size = (range1.
max_d - range1.
min_d) / ((
float)n_d);
54 float k_size = (range1.
max_k - range1.
min_k) / ((
float)n_k);
56 unsigned int z0_bin = 0;
57 unsigned int dzdl_bin = 0;
58 unsigned int k_bin = 0;
59 unsigned int d_bin = 0;
60 unsigned int phi_bin = 0;
62 bp.
bin5D(n_d, n_k, n_dzdl, n_z0, phi_bin, d_bin, k_bin, dzdl_bin, z0_bin);
65 range2.
min_d = range1.
min_d + d_size * ((float)(d_bin));
67 range2.
min_k = range1.
min_k + k_size * ((float)(k_bin));
71 range2.
min_z0 = range1.
min_z0 + z0_size * ((float)(z0_bin));
76 ParRange& prange,
unsigned int n_phi,
77 unsigned int n_d,
unsigned int n_k,
78 unsigned int n_dzdl,
unsigned int n_z0) {
80 float z0_size = (range1.
max_z0 - range1.
min_z0) / ((
float)n_z0);
81 float phi_size = (range1.
max_phi - range1.
min_phi) / ((
float)n_phi);
82 float d_size = (range1.
max_d - range1.
min_d) / ((
float)n_d);
83 float k_size = (range1.
max_k - range1.
min_k) / ((
float)n_k);
99 unsigned int max_hits,
100 vector<SimpleTrack3D>& tracks,
101 unsigned int maxtracks) {
107 vector<SimpleHit3D>
hits;
109 if (layer_end > 0 && layer_start > 0){
110 for (
unsigned int in = 0;
in < hits_init.size();
in++) {
111 if(hits_init[
in].get_layer() > layer_start-1 && hits_init[
in].get_layer() < layer_end+1 )
112 hits.push_back(hits_init[
in]);
114 }
else if ( layer_end < 0 && layer_start > 0) {
115 for (
unsigned int in = 0;
in < hits_init.size();
in++) {
116 if(hits_init[
in].get_layer() > layer_start-1)
117 hits.push_back(hits_init[
in]);
119 }
else if ( layer_end > 0 && layer_start < 0) {
120 for (
unsigned int in = 0;
in < hits_init.size();
in++) {
121 if( hits_init[
in].get_layer() < layer_end+1 )
122 hits.push_back(hits_init[
in]);
125 for (
unsigned int in = 0;
in < hits_init.size();
in++) {
126 hits.push_back(hits_init[
in]);
130 index_mapping.clear();
131 index_mapping.resize(hits.size(), 0);
133 for (
unsigned int i = 0; i < hits.size(); i++) {
134 index_mapping[i] = hits[i].get_id();
137 (*hit_used).assign(hits.size(),
false);
139 initEvent(hits, min_hits);
141 max_tracks = maxtracks;
143 (*(hits_vec[start_zoom])) =
hits;
144 current_range = top_range;
146 for (
unsigned int z = 0;
z <= max_zoom;
z++) {
147 zoomranges.push_back(top_range);
149 vector<SimpleTrack3D> temp_tracks;
150 if ((separate_by_helicity ==
true) && (only_one_helicity ==
false)) {
152 findHelices(min_hits, max_hits, temp_tracks, maxtracks, start_zoom);
154 findHelices(min_hits, max_hits, temp_tracks, maxtracks, start_zoom);
156 findHelices(min_hits, max_hits, temp_tracks, maxtracks, start_zoom);
159 vector<SimpleHit3D> tr_hits;
160 if (cull_input_hits ==
true) {
161 for (
unsigned int i = 0; i < hits.size(); i++) {
162 if ((*hit_used)[hits[i].get_id()] ==
false) {
163 tr_hits.push_back(hits[i]);
164 tr_hits.back().set_id(index_mapping[i]);
169 for (
unsigned int i = 0; i < hits.size(); i++) {
170 hits[i].set_id(index_mapping[i]);
172 for (
unsigned int t = 0;
t < temp_tracks.size();
t++) {
173 for (
unsigned int h = 0;
h < temp_tracks[
t].hits.size();
h++) {
174 if (temp_tracks[
t].hits[
h].get_id() != (unsigned)-1) {
175 temp_tracks[
t].hits[
h].set_id(index_mapping[temp_tracks[
t].hits[
h].get_id()]);
180 finalize(temp_tracks, tracks);
182 if (cull_input_hits ==
true) {
186 if (print_timings ==
true) {
187 cout <<
"vote time = " << vote_time << endl;
188 cout <<
"xy vote time = " << xy_vote_time << endl;
189 cout <<
"z vote time = " << z_vote_time << endl;
190 cout <<
"cluster time = " << cluster_time << endl;
197 vector<unsigned int> old_indexes = clus1.
hit_indexes;
198 for (
unsigned int i = 0; i < clus2.
hit_indexes.size(); ++i) {
201 vector<unsigned int>
C;
202 C.assign(MAX + 1, 0);
205 float size_diff_1 = ((float)(clus1.
hit_indexes.size() - old_indexes.size()) /
206 ((
float)(old_indexes.size())));
211 if ((size_diff_1 < overlap_cut) || (size_diff_2 < overlap_cut)) {
221 unsigned int zoomlevel,
unsigned int MAX,
unsigned int ca,
unsigned int d,
222 unsigned int r,
unsigned int th,
unsigned int zz0,
unsigned int bin,
223 unsigned int newbin, vector<unsigned char>& good_bins,
unsigned int volume,
224 float cluster_size_cut,
float overlap_cut,
225 vector<ParameterCluster>& clusters,
unsigned int* bins_start,
226 unsigned int* bins_end, vector<unsigned int>& map_clus,
227 vector<unsigned char>& too_big, vector<unsigned int>& ,
228 vector<unsigned int>&
C) {
229 if ((good_bins[newbin] == 1) && (map_clus[newbin] != 4294967295)) {
230 if (too_big[map_clus[newbin]] == 0) {
231 unsigned int tempsize = clusters[map_clus[newbin]].hit_indexes.size();
232 for (
unsigned int ind = bins_start[bin]; ind <= bins_end[
bin]; ++ind) {
233 clusters[map_clus[newbin]].hit_indexes.push_back(
234 (*(bins_vec[zoomlevel]))[ind].entry);
237 C.assign(MAX + 1, 0);
239 unsigned int size_diff =
240 clusters[map_clus[newbin]].hit_indexes.size() - tempsize;
241 unsigned int size_diff_2 = clusters[map_clus[newbin]].hit_indexes.size() -
242 (1 + bins_end[
bin] - bins_start[
bin]);
243 if ((((
float)size_diff) / ((
float)(1 + bins_end[bin] - bins_start[bin])) <
245 (((float)size_diff_2) / ((float)(tempsize)) < overlap_cut)) {
246 clusters[map_clus[newbin]].range.mergeRange(ca, d, r, th, zz0);
248 unsigned int cluster_volume =
254 if (((
float)cluster_volume) / ((
float)
volume) > cluster_size_cut) {
255 too_big[map_clus[newbin]] = 1;
257 map_clus[
bin] = map_clus[newbin];
260 clusters[map_clus[newbin]].hit_indexes.resize(tempsize);
268 unsigned int n_phi,
unsigned int n_d,
269 unsigned int n_k,
unsigned int n_dzdl,
270 unsigned int n_z0,
unsigned int min_hits,
271 vector<ParameterCluster>& clusters,
272 bool& ,
bool& is_super_bin) {
273 unsigned int volume = n_phi * n_d * n_k * n_dzdl * n_z0;
274 float cluster_size_cut = 1.0;
276 float overlap_cut = 0.1;
277 is_super_bin =
false;
279 vector<unsigned int> map_clus(volume, 4294967295);
280 vector<unsigned char> good_bins(volume, 0);
281 vector<unsigned char> too_big(volume, 0);
283 for (
unsigned int ca = 0; ca < n_phi; ++ca) {
284 for (
unsigned int d = 0;
d < n_d; ++
d) {
285 for (
unsigned int r = 0;
r < n_k; ++
r) {
286 for (
unsigned int th = 0; th < n_dzdl; ++th) {
287 for (
unsigned int zz0 = 0; zz0 < n_z0; ++zz0) {
291 if (bins_end[bin] == 4294967295) {
294 if ((1 + bins_end[bin] - bins_start[bin]) >= min_hits) {
295 if (check_layers ==
true) {
296 unsigned int layer_mask[4] = {0, 0, 0, 0};
297 for (
unsigned int i = bins_start[bin]; i <= bins_end[
bin];
299 if ((*(hits_vec[zoomlevel]))[(*(bins_vec[zoomlevel]))[i]
305 hits_vec[zoomlevel]))[(*(bins_vec[zoomlevel]))[i]
308 }
else if ((*(hits_vec[zoomlevel]))
309 [(*(bins_vec[zoomlevel]))[i].entry]
313 (1 << ((*(hits_vec[zoomlevel]))
314 [(*(bins_vec[zoomlevel]))[i].entry]
317 }
else if ((*(hits_vec[zoomlevel]))
318 [(*(bins_vec[zoomlevel]))[i].entry]
322 (1 << ((*(hits_vec[zoomlevel]))
323 [(*(bins_vec[zoomlevel]))[i].entry]
326 }
else if ((*(hits_vec[zoomlevel]))
327 [(*(bins_vec[zoomlevel]))[i].entry]
328 .get_layer() < 128) {
331 (1 << ((*(hits_vec[zoomlevel]))
332 [(*(bins_vec[zoomlevel]))[i].entry]
337 unsigned int nlayers = __builtin_popcount(layer_mask[0]) +
338 __builtin_popcount(layer_mask[1]) +
339 __builtin_popcount(layer_mask[2]) +
340 __builtin_popcount(layer_mask[3]);
341 if (nlayers >= req_layers) {
351 if (good_bins[bin] == 0) {
356 unsigned int newbin = bin - n_d * n_k * n_dzdl * n_z0;
357 if (attemptClusterMerge(
358 zoomlevel, MAX, ca,
d,
r, th, zz0, bin, newbin, good_bins,
359 volume, cluster_size_cut, overlap_cut, clusters,
360 (
unsigned int*)(bins_start), (
unsigned int*)(bins_end),
361 map_clus, too_big, temp_merged_clus, C_clus) ==
true) {
367 unsigned int newbin = bin - n_k * n_dzdl * n_z0;
368 if (attemptClusterMerge(
369 zoomlevel, MAX, ca,
d,
r, th, zz0, bin, newbin, good_bins,
370 volume, cluster_size_cut, overlap_cut, clusters,
371 (
unsigned int*)(bins_start), (
unsigned int*)(bins_end),
372 map_clus, too_big, temp_merged_clus, C_clus) ==
true) {
378 unsigned int newbin = bin - n_dzdl * n_z0;
379 if (attemptClusterMerge(
380 zoomlevel, MAX, ca,
d,
r, th, zz0, bin, newbin, good_bins,
381 volume, cluster_size_cut, overlap_cut, clusters,
382 (
unsigned int*)(bins_start), (
unsigned int*)(bins_end),
383 map_clus, too_big, temp_merged_clus, C_clus) ==
true) {
389 unsigned int newbin = bin - n_z0;
390 if (attemptClusterMerge(
391 zoomlevel, MAX, ca,
d,
r, th, zz0, bin, newbin, good_bins,
392 volume, cluster_size_cut, overlap_cut, clusters,
393 (
unsigned int*)(bins_start), (
unsigned int*)(bins_end),
394 map_clus, too_big, temp_merged_clus, C_clus) ==
true) {
400 unsigned int newbin = bin - 1;
401 if (attemptClusterMerge(
402 zoomlevel, MAX, ca,
d,
r, th, zz0, bin, newbin, good_bins,
403 volume, cluster_size_cut, overlap_cut, clusters,
404 (
unsigned int*)(bins_start), (
unsigned int*)(bins_end),
405 map_clus, too_big, temp_merged_clus, C_clus) ==
true) {
409 if (num_clusters[zoomlevel] >= clusters.size()) {
412 clusters[num_clusters[zoomlevel]].range =
414 map_clus[
bin] = num_clusters[zoomlevel];
415 clusters[num_clusters[zoomlevel]].hit_indexes.clear();
416 for (
unsigned int ind = bins_start[bin]; ind <= bins_end[
bin];
418 clusters[num_clusters[zoomlevel]].hit_indexes.push_back(
419 (*(bins_vec[zoomlevel]))[ind].entry);
422 num_clusters[zoomlevel] += 1;
428 if (iterate_clustering ==
false) {
431 if (num_clusters[zoomlevel] == 0) {
434 vector<ParameterCluster> in_clusters;
435 for (
unsigned int i = 0; i < num_clusters[zoomlevel]; ++i) {
436 in_clusters.push_back(clusters[i]);
438 vector<ParameterCluster> out_clusters;
439 out_clusters.push_back(in_clusters[0]);
440 for (
unsigned int i = 1; i < in_clusters.size(); ++i) {
442 for (
unsigned int j = 0; j < out_clusters.size(); ++j) {
443 merged =
mergeClusters(out_clusters[j], in_clusters[i], MAX, overlap_cut);
444 if (merged ==
true) {
449 if (merged ==
false) {
450 out_clusters.push_back(in_clusters[i]);
453 clusters = out_clusters;
454 num_clusters[zoomlevel] = out_clusters.size();
458 vector<SimpleTrack3D>& tracks,
459 unsigned int maxtracks,
unsigned int zoomlevel) {
460 unsigned int tracks_at_start = tracks.size();
462 if ((maxtracks != 0) && (tracks.size() >= max_tracks)) {
469 if (print_timings ==
true) {
470 gettimeofday(&t1, NULL);
473 if (print_timings ==
true) {
474 gettimeofday(&t2, NULL);
475 time1 = ((double)(t1.tv_sec) + (double)(t1.tv_usec) / 1000000.);
476 time2 = ((double)(t2.tv_sec) + (double)(t2.tv_usec) / 1000000.);
477 vote_time += (time2 - time1);
480 unsigned int n_entries = bins_vec[zoomlevel]->size();
481 if (n_entries == 0) {
485 unsigned int n_phi = n_phi_bins[zoomlevel];
486 unsigned int n_d = n_d_bins[zoomlevel];
487 unsigned int n_k = n_k_bins[zoomlevel];
488 unsigned int n_dzdl = n_dzdl_bins[zoomlevel];
489 unsigned int n_z0 = n_z0_bins[zoomlevel];
491 num_clusters[zoomlevel] = 0;
492 bool use_clusters =
true;
493 bool is_super_bin =
false;
494 if (zoomlevel >= cluster_start_bin) {
495 if (print_timings ==
true) {
496 gettimeofday(&t1, NULL);
498 makeClusters(zoomlevel, hits_vec[zoomlevel]->size(), n_phi, n_d, n_k,
499 n_dzdl, n_z0, min_hits, *(clusters_vec[zoomlevel]),
500 use_clusters, is_super_bin);
503 gettimeofday(&t2, NULL);
504 time1 = ((double)(t1.tv_sec) + (double)(t1.tv_usec) / 1000000.);
505 time2 = ((double)(t2.tv_sec) + (double)(t2.tv_usec) / 1000000.);
506 cluster_time += (time2 - time1);
509 use_clusters =
false;
512 if (use_clusters ==
false) {
513 unsigned int count = 0;
514 hits_vec[zoomlevel + 1]->clear();
515 setRange((*(bins_vec[zoomlevel]))[count], zoomranges[zoomlevel],
516 zoomranges[zoomlevel + 1], n_phi, n_d, n_k, n_dzdl, n_z0);
518 while (count < n_entries) {
519 if (remove_hits ==
true) {
521 hits_vec[zoomlevel]))[(*(bins_vec[zoomlevel]))[count].entry]
522 .get_id()] ==
false) {
523 hits_vec[zoomlevel + 1]->push_back(
524 (*(hits_vec[zoomlevel]))[(*(bins_vec[zoomlevel]))[count].entry]);
527 hits_vec[zoomlevel + 1]->push_back(
528 (*(hits_vec[zoomlevel]))[(*(bins_vec[zoomlevel]))[count].entry]);
534 if ((count == n_entries) || ((*(bins_vec[zoomlevel]))[count].
bin !=
535 (*(bins_vec[zoomlevel]))[count - 1].
bin)) {
536 if (hits_vec[zoomlevel + 1]->size() >= min_hits) {
537 if (breakRecursion(*(hits_vec[zoomlevel + 1]),
538 zoomranges[zoomlevel + 1]) ==
true) {
539 }
else if (((zoomlevel + 1) == max_zoom)) {
540 findTracks(*(hits_vec[zoomlevel + 1]), tracks,
541 zoomranges[zoomlevel + 1]);
542 }
else if (((zoomlevel + 1) >= min_zoom) &&
543 ((hits_vec[zoomlevel + 1]->size() <= max_hits))) {
544 findTracks(*(hits_vec[zoomlevel + 1]), tracks,
545 zoomranges[zoomlevel + 1]);
546 }
else if ((hits_vec[zoomlevel + 1]->size() <= max_hits_pairs) &&
547 ((zoomlevel + 1) >= min_zoom)) {
548 findHelicesByPairsBegin(min_hits, max_hits, tracks, maxtracks,
551 findHelices(min_hits, max_hits, tracks, maxtracks, zoomlevel + 1);
553 if (maxtracks != 0) {
554 double phi_proportion =
555 (zoomranges[zoomlevel].max_phi -
556 zoomranges[zoomlevel].min_phi) /
557 ((zoomranges[0].max_phi - zoomranges[0].min_phi));
558 double d_proportion =
559 (zoomranges[zoomlevel].max_d - zoomranges[zoomlevel].min_d) /
560 ((zoomranges[0].max_d - zoomranges[0].min_d));
561 double k_proportion =
562 (zoomranges[zoomlevel].max_k - zoomranges[zoomlevel].min_k) /
563 ((zoomranges[0].max_k - zoomranges[0].min_k));
564 unsigned int expected_tracks =
565 (
unsigned int)(fabs(((
double)maxtracks) * phi_proportion *
566 d_proportion * k_proportion)) +
569 if ((tracks.size() - tracks_at_start) > expected_tracks) {
574 if (count == n_entries) {
577 hits_vec[zoomlevel + 1]->clear();
578 setRange((*(bins_vec[zoomlevel]))[count], zoomranges[zoomlevel],
579 zoomranges[zoomlevel + 1], n_phi, n_d, n_k, n_dzdl, n_z0);
583 if (clusters_vec[zoomlevel]->size() == 0) {
588 for (
unsigned int i = 0, size = num_clusters[zoomlevel]; i < size; ++i) {
589 hits_vec[zoomlevel + 1]->clear();
590 vector<unsigned int>::iterator index_iter;
591 for (index_iter = (*(clusters_vec[zoomlevel]))[i].hit_indexes.begin();
592 index_iter != (*(clusters_vec[zoomlevel]))[i].hit_indexes.end();
594 if (remove_hits ==
true) {
595 if ((*hit_used)[(*(hits_vec[zoomlevel]))[*index_iter].get_id()] ==
597 hits_vec[zoomlevel + 1]->push_back(
598 (*(hits_vec[zoomlevel]))[*index_iter]);
601 hits_vec[zoomlevel + 1]->push_back(
602 (*(hits_vec[zoomlevel]))[*index_iter]);
607 (*(clusters_vec[zoomlevel]))[i].range, n_phi, n_d, n_k,
609 if ((breakRecursion(*(hits_vec[zoomlevel + 1]),
610 zoomranges[zoomlevel + 1]) ==
true)) {
611 }
else if ((zoomlevel + 1) == max_zoom) {
612 findTracks(*(hits_vec[zoomlevel + 1]), tracks,
613 zoomranges[zoomlevel + 1]);
614 }
else if (((zoomlevel + 1) >= min_zoom) &&
615 (hits_vec[zoomlevel + 1]->size() <= max_hits)) {
616 findTracks(*(hits_vec[zoomlevel + 1]), tracks,
617 zoomranges[zoomlevel + 1]);
618 }
else if ((hits_vec[zoomlevel + 1]->size() <= max_hits_pairs) &&
619 ((zoomlevel + 1) >= min_zoom)) {
620 findHelicesByPairsBegin(min_hits, max_hits, tracks, maxtracks,
623 findHelices(min_hits, max_hits, tracks, maxtracks, zoomlevel + 1);
630 vector<SimpleHit3D>&
hits,
631 unsigned int min_hits,
unsigned int max_hits,
632 vector<SimpleTrack3D>& tracks,
633 unsigned int maxtracks) {
634 unsigned int n_layers_orig = n_layers;
635 unsigned int req_layers_orig = req_layers;
637 bool smooth_back_orig = smooth_back;
640 int min_layer = 999999;
641 for (
unsigned int i = 0; i < hits.size(); ++i) {
642 if (hits[i].get_layer() > max_layer) {
643 max_layer = hits[i].get_layer();
645 if (hits[i].get_layer() < min_layer) {
646 min_layer = hits[i].get_layer();
648 if (hits[i].get_layer() < 0) {
649 findSeededHelices_run(seeds, hits, min_hits, max_hits, tracks, maxtracks);
654 if (min_layer == 999999) {
655 findSeededHelices_run(seeds, hits, min_hits, max_hits, tracks, maxtracks);
659 n_layers = min_layer;
661 float req_prop = ((float)req_layers) / ((float)((1 + max_layer - min_layer)));
663 int layers_left = (1 + max_layer - min_layer);
664 unsigned int cur_layer = min_layer;
665 vector<SimpleTrack3D> cur_seeds = seeds;
666 vector<SimpleTrack3D> cur_tracks;
667 vector<SimpleHit3D> cur_hits;
669 while (layers_left > 0) {
670 unsigned int layers_iter = layers_at_a_time;
671 if ((cur_layer + layers_at_a_time - 1) > (
unsigned int) max_layer) {
672 layers_iter = (max_layer - (cur_layer - 1));
677 for (
unsigned int i = 0; i < hits.size(); ++i) {
678 if ((hits[i].get_layer() >= (
int) cur_layer) &&
679 (hits[i].get_layer() <= (
int) (cur_layer + layers_at_a_time - 1))) {
680 cur_hits.push_back(hits[i]);
684 n_layers += layers_iter;
686 unsigned int nreq = floor(((
float)(layers_iter)) * req_prop);
695 if (layers_left <= (
int) layers_at_a_time) {
696 smooth_back = smooth_back_orig;
699 findSeededHelices_run(cur_seeds, cur_hits, nreq, max_hits, cur_tracks,
701 setSeedStates(track_states);
702 cur_seeds = cur_tracks;
703 layers_left -= layers_at_a_time;
704 cur_layer += layers_at_a_time;
708 n_layers = n_layers_orig;
709 req_layers = req_layers_orig;
713 vector<SimpleHit3D>&
hits,
714 unsigned int min_hits,
715 unsigned int max_hits,
716 vector<SimpleTrack3D>& tracks,
717 unsigned int maxtracks) {
722 initEvent(hits, min_hits);
726 index_mapping.clear();
727 index_mapping.resize(hits.size(), 0);
728 for (
unsigned int i = 0; i < hits.size(); i++) {
729 index_mapping[i] = hits[i].get_id();
732 (*hit_used).assign(hits.size(),
false);
734 max_tracks = maxtracks;
736 (*(hits_vec[0])) =
hits;
737 (*(seeds_vec[0])) = seeds;
738 current_range = top_range;
740 for (
unsigned int z = 0;
z <= max_zoom;
z++) {
741 zoomranges.push_back(top_range);
743 vector<SimpleTrack3D> temp_tracks;
745 if (separate_by_helicity ==
true) {
747 findSeededHelices(min_hits, max_hits, temp_tracks, maxtracks, 0);
749 findSeededHelices(min_hits, max_hits, temp_tracks, maxtracks, 0);
751 findSeededHelices(min_hits, max_hits, temp_tracks, maxtracks, 0);
754 for (
unsigned int i = 0; i < hits.size(); i++) {
755 hits[i].set_id(index_mapping[i]);
758 finalize(temp_tracks, tracks);
760 if (print_timings ==
true) {
761 cout <<
"vote time = " << vote_time << endl;
762 cout <<
"xy vote time = " << xy_vote_time << endl;
763 cout <<
"z vote time = " << z_vote_time << endl;
779 static int seed_bin(vector<floatBin>& bins,
float val) {
781 if ((bin < bins[0]) || (bins.back() <
bin)) {
784 pair<vector<floatBin>::iterator, vector<floatBin>::iterator> bounds;
785 bounds = equal_range(bins.begin(), bins.end(),
bin);
786 return ((
int)(bounds.first - bins.begin()));
798 range.
min_d = seed.
d - 0.001;
799 range.
max_d = seed.
d + 0.001;
802 if (range.
min_k < 0.) {
816 vector<SimpleTrack3D>& tracks,
817 unsigned int maxtracks,
818 unsigned int zoomlevel) {
819 unsigned int tracks_at_start = tracks.size();
821 if ((maxtracks != 0) && (tracks.size() >= max_tracks)) {
828 gettimeofday(&t1, NULL);
830 gettimeofday(&t2, NULL);
831 time1 = ((double)(t1.tv_sec) + (double)(t1.tv_usec) / 1000000.);
832 time2 = ((double)(t2.tv_sec) + (double)(t2.tv_usec) / 1000000.);
833 vote_time += (time2 - time1);
835 unsigned int n_phi = n_phi_bins[zoomlevel];
836 unsigned int n_d = n_d_bins[zoomlevel];
837 unsigned int n_k = n_k_bins[zoomlevel];
838 unsigned int n_dzdl = n_dzdl_bins[zoomlevel];
839 unsigned int n_z0 = n_z0_bins[zoomlevel];
846 low = zoomranges[zoomlevel].min_phi;
847 delta = (zoomranges[zoomlevel].max_phi - zoomranges[zoomlevel].min_phi) /
850 for (
unsigned int i = 0; i < n_phi; ++i) {
851 phibins.push_back(
floatBin(low, high));
856 vector<floatBin> dbins;
857 low = zoomranges[zoomlevel].min_d;
858 delta = (zoomranges[zoomlevel].max_d - zoomranges[zoomlevel].min_d) /
861 for (
unsigned int i = 0; i < n_d; ++i) {
862 dbins.push_back(
floatBin(low, high));
867 vector<floatBin> kbins;
868 low = zoomranges[zoomlevel].min_k;
869 delta = (zoomranges[zoomlevel].max_k - zoomranges[zoomlevel].min_k) /
872 for (
unsigned int i = 0; i < n_k; ++i) {
873 kbins.push_back(
floatBin(low, high));
878 vector<floatBin> z0bins;
879 low = zoomranges[zoomlevel].min_z0;
880 delta = (zoomranges[zoomlevel].max_z0 - zoomranges[zoomlevel].min_z0) /
883 for (
unsigned int i = 0; i < n_z0; ++i) {
884 z0bins.push_back(
floatBin(low, high));
889 vector<floatBin> dzdlbins;
890 low = zoomranges[zoomlevel].min_dzdl;
891 delta = (zoomranges[zoomlevel].max_dzdl - zoomranges[zoomlevel].min_dzdl) /
894 for (
unsigned int i = 0; i < n_dzdl; ++i) {
895 dzdlbins.push_back(
floatBin(low, high));
901 for (
unsigned int i = 0; i < seeds_vec[zoomlevel]->size(); ++i) {
902 float d = (*(seeds_vec[zoomlevel]))[i].d;
903 float phi = (*(seeds_vec[zoomlevel]))[i].phi;
904 float kappa = pow((*(seeds_vec[zoomlevel]))[i].kappa, 1.);
905 float z0 = (*(seeds_vec[zoomlevel]))[i].z0;
906 float dzdl = (*(seeds_vec[zoomlevel]))[i].dzdl;
910 unsigned int phi_bin = 0;
915 phi_bin = (
unsigned int)(tempbin);
918 unsigned int d_bin = 0;
923 d_bin = (
unsigned int)(tempbin);
926 unsigned int kappa_bin = 0;
931 kappa_bin = (
unsigned int)(tempbin);
934 unsigned int z0_bin = 0;
939 z0_bin = (
unsigned int)(tempbin);
942 unsigned int dzdl_bin = 0;
947 dzdl_bin = (
unsigned int)(tempbin);
951 n_d, n_k, n_dzdl, n_z0, phi_bin, d_bin, kappa_bin, dzdl_bin, z0_bin);
955 sort(bins_vec[zoomlevel]->begin(), bins_vec[zoomlevel]->end());
957 unsigned int n_entries = bins_vec[zoomlevel]->size();
958 if (n_entries == 0) {
962 unsigned int count = 0;
963 hits_vec[zoomlevel + 1]->clear();
964 seeds_vec[zoomlevel + 1]->clear();
965 setRange((*(bins_vec[zoomlevel]))[count], zoomranges[zoomlevel],
966 zoomranges[zoomlevel + 1], n_phi, n_d, n_k, n_dzdl, n_z0);
968 while (count < n_entries) {
969 if ((*(bins_vec[zoomlevel]))[count].is_seed ==
false) {
970 if (remove_hits ==
true) {
972 hits_vec[zoomlevel]))[(*(bins_vec[zoomlevel]))[count].entry]
973 .get_id()] ==
false) {
974 hits_vec[zoomlevel + 1]->push_back(
975 (*(hits_vec[zoomlevel]))[(*(bins_vec[zoomlevel]))[count].entry]);
978 hits_vec[zoomlevel + 1]->push_back(
979 (*(hits_vec[zoomlevel]))[(*(bins_vec[zoomlevel]))[count].entry]);
982 seeds_vec[zoomlevel + 1]->push_back(
983 (*(seeds_vec[zoomlevel]))[(*(bins_vec[zoomlevel]))[count].entry]);
990 if ((count == n_entries) || ((*(bins_vec[zoomlevel]))[count].
bin !=
991 (*(bins_vec[zoomlevel]))[count - 1].
bin)) {
992 if ((hits_vec[zoomlevel + 1]->size() >= min_hits) &&
993 (seeds_vec[zoomlevel + 1]->size() != 0)) {
994 if (breakRecursion(*(hits_vec[zoomlevel + 1]),
995 zoomranges[zoomlevel + 1]) ==
true) {
996 findSeededTracks(*(seeds_vec[zoomlevel + 1]),
997 *(hits_vec[zoomlevel + 1]), tracks,
998 zoomranges[zoomlevel + 1]);
999 }
else if (((zoomlevel + 1) >= min_zoom) &&
1000 ((hits_vec[zoomlevel + 1]->size() <= max_hits) ||
1001 (zoomlevel + 1) == max_zoom)) {
1002 findSeededTracks(*(seeds_vec[zoomlevel + 1]),
1003 *(hits_vec[zoomlevel + 1]), tracks,
1004 zoomranges[zoomlevel + 1]);
1006 findSeededHelices(min_hits, max_hits, tracks, maxtracks,
1009 if (maxtracks != 0) {
1010 double phi_proportion =
1011 (zoomranges[zoomlevel].max_phi - zoomranges[zoomlevel].min_phi) /
1012 ((zoomranges[0].max_phi - zoomranges[0].min_phi));
1013 double d_proportion =
1014 (zoomranges[zoomlevel].max_d - zoomranges[zoomlevel].min_d) /
1015 ((zoomranges[0].max_d - zoomranges[0].min_d));
1016 double k_proportion =
1017 (zoomranges[zoomlevel].max_k - zoomranges[zoomlevel].min_k) /
1018 ((zoomranges[0].max_k - zoomranges[0].min_k));
1019 unsigned int expected_tracks =
1020 (
unsigned int)(fabs(((
double)maxtracks) * phi_proportion *
1021 d_proportion * k_proportion)) +
1024 if ((tracks.size() - tracks_at_start) > expected_tracks) {
1029 if (count == n_entries) {
1032 hits_vec[zoomlevel + 1]->clear();
1033 seeds_vec[zoomlevel + 1]->clear();
1034 setRange((*(bins_vec[zoomlevel]))[count], zoomranges[zoomlevel],
1035 zoomranges[zoomlevel + 1], n_phi, n_d, n_k, n_dzdl, n_z0);