CbmRoot
CbmRichRingFinderHoughImpl.cxx
Go to the documentation of this file.
1 
9 #include "CbmRichRingLight.h"
10 
11 #include "../../littrack/std/utils/CbmLitMemoryManagment.h"
12 #include "CbmRichRingFitterCOP.h"
13 #include "CbmRichRingSelectAnn.h"
14 
15 #include "TSystem.h"
16 
17 #include <algorithm>
18 #include <cmath>
19 #include <iostream>
20 #include <set>
21 
22 using std::cout;
23 using std::endl;
24 using std::set;
25 using std::sort;
26 using std::vector;
27 
29  : fNofParts(0)
30  ,
31 
32  fMaxDistance(0.f)
33  , fMinDistance(0.f)
34  , fMinDistanceSq(0.f)
35  , fMaxDistanceSq(0.f)
36  ,
37 
38  fMinRadius(0.f)
39  , fMaxRadius(0.f)
40  ,
41 
42  fDx(0.f)
43  , fDy(0.f)
44  , fDr(0.f)
45  , fNofBinsX(0)
46  , fNofBinsY(0)
47  , fNofBinsXY(0)
48  ,
49 
50  fHTCut(0)
51  ,
52 
53  fNofBinsR(0)
54  , fHTCutR(0)
55  ,
56 
57  fMinNofHitsInArea(0)
58  ,
59 
60  fRmsCoeffEl(0.f)
61  , fMaxCutEl(0.f)
62  , fRmsCoeffCOP(0.f)
63  , fMaxCutCOP(0.f)
64  ,
65 
66  fAnnCut(0.f)
67  , fUsedHitsAllCut(0.f)
68  ,
69 
70  fTimeCut(0.)
71  ,
72 
73  fCurMinX(0.f)
74  , fCurMinY(0.f)
75  ,
76 
77  fUseAnnSelect(true)
78  ,
79 
80  fData()
81  , fHist()
82  , fHistR()
83  , fHitInd()
84  , fFoundRings()
85  , fFitCOP(NULL)
86  , fANNSelect(NULL)
87  ,
88 
89  fCurTime(0.)
90 
91 {}
92 
94  if (NULL != fFitCOP) delete fFitCOP;
95  if (NULL != fANNSelect) delete fANNSelect;
96 }
97 
99  SetParameters();
100 
101  fHist.resize(fNofBinsXY);
102  fHistR.resize(fNofBinsR);
103  fHitInd.resize(fNofParts);
104 
106  if (fUseAnnSelect) {
108  fANNSelect->Init();
109  }
110 }
111 
113  // if (fData.size() > MAX_NOF_HITS) {
114  // cout << endl << endl << "-E- CbmRichRingFinderHoughImpl::DoFind" << ". Number of hits is more than " << MAX_NOF_HITS << endl << endl;
115  // return;
116  // }
117 
118  for_each(fFoundRings.begin(), fFoundRings.end(), DeleteObject());
119  fFoundRings.clear();
120  fFoundRings.reserve(100);
121 
122  std::sort(fData.begin(), fData.end(), CbmRichHoughHitCmpUp());
124  RingSelection();
125  fData.clear();
126 }
127 
129  fMaxDistance = 11.5;
130  fMinDistance = 2.5;
133 
134  fMinRadius = 3.5;
135  fMaxRadius = 5.7;
136 
137  fHTCut = 20;
138  fHTCutR = 12;
139  fMinNofHitsInArea = 4;
140 
141  fNofBinsX = 25;
142  fNofBinsY = 25;
143  fNofBinsR = 32;
144 
145  fAnnCut = 0.6;
146  fUsedHitsAllCut = 0.4;
147 
148  fRmsCoeffEl = 2.5;
149  fMaxCutEl = 1.0;
150  fRmsCoeffCOP = 3.;
151  fMaxCutCOP = 1.0;
152 
153  fTimeCut = 3.;
154 
155  fNofParts = 1;
156  fDx = 2.f * fMaxDistance / (float) fNofBinsX;
157  fDy = 2.f * fMaxDistance / (float) fNofBinsY;
158  fDr = fMaxRadius / (float) fNofBinsR;
160 }
161 
163  int indmin, indmax;
164  unsigned int size = fData.size();
165  for (unsigned int iHit = 0; iHit < size; iHit++) {
166  if (fData[iHit].fIsUsed == true) continue;
167 
168  fCurMinX = fData[iHit].fHit.fX - fMaxDistance;
169  fCurMinY = fData[iHit].fHit.fY - fMaxDistance;
170  fCurTime = fData[iHit].fTime;
171 
173  fData[iHit].fHit.fX, fData[iHit].fHit.fY, &indmin, &indmax);
174  HoughTransform(indmin, indmax);
175  FindPeak(indmin, indmax);
176  }
177 }
178 
180  float y0,
181  int* indmin,
182  int* indmax) {
183  CbmRichHoughHit mpnt;
184  vector<CbmRichHoughHit>::iterator itmin, itmax;
185 
186  //find all hits which are in the corridor
187  mpnt.fHit.fX = x0 - fMaxDistance;
188  itmin =
189  std::lower_bound(fData.begin(), fData.end(), mpnt, CbmRichHoughHitCmpUp());
190 
191  mpnt.fHit.fX = x0 + fMaxDistance;
192  itmax =
193  std::lower_bound(fData.begin(), fData.end(), mpnt, CbmRichHoughHitCmpUp())
194  - 1;
195 
196  *indmin = itmin - fData.begin();
197  *indmax = itmax - fData.begin();
198 
199  int arSize = *indmax - *indmin + 1;
200  if (arSize <= fMinNofHitsInArea) return;
201 
202  for (unsigned short i = 0; i < fNofParts; i++) {
203  fHitInd[i].clear();
204  fHitInd[i].reserve((*indmax - *indmin) / fNofParts);
205  }
206 
207  for (int i = *indmin; i <= *indmax; i++) {
208  if (fData[i].fIsUsed == true) continue;
209  float ry = y0 - fData[i].fHit.fY;
210  if (fabs(ry) > fMaxDistance) continue;
211  float rx = x0 - fData[i].fHit.fX;
212  float d = rx * rx + ry * ry;
213  if (d > fMaxDistanceSq) continue;
214  if (std::fabs(fCurTime - fData[i].fTime) > fTimeCut) continue;
215 
216  fHitInd[i % fNofParts].push_back(i);
217  }
218 
219  for (unsigned short j = 0; j < fNofBinsXY; j++) {
220  fHist[j] = 0;
221  }
222 
223  for (unsigned short k = 0; k < fNofBinsR; k++) {
224  fHistR[k] = 0;
225  }
226 }
227 
229  unsigned int indmax) {
230  for (int iPart = 0; iPart < fNofParts; iPart++) {
231  HoughTransformGroup(indmin, indmax, iPart);
232  } //iPart
233 }
234 
236  unsigned int /*indmax*/,
237  int iPart) {
238  unsigned int nofHits = fHitInd[iPart].size();
239  float xcs, ycs; // xcs = xc - fCurMinX
240  float dx = 1.0f / fDx, dy = 1.0f / fDy, dr = 1.0f / fDr;
241 
242  vector<CbmRichHoughHit> data;
243  data.resize(nofHits);
244  for (unsigned int i = 0; i < nofHits; i++) {
245  data[i] = fData[fHitInd[iPart][i]];
246  }
247 
248  typedef vector<CbmRichHoughHit>::iterator iH;
249 
250  for (iH iHit1 = data.begin(); iHit1 != data.end(); iHit1++) {
251  float iH1X = iHit1->fHit.fX;
252  float iH1Y = iHit1->fHit.fY;
253  double time1 = iHit1->fTime;
254 
255  for (iH iHit2 = iHit1 + 1; iHit2 != data.end(); iHit2++) {
256  float iH2X = iHit2->fHit.fX;
257  float iH2Y = iHit2->fHit.fY;
258  double time2 = iHit2->fTime;
259  if (std::fabs(time1 - time2) > fTimeCut) continue;
260 
261  float rx0 = iH1X - iH2X; //rx12
262  float ry0 = iH1Y - iH2Y; //ry12
263  float r12 = rx0 * rx0 + ry0 * ry0;
264 
265  if (r12 < fMinDistanceSq || r12 > fMaxDistanceSq) continue;
266 
267  float t10 = iHit1->fX2plusY2 - iHit2->fX2plusY2;
268  for (iH iHit3 = iHit2 + 1; iHit3 != data.end(); iHit3++) {
269  float iH3X = iHit3->fHit.fX;
270  float iH3Y = iHit3->fHit.fY;
271  double time3 = iHit3->fTime;
272 
273  if (std::fabs(time1 - time3) > fTimeCut) continue;
274  if (std::fabs(time2 - time3) > fTimeCut) continue;
275 
276  float rx1 = iH1X - iH3X; //rx13
277  float ry1 = iH1Y - iH3Y; //ry13
278  float r13 = rx1 * rx1 + ry1 * ry1;
279  if (r13 < fMinDistanceSq || r13 > fMaxDistanceSq) continue;
280 
281  float rx2 = iH2X - iH3X; //rx23
282  float ry2 = iH2Y - iH3Y; //ry23
283  float r23 = rx2 * rx2 + ry2 * ry2;
284  if (r23 < fMinDistanceSq || r23 > fMaxDistanceSq) continue;
285 
286  float det = rx2 * ry0 - rx0 * ry2;
287  if (det == 0.0f) continue;
288  float t19 = 0.5f / det;
289  float t5 = iHit2->fX2plusY2 - iHit3->fX2plusY2;
290 
291  float xc = (t5 * ry0 - t10 * ry2) * t19;
292  xcs = xc - fCurMinX;
293  int intX = int(xcs * dx);
294  if (intX < 0 || intX >= fNofBinsX) continue;
295 
296  float yc = (t10 * rx2 - t5 * rx0) * t19;
297  ycs = yc - fCurMinY;
298  int intY = int(ycs * dy);
299  if (intY < 0 || intY >= fNofBinsY) continue;
300 
301  //radius calculation
302  float t6 = iH1X - xc;
303  float t7 = iH1Y - yc;
304  //if (t6 > fMaxRadius || t7 > fMaxRadius) continue;
305  float r = sqrt(t6 * t6 + t7 * t7);
306  //if (r < fMinRadius) continue;
307  int intR = int(r * dr);
308  if (intR < 0 || intR >= fNofBinsR) continue;
309 
310  fHist[intX * fNofBinsX + intY]++;
311  fHistR[intR]++;
312  } // iHit1
313  } // iHit2
314  } // iHit3
315  //hitIndPart.clear();
316 }
317 
318 void CbmRichRingFinderHoughImpl::FindPeak(int indmin, int indmax) {
319  // Find MAX bin R and compare it with CUT
320  int maxBinR = -1, maxR = -1;
321  unsigned int size = fHistR.size();
322  for (unsigned int k = 0; k < size; k++) {
323  if (fHistR[k] > maxBinR) {
324  maxBinR = fHistR[k];
325  maxR = k;
326  }
327  }
328  //cout << "maxBinR:" << maxBinR << endl;
329  if (maxBinR < fHTCutR) return;
330 
331  // Find MAX bin XY and compare it with CUT
332  int maxBinXY = -1, maxXY = -1;
333  size = fHist.size();
334  for (unsigned int k = 0; k < size; k++) {
335  if (fHist[k] > maxBinXY) {
336  maxBinXY = fHist[k];
337  maxXY = k;
338  }
339  }
340  //cout << "maxBinXY:" << maxBinXY << endl;
341  if (maxBinXY < fHTCut) return;
342 
343  CbmRichRingLight* ring1 = new CbmRichRingLight();
344 
345  // Find Preliminary Xc, Yc, R
346  float xc, yc, r;
347  float rx, ry, dr;
348  xc = (maxXY / fNofBinsX + 0.5f) * fDx + fCurMinX;
349  yc = (maxXY % fNofBinsX + 0.5f) * fDy + fCurMinY;
350  r = (maxR + 0.5f) * fDr;
351  for (int j = indmin; j < indmax + 1; j++) {
352  if (std::fabs(fCurTime - fData[j].fTime) > fTimeCut) continue;
353  rx = fData[j].fHit.fX - xc;
354  ry = fData[j].fHit.fY - yc;
355 
356  dr = fabs(sqrt(rx * rx + ry * ry) - r);
357  if (dr > 0.9f) continue;
358  ring1->AddHit(fData[j].fHit);
359  }
360 
361  if (ring1->GetNofHits() < 7) {
362  delete ring1;
363  return;
364  }
365 
366  fFitCOP->DoFit(ring1);
367  float drCOPCut = fRmsCoeffCOP * sqrt(ring1->GetChi2() / ring1->GetNofHits());
368  if (drCOPCut > fMaxCutCOP) drCOPCut = fMaxCutCOP;
369 
370  xc = ring1->GetCenterX();
371  yc = ring1->GetCenterY();
372  r = ring1->GetRadius();
373 
374  delete ring1;
375 
376  CbmRichRingLight* ring2 = new CbmRichRingLight();
377  for (int j = indmin; j < indmax + 1; j++) {
378  if (std::fabs(fCurTime - fData[j].fTime) > fTimeCut) continue;
379 
380  rx = fData[j].fHit.fX - xc;
381  ry = fData[j].fHit.fY - yc;
382 
383 
384  dr = fabs(sqrt(rx * rx + ry * ry) - r);
385  if (dr > drCOPCut) continue;
386  //fData[j+indmin].fIsUsed = true;
387  ring2->AddHit(fData[j].fHit);
388  }
389 
390  if (ring2->GetNofHits() < 7) {
391  delete ring2;
392  return;
393  }
394 
395  fFitCOP->DoFit(ring2);
396 
397  if (fUseAnnSelect) fANNSelect->DoSelect(ring2);
398  float select = (fUseAnnSelect) ? ring2->GetSelectionNN() : 1.;
399  // cout << select << endl;
400 
401  // Remove found hits only for good quality rings
402  if (select > fAnnCut) { RemoveHitsAroundRing(indmin, indmax, ring2); }
403 
404  if (select > -0.7) {
405  fFoundRings.push_back(ring2);
406  ring2 = NULL;
407  }
408  delete ring2;
409 }
410 
412  int indmax,
413  CbmRichRingLight* ring) {
414  float rms = sqrt(ring->GetChi2() / ring->GetNofHits());
415  float dCut = fRmsCoeffEl * rms;
416  if (dCut > fMaxCutEl) dCut = fMaxCutEl;
417 
418  for (int j = indmin; j < indmax + 1; j++) {
419  if (std::fabs(fCurTime - fData[j].fTime) > fTimeCut) continue;
420  float rx = fData[j].fHit.fX - ring->GetCenterX();
421  float ry = fData[j].fHit.fY - ring->GetCenterY();
422 
423  float dr = fabs(sqrt(rx * rx + ry * ry) - ring->GetRadius());
424  if (dr < dCut) { fData[j].fIsUsed = true; }
425  }
426 }
427 
429  int nofRings = fFoundRings.size();
430  sort(fFoundRings.begin(), fFoundRings.end(), CbmRichRingComparatorMore());
431  set<unsigned int> usedHitsAll;
432  vector<unsigned int> goodRingIndex;
433  goodRingIndex.reserve(nofRings);
434  // CbmRichRingLight* ring2;
435 
436  for (int iRing = 0; iRing < nofRings; iRing++) {
437  CbmRichRingLight* ring = fFoundRings[iRing];
438  ring->SetRecFlag(-1);
439  int nofHits = ring->GetNofHits();
440  bool isGoodRingAll = true;
441  int nofUsedHitsAll = 0;
442  for (int iHit = 0; iHit < nofHits; iHit++) {
443  set<unsigned int>::iterator it = usedHitsAll.find(ring->GetHitId(iHit));
444  if (it != usedHitsAll.end()) { nofUsedHitsAll++; }
445  }
446  if ((float) nofUsedHitsAll / (float) nofHits > fUsedHitsAllCut) {
447  isGoodRingAll = false;
448  }
449 
450  if (isGoodRingAll) {
451  fFoundRings[iRing]->SetRecFlag(1);
452  for (unsigned int iRSet = 0; iRSet < goodRingIndex.size(); iRSet++) {
453  ReAssignSharedHits(goodRingIndex[iRSet], iRing);
454  }
455  goodRingIndex.push_back(iRing);
456  for (int iHit = 0; iHit < nofHits; iHit++) {
457  usedHitsAll.insert(ring->GetHitId(iHit));
458  }
459  } // isGoodRing
460  } // iRing
461 
462  // usedHits.clear();
463  usedHitsAll.clear();
464  goodRingIndex.clear();
465 }
466 
468  int ringInd2) {
469  CbmRichRingLight* ring1 = fFoundRings[ringInd1];
470  CbmRichRingLight* ring2 = fFoundRings[ringInd2];
471  int nofHits1 = ring1->GetNofHits();
472  int nofHits2 = ring2->GetNofHits();
473 
474  for (int iHit1 = 0; iHit1 < nofHits1; iHit1++) {
475  unsigned int hitId1 = ring1->GetHitId(iHit1);
476  for (int iHit2 = 0; iHit2 < nofHits2; iHit2++) {
477  unsigned int hitId2 = ring2->GetHitId(iHit2);
478  if (hitId1 != hitId2) continue;
479  int hitIndData = GetHitIndexById(hitId1);
480  float hitX = fData[hitIndData].fHit.fX;
481  float hitY = fData[hitIndData].fHit.fY;
482  float rx1 = hitX - ring1->GetCenterX();
483  float ry1 = hitY - ring1->GetCenterY();
484  float dr1 = fabs(sqrt(rx1 * rx1 + ry1 * ry1) - ring1->GetRadius());
485 
486  float rx2 = hitX - ring2->GetCenterX();
487  float ry2 = hitY - ring2->GetCenterY();
488  float dr2 = fabs(sqrt(rx2 * rx2 + ry2 * ry2) - ring2->GetRadius());
489 
490  if (dr1 > dr2) {
491  ring1->RemoveHit(hitId1);
492  } else {
493  ring2->RemoveHit(hitId2);
494  }
495  } //iHit2
496  } //iHit1
497 }
498 
500  unsigned int size = fData.size();
501  for (unsigned int i = 0; i < size; i++) {
502  if (fData[i].fHit.fId == hitId) return i;
503  }
504  return -1;
505 }
506 
508  float y[],
509  float* xc,
510  float* yc,
511  float* r) {
512  float t1, t2, t3, t4, t5, t6, t8, t9, t10, t11, t14, t16, t19, t21, t41;
513 
514  t1 = x[1] * x[1];
515  t2 = x[2] * x[2];
516  t3 = y[1] * y[1];
517  t4 = y[2] * y[2];
518  t5 = t1 - t2 + t3 - t4;
519  t6 = y[0] - y[1];
520  t8 = x[0] * x[0];
521  t9 = y[0] * y[0];
522  t10 = t8 - t1 + t9 - t3;
523  t11 = y[1] - y[2];
524  t14 = x[1] - x[2];
525  t16 = x[0] - x[1];
526  t19 = 1.0f / (t14 * t6 - t16 * t11);
527 
528  *xc = 0.5e0 * (t5 * t6 - t10 * t11) * t19;
529  *yc = 0.5e0 * (t10 * t14 - t5 * t16) * t19;
530 
531  t21 = (x[0] - *xc) * (x[0] - *xc);
532  t41 = (y[0] - *yc) * (y[0] - *yc);
533  *r = sqrt(t21 + t41);
534 }
CbmRichRingFinderHoughImpl::GetHitIndexById
int GetHitIndexById(unsigned int hitId)
Return hit indez in the internal Array.
Definition: CbmRichRingFinderHoughImpl.cxx:499
CbmRichRingFinderHoughImpl::fHitInd
vector< vector< unsigned int > > fHitInd
Definition: CbmRichRingFinderHoughImpl.h:85
CbmRichRingFinderHoughImpl.h
Ring finder implementation based on Hough Transform method.
f
float f
Definition: L1/vectors/P4_F32vec4.h:24
CbmRichRingSelectAnn
Implementation for concrete RICH ring selection algorithm: reject rings using a trained neural net (i...
Definition: CbmRichRingSelectAnn.h:39
CbmRichRingFinderHoughImpl::fCurMinY
float fCurMinY
Definition: CbmRichRingFinderHoughImpl.h:77
CbmRichHoughHitCmpUp
CbmRichHoughHit comparator for hits sorting by X coordinate.
Definition: CbmRichRingFinderData.h:46
CbmRichRingFinderHoughImpl::fMinDistance
float fMinDistance
Definition: CbmRichRingFinderHoughImpl.h:43
sqrt
friend F32vec4 sqrt(const F32vec4 &a)
Definition: L1/vectors/P4_F32vec4.h:41
CbmRichRingFinderHoughImpl::fNofParts
unsigned short fNofParts
Definition: CbmRichRingFinderHoughImpl.h:40
CbmRichRingFinderHoughImpl::HoughTransformReconstruction
virtual void HoughTransformReconstruction()
Run HT for each hit.
Definition: CbmRichRingFinderHoughImpl.cxx:162
CbmRichRingSelectAnn::Init
virtual void Init()
Initialize ANN.
Definition: CbmRichRingSelectAnn.cxx:32
CbmRichRingFinderHoughImpl::HoughTransform
virtual void HoughTransform(unsigned int indmin, unsigned int indmax)
Run HoughTransformGroup for each group of hits.
Definition: CbmRichRingFinderHoughImpl.cxx:228
CbmRichRingLight::GetChi2
float GetChi2() const
Definition: CbmRichRingLight.h:242
CbmRichRingFinderHoughImpl::fNofBinsY
unsigned short fNofBinsY
Definition: CbmRichRingFinderHoughImpl.h:54
CbmRichRingFinderHoughImpl::fCurTime
double fCurTime
Definition: CbmRichRingFinderHoughImpl.h:90
CbmRichHoughHit
Implementation of RICH hit for ring finder algorithm.
Definition: CbmRichRingFinderData.h:21
CbmRichRingFinderHoughImpl::Init
void Init()
Definition: CbmRichRingFinderHoughImpl.cxx:98
CbmRichRingFinderHoughImpl::fANNSelect
CbmRichRingSelectAnn * fANNSelect
Definition: CbmRichRingFinderHoughImpl.h:88
i
int i
Definition: L1/vectors/P4_F32vec4.h:25
CbmRichRingFinderHoughImpl::RemoveHitsAroundRing
void RemoveHitsAroundRing(int indmin, int indmax, CbmRichRingLight *ring)
Set fIsUsed flag to true for hits attached to the ring.
Definition: CbmRichRingFinderHoughImpl.cxx:411
CbmRichRingFinderHoughImpl::fMinRadius
float fMinRadius
Definition: CbmRichRingFinderHoughImpl.h:47
CbmRichRingFinderHoughImpl::ReAssignSharedHits
void ReAssignSharedHits(int ringInd1, int ringInd2)
Reassign shared hits from two rings to only one of the rings.
Definition: CbmRichRingFinderHoughImpl.cxx:467
CbmRichRingFinderHoughImpl::fDy
float fDy
Definition: CbmRichRingFinderHoughImpl.h:51
CbmRichRingFinderHoughImpl::SetParameters
void SetParameters()
Set parameters of the algorithm.
Definition: CbmRichRingFinderHoughImpl.cxx:128
CbmRichRingFinderHoughImpl::fDr
float fDr
Definition: CbmRichRingFinderHoughImpl.h:52
CbmRichHitLight::fX
float fX
Definition: CbmRichRingLight.h:34
CbmRichRingLight::RemoveHit
bool RemoveHit(int hitId)
Remove hit from the ring.
Definition: CbmRichRingLight.h:94
CbmRichRingFinderHoughImpl::fAnnCut
float fAnnCut
Definition: CbmRichRingFinderHoughImpl.h:71
CbmRichRingFinderHoughImpl::CbmRichRingFinderHoughImpl
CbmRichRingFinderHoughImpl()
Standard constructor.
Definition: CbmRichRingFinderHoughImpl.cxx:28
CbmRichRingLight::GetHitId
unsigned int GetHitId(int ind)
Return hit index in TClonesArray.
Definition: CbmRichRingLight.h:120
DeleteObject
Definition: CbmLitMemoryManagment.h:5
CbmRichRingSelectAnn::DoSelect
void DoSelect(CbmRichRingLight *ring)
Definition: CbmRichRingSelectAnn.cxx:53
CbmRichRingFinderHoughImpl::fRmsCoeffCOP
float fRmsCoeffCOP
Definition: CbmRichRingFinderHoughImpl.h:68
CbmRichRingFitterCOP
Here the ring is fitted with the COP algorithm from A. Ayriyan/G. Ososkov.
Definition: CbmRichRingFitterCOP.h:28
CbmRichRingFinderHoughImpl::fMaxDistance
float fMaxDistance
Definition: CbmRichRingFinderHoughImpl.h:42
CbmRichRingFinderHoughImpl::fMinNofHitsInArea
unsigned short fMinNofHitsInArea
Definition: CbmRichRingFinderHoughImpl.h:64
CbmRichRingLight::GetNofHits
int GetNofHits() const
Return number of hits in ring.
Definition: CbmRichRingLight.h:108
CbmRichRingSelectAnn.h
Implementation for concrete RICH ring selection algorithm: reject rings using a trained neural net (i...
CbmRichRingFinderHoughImpl::fNofBinsXY
unsigned short fNofBinsXY
Definition: CbmRichRingFinderHoughImpl.h:55
CbmRichRingFinderHoughImpl::DefineLocalAreaAndHits
virtual void DefineLocalAreaAndHits(float x0, float y0, int *indmin, int *indmax)
Find hits in a local area.
Definition: CbmRichRingFinderHoughImpl.cxx:179
CbmRichRingFinderHoughImpl::HoughTransformGroup
virtual void HoughTransformGroup(unsigned int indmin, unsigned int indmax, int iPart)
Definition: CbmRichRingFinderHoughImpl.cxx:235
CbmRichRingFinderHoughImpl::fUseAnnSelect
bool fUseAnnSelect
Definition: CbmRichRingFinderHoughImpl.h:79
d
double d
Definition: P4_F64vec2.h:24
CbmRichRingFinderHoughImpl::~CbmRichRingFinderHoughImpl
virtual ~CbmRichRingFinderHoughImpl()
Distructor.
Definition: CbmRichRingFinderHoughImpl.cxx:93
CbmRichRingFinderHoughImpl::fHTCut
unsigned short fHTCut
Definition: CbmRichRingFinderHoughImpl.h:58
CbmRichRingFinderHoughImpl::fHist
vector< unsigned short > fHist
Definition: CbmRichRingFinderHoughImpl.h:82
CbmRichRingFinderHoughImpl::FindPeak
void FindPeak(int indmin, int indmax)
Definition: CbmRichRingFinderHoughImpl.cxx:318
CbmRichRingLight::GetCenterY
float GetCenterY() const
Definition: CbmRichRingLight.h:160
CbmRichRingLight.h
CbmRichRingFinderHoughImpl::fMaxRadius
float fMaxRadius
Definition: CbmRichRingFinderHoughImpl.h:48
CbmRichRingFinderHoughImpl::fMinDistanceSq
float fMinDistanceSq
Definition: CbmRichRingFinderHoughImpl.h:44
CbmRichRingFinderHoughImpl::fHistR
vector< unsigned short > fHistR
Definition: CbmRichRingFinderHoughImpl.h:83
CbmRichRingLight::AddHit
void AddHit(CbmRichHitLight hit)
Add new hit to the ring.
Definition: CbmRichRingLight.h:87
CbmRichRingFinderHoughImpl::fRmsCoeffEl
float fRmsCoeffEl
Definition: CbmRichRingFinderHoughImpl.h:66
CbmRichRingLight::SetRecFlag
void SetRecFlag(int r)
Definition: CbmRichRingLight.h:263
CbmRichRingFinderHoughImpl::fUsedHitsAllCut
float fUsedHitsAllCut
Definition: CbmRichRingFinderHoughImpl.h:72
CbmRichRingFinderHoughImpl::fTimeCut
double fTimeCut
Definition: CbmRichRingFinderHoughImpl.h:74
CbmRichRingFinderHoughImpl::fFoundRings
vector< CbmRichRingLight * > fFoundRings
Definition: CbmRichRingFinderHoughImpl.h:86
CbmRichRingFinderHoughImpl::DoFind
void DoFind()
Start point to run algorithm.
Definition: CbmRichRingFinderHoughImpl.cxx:112
CbmRichRingFinderHoughImpl::CalculateRingParameters
void CalculateRingParameters(float x[], float y[], float *xc, float *yc, float *r)
Calculate circle center and radius.
Definition: CbmRichRingFinderHoughImpl.cxx:507
CbmRichRingFinderHoughImpl::fNofBinsR
unsigned short fNofBinsR
Definition: CbmRichRingFinderHoughImpl.h:60
CbmRichRingFinderHoughImpl::fCurMinX
float fCurMinX
Definition: CbmRichRingFinderHoughImpl.h:76
CbmRichRingFitterCOP.h
Here the ring is fitted with the COP algorithm from A. Ayriyan/G. Ososkov.
CbmRichRingFinderHoughImpl::fDx
float fDx
Definition: CbmRichRingFinderHoughImpl.h:50
CbmRichRingLight::GetRadius
float GetRadius() const
Definition: CbmRichRingLight.h:161
CbmRichRingFinderHoughImpl::fFitCOP
CbmRichRingFitterCOP * fFitCOP
Definition: CbmRichRingFinderHoughImpl.h:87
CbmRichRingComparatorMore
CbmRichRingLight comparator based on the selection ANN criterion.
Definition: CbmRichRingFinderData.h:66
CbmRichRingLight::GetSelectionNN
float GetSelectionNN() const
Definition: CbmRichRingLight.h:241
CbmRichRingFitterCOP::DoFit
virtual void DoFit(CbmRichRingLight *ring)
Inherited from CbmRichRingFitterBase.
Definition: CbmRichRingFitterCOP.cxx:19
x
Double_t x
Definition: CbmMvdSensorDigiToHitTask.cxx:68
fabs
friend F32vec4 fabs(const F32vec4 &a)
Definition: L1/vectors/P4_F32vec4.h:60
CbmRichRingFinderHoughImpl::fHTCutR
unsigned short fHTCutR
Definition: CbmRichRingFinderHoughImpl.h:62
CbmRichRingFinderHoughImpl::fMaxCutCOP
float fMaxCutCOP
Definition: CbmRichRingFinderHoughImpl.h:69
CbmRichRingFinderHoughImpl::RingSelection
void RingSelection()
Ring selection procedure.
Definition: CbmRichRingFinderHoughImpl.cxx:428
y
Double_t y
Definition: CbmMvdSensorDigiToHitTask.cxx:68
CbmRichRingFinderHoughImpl::fNofBinsX
unsigned short fNofBinsX
Definition: CbmRichRingFinderHoughImpl.h:53
CbmRichRingLight::GetCenterX
float GetCenterX() const
Definition: CbmRichRingLight.h:159
CbmRichRingFinderHoughImpl::fMaxCutEl
float fMaxCutEl
Definition: CbmRichRingFinderHoughImpl.h:67
CbmRichRingFinderHoughImpl::fData
vector< CbmRichHoughHit > fData
Definition: CbmRichRingFinderHoughImpl.h:81
CbmRichRingLight
Definition: CbmRichRingLight.h:39
CbmRichRingFinderHoughImpl::fMaxDistanceSq
float fMaxDistanceSq
Definition: CbmRichRingFinderHoughImpl.h:45
CbmRichHoughHit::fHit
CbmRichHitLight fHit
Definition: CbmRichRingFinderData.h:30