ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
PHTpcSeedFinder.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file PHTpcSeedFinder.cc
1 
7 #include "PHTpcSeedFinder.h"
8 #include "PHTpcTrackerUtil.h"
9 
10 #include "externals/kdfinder.hpp" // for TrackCandidate, Circle, find_tr...
11 
12 #include <phool/PHLog.h>
13 
14 #include <log4cpp/CategoryStream.hh> // for CategoryStream
15 
18 
19 float round(float var)
20 {
21  return (int(var * 100.0) / 100.0);
22 }
23 
25  : mMaxDistance1(3.0)
26  , mTripletAngle1(M_PI / 8)
27  , mMinHits1(10)
28  , mMaxDistance2(6.0)
29  , mTripletAngle2(M_PI / 8)
30  , mMinHits2(5)
31  , mNThreads(1)
32  , mRemoveLoopers(false)
33  , mMinLooperRadius(10.0)
34  , mMaxLooperRadius(70.0)
35 {
36 }
37 
38 std::vector<kdfinder::TrackCandidate<double>*> PHTpcSeedFinder::findSeeds(TrkrClusterContainer* cluster_map, TrkrHitSetContainer* hitsets, double B)
39 {
40  LOG_WARN_IF("tracking.PHTpcSeedFnder.findSeeds", !cluster_map) << __FILE__ << "," << __LINE__ << " Can't find node TRKR_CLUSTER";
41  std::vector<kdfinder::TrackCandidate<double>*> result;
42  if (!cluster_map)
43  {
44  return result;
45  }
46  LOG_WARN_IF("tracking.PHTpcSeedFnder.findSeeds", !hitsets) << __FILE__ << "," << __LINE__ << " Can't find node TRKR_HITSET";
47  if (!hitsets)
48  {
49  return result;
50  }
51 
52  //----- convert clusters to kdhits -----
53  std::vector<std::vector<double> > kdhits(PHTpcTrackerUtil::convert_clusters_to_hits(cluster_map, hitsets));
54  LOG_DEBUG("tracking.PHTpcSeedFinder.findSeeds") << "imported " << kdhits.size() << " clusters";
55 
56  //----- find unfitted kdtracks -----
57  std::vector<std::vector<double> > unused_hits;
58  std::vector<std::vector<std::vector<double> > > kdtracks;
59 
60  kdtracks = kdfinder::find_tracks_iterative<double>(kdhits, unused_hits,
61  mMaxDistance1 /* max distance in cm*/, mTripletAngle1 /* triplet angle */, mMinHits1 /* min hits to keep track */, // first iteration
62  mMaxDistance2, mTripletAngle2, mMinHits2, // second iteration params
63  mNThreads /* nthreads */,
64  false); // print stats
65 
66  LOG_DEBUG("tracking.PHTpcSeedFinder.findSeeds") << "kdtracks: " << kdtracks.size();
67 
68  //----- create fitted track candidates -----
69  result = kdfinder::get_track_candidates<double>(kdtracks, B);
70  LOG_DEBUG("tracking.PHTpcSeedFinder.findSeeds") << "initial kdtrack candidates: " << result.size();
71 
72  //----- filter out junk = bad tracks, <25 MeV Pt tracks etc.. -----
73  for (auto it = result.begin(); it != result.end();)
74  {
75  if (!(*it)->isFitted() || (*it)->Pt() < 0.025 || (*it)->Pt() > 200.0)
76  {
77  (*it)->deleteHits();
78  delete (*it);
79  it = result.erase(it);
80  }
81  else
82  {
83  ++it;
84  }
85  }
86  LOG_DEBUG("tracking.PHTpcSeedFinder.findSeeds") << "kdtrack candidates after filtering: " << result.size();
87 
88  //----- filter out loopers, if enabled -----
89  if (mRemoveLoopers)
90  {
91  for (auto it = result.begin(); it != result.end();)
92  {
93  kdfinder::Circle<double>* c = (*it)->getCircleFit();
94  double x = c->a, y = c->b, r = c->r, xyr = std::sqrt(x * x + y * y);
95  if (((xyr - r) > mMinLooperRadius) && ((xyr + r) < mMaxLooperRadius))
96  {
97  (*it)->deleteHits();
98  delete (*it);
99  it = result.erase(it);
100  }
101  else
102  {
103  ++it;
104  }
105  }
106  }
107 
108  /*
109  //----- primitive yet very fast vertex seed finding -----
110  // FIXME: use it instead of Rave?
111  std::vector< std::pair< std::vector<double>, std::vector<size_t> > > vertices =
112  kdfinder::find_vertex_seeds<double>( result, 0, 0, 0.5, 2.0, 3, false );
113  LOG_DEBUG("tracking.PHTpcSeedFinder.findSeeds") << "track candidate vertex seeds: " << vertices.size();
114  for( auto vertex: vertices ) {
115  std::vector<double> pos = vertex.first;
116  std::vector<size_t> ids = vertex.second;
117  LOG_DEBUG("tracking.PHTpcSeedFinder.findSeeds") << "vertex tracks: " << ids.size() << ", pos: " << pos[0] << ", " << pos[1] << ", " << pos[2];
118  }
119  */
120 
121  return result;
122 }