CbmRichRingFinderHough.cxx 4.19 KB
Newer Older
1
2
3
4
5
6
7
8
/**
* \file CbmRichRingFinderHough.cxx
*
* \author Semen Lebedev
* \date 2008
**/

#include "CbmRichRingFinderHough.h"
9

10
11
12
#include "CbmRichRingFinderHoughImpl.h"
//#include "CbmRichRingFinderHoughSimd.h"
//#include "../../littrack/utils/CbmLitMemoryManagment.h"
13
#include "CbmEvent.h"
14
15
16
#include "CbmRichHit.h"
#include "CbmRichRing.h"

17
#include <Logger.h>
Administrator's avatar
Administrator committed
18

Administrator's avatar
Administrator committed
19
#include "TClonesArray.h"
20
21
22
23
24
25
26
27
#include "TStopwatch.h"

#include <iostream>

using std::cout;
using std::endl;
using std::vector;

Administrator's avatar
Administrator committed
28
CbmRichRingFinderHough::CbmRichRingFinderHough()
29
{
30
#ifdef HOUGH_SERIAL
Administrator's avatar
Administrator committed
31
  fHTImpl = new CbmRichRingFinderHoughImpl();
32
33
34
#endif

#ifdef HOUGH_SIMD
Administrator's avatar
Administrator committed
35
  fHTImpl = new CbmRichRingFinderHoughSimd();
36
37
38
#endif
}

39
40
void CbmRichRingFinderHough::Init()
{
Administrator's avatar
Administrator committed
41
42
  fHTImpl->SetUseAnnSelect(fUseAnnSelect);
  fHTImpl->Init();
43
44
}

Administrator's avatar
Administrator committed
45
46
CbmRichRingFinderHough::~CbmRichRingFinderHough() { delete fHTImpl; }

47
48
49
Int_t CbmRichRingFinderHough::DoFind(CbmEvent* event, TClonesArray* rHitArray, TClonesArray* /*rProjArray*/,
                                     TClonesArray* rRingArray)
{
Administrator's avatar
Administrator committed
50
51
  TStopwatch timer;
  timer.Start();
52
  fEventNum++;
Administrator's avatar
Administrator committed
53
54
55
56

  vector<CbmRichHoughHit> UpH;
  vector<CbmRichHoughHit> DownH;

57
58
  if (rHitArray == nullptr) {
    LOG(error) << "CbmRichRingFinderHough::DoFind(): Hit array is nullptr.";
Administrator's avatar
Administrator committed
59
60
    return -1;
  }
61
62
63
64

  const Int_t nofRichHits = event ? event->GetNofData(ECbmDataType::kRichHit) : rHitArray->GetEntriesFast();
  if (nofRichHits <= 0) {
    LOG(debug) << "CbmRichRingFinderHough::DoFind(): No RICH hits in this event.";
Administrator's avatar
Administrator committed
65
66
67
    return -1;
  }

68
69
  UpH.reserve(nofRichHits / 2);
  DownH.reserve(nofRichHits / 2);
Administrator's avatar
Administrator committed
70
71
72

  // convert CbmRichHit to CbmRichHoughHit and
  // sort hits according to the photodetector (up or down)
73
74
75
76
  for (Int_t iH0 = 0; iH0 < nofRichHits; iH0++) {
    Int_t iH        = event ? event->GetIndex(ECbmDataType::kRichHit, iH0) : iH0;
    CbmRichHit* hit = static_cast<CbmRichHit*>(rHitArray->At(iH));
    if (hit != nullptr) {
Administrator's avatar
Administrator committed
77
      CbmRichHoughHit tempPoint;
78
79
80
81
82
83
84
      tempPoint.fHit.fX   = hit->GetX();
      tempPoint.fHit.fY   = hit->GetY();
      tempPoint.fHit.fId  = iH;
      tempPoint.fX2plusY2 = hit->GetX() * hit->GetX() + hit->GetY() * hit->GetY();
      tempPoint.fTime     = hit->GetTime();
      tempPoint.fIsUsed   = false;
      if (hit->GetY() >= 0) UpH.push_back(tempPoint);
Administrator's avatar
Administrator committed
85
86
87
88
89
90
91
92
93
94
95
96
      else
        DownH.push_back(tempPoint);
    }
  }

  timer.Stop();
  Double_t dt1 = timer.RealTime();

  timer.Start();

  fHTImpl->SetData(UpH);
  fHTImpl->DoFind();
97
  if (rRingArray != nullptr) AddRingsToOutputArray(event, rRingArray, rHitArray, fHTImpl->GetFoundRings());
Administrator's avatar
Administrator committed
98
99
100
101
102
103
104
105
106
  //for_each(UpH.begin(), UpH.end(), DeleteObject());
  UpH.clear();

  timer.Stop();
  Double_t dt2 = timer.RealTime();

  timer.Start();
  fHTImpl->SetData(DownH);
  fHTImpl->DoFind();
107
  if (rRingArray != nullptr) AddRingsToOutputArray(event, rRingArray, rHitArray, fHTImpl->GetFoundRings());
Administrator's avatar
Administrator committed
108
109
110
111
112
  //for_each(DownH.begin(), DownH.end(), DeleteObject());
  DownH.clear();
  timer.Stop();
  Double_t dt3 = timer.RealTime();

113
  int nofFoundRings = event ? event->GetNofData(ECbmDataType::kRichRing) : rRingArray->GetEntriesFast();
114
115
116
  LOG(debug) << "CbmRichRingFinderHough::DoFind(): Event:" << fEventNum << " hits:" << nofRichHits
             << " rings:" << nofFoundRings << " ringsInTS:" << rRingArray->GetEntriesFast()
             << " Time:" << dt1 + dt2 + dt3;
Administrator's avatar
Administrator committed
117
118

  return 1;
119
120
}

121
122
123
124
void CbmRichRingFinderHough::AddRingsToOutputArray(CbmEvent* event, TClonesArray* rRingArray, TClonesArray* rHitArray,
                                                   const vector<CbmRichRingLight*>& rings)
{

Administrator's avatar
Administrator committed
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
  for (UInt_t iRing = 0; iRing < rings.size(); iRing++) {
    if (rings[iRing]->GetRecFlag() == -1) continue;
    CbmRichRing* r    = new CbmRichRing();
    double ringTime   = 0.;
    Int_t ringCounter = 0;

    for (Int_t iH = 0; iH < rings[iRing]->GetNofHits(); iH++) {
      Int_t hitId = rings[iRing]->GetHitId(iH);
      r->AddHit(hitId);
      CbmRichHit* hit = static_cast<CbmRichHit*>(rHitArray->At(hitId));
      if (hit != nullptr) {
        ringCounter++;
        ringTime += hit->GetTime();
      }
    }
    r->SetTime(ringTime / (double) ringCounter);
141
142
143
    int nofRings = rRingArray->GetEntriesFast();
    new ((*rRingArray)[nofRings]) CbmRichRing(*r);
    if (event != nullptr) event->AddData(ECbmDataType::kRichRing, nofRings);
Administrator's avatar
Administrator committed
144
  }
145
}