]> git.uio.no Git - u/mrichter/AliRoot.git/blame - HLT/hough/AliHLTHoughTransformerRow.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
compilation warnings corrected
[u/mrichter/AliRoot.git] / HLT / hough / AliHLTHoughTransformerRow.cxx
CommitLineData
0bd0c1ef 1// @(#) $Id$
2
de3c3890 3
ce4c28d5 4// Author: Cvetan Cheshkov <mailto:cvetan.cheshkov@cern.ch>
0bd0c1ef 5
4aa41877 6#include "AliHLTStandardIncludes.h"
0bd0c1ef 7
4aa41877 8#include "AliHLTLogging.h"
9#include "AliHLTMemHandler.h"
10#include "AliHLTTransform.h"
11#include "AliHLTDigitData.h"
12#include "AliHLTHistogramAdaptive.h"
13#include "AliHLTHoughTrack.h"
14#include "AliHLTHoughTransformerRow.h"
0bd0c1ef 15
c6747af6 16#if __GNUC__ >= 3
0bd0c1ef 17using namespace std;
18#endif
19
20//_____________________________________________________________
4aa41877 21// AliHLTHoughTransformerRow
0bd0c1ef 22//
ce4c28d5 23// Impelementation of the so called "TPC rows Hough transformation" class
0bd0c1ef 24//
ce4c28d5 25// Transforms the TPC data into the hough space and counts the missed TPC
26// rows corresponding to each track cnadidate - hough space bin
0bd0c1ef 27
4aa41877 28ClassImp(AliHLTHoughTransformerRow)
0bd0c1ef 29
4aa41877 30Float_t AliHLTHoughTransformerRow::fgBeta1 = 1.0/AliHLTTransform::Row2X(84);
31Float_t AliHLTHoughTransformerRow::fgBeta2 = 1.0/(AliHLTTransform::Row2X(158)*(1.0+tan(AliHLTTransform::Pi()*10/180)*tan(AliHLTTransform::Pi()*10/180)));
32Float_t AliHLTHoughTransformerRow::fgDAlpha = 0.22;
33Float_t AliHLTHoughTransformerRow::fgDEta = 0.40;
34Double_t AliHLTHoughTransformerRow::fgEtaCalcParam1 = 1.0289;
35Double_t AliHLTHoughTransformerRow::fgEtaCalcParam2 = 0.15192;
36Double_t AliHLTHoughTransformerRow::fgEtaCalcParam3 = 1./(32.*600.*600.);
0bd0c1ef 37
4aa41877 38AliHLTHoughTransformerRow::AliHLTHoughTransformerRow()
0bd0c1ef 39{
40 //Default constructor
41 fParamSpace = 0;
0bd0c1ef 42
a8ffd46b 43 fGapCount = 0;
44 fCurrentRowCount = 0;
45#ifdef do_mc
46 fTrackID = 0;
47#endif
48 fTrackNRows = 0;
49 fTrackFirstRow = 0;
50 fTrackLastRow = 0;
51 fInitialGapCount = 0;
52
53 fPrevBin = 0;
54 fNextBin = 0;
55 fNextRow = 0;
56
57 fStartPadParams = 0;
58 fEndPadParams = 0;
f644512a 59 fLUTr = 0;
a8ffd46b 60 fLUTforwardZ = 0;
a8ffd46b 61 fLUTbackwardZ = 0;
a8ffd46b 62
0bd0c1ef 63}
64
4aa41877 65AliHLTHoughTransformerRow::AliHLTHoughTransformerRow(Int_t slice,Int_t patch,Int_t netasegments,Bool_t /*DoMC*/,Float_t zvertex) : AliHLTHoughBaseTransformer(slice,patch,netasegments,zvertex)
0bd0c1ef 66{
67 //Normal constructor
68 fParamSpace = 0;
a8ffd46b 69
70 fGapCount = 0;
71 fCurrentRowCount = 0;
0bd0c1ef 72#ifdef do_mc
a8ffd46b 73 fTrackID = 0;
0bd0c1ef 74#endif
75
a8ffd46b 76 fTrackNRows = 0;
77 fTrackFirstRow = 0;
78 fTrackLastRow = 0;
79 fInitialGapCount = 0;
80
81 fPrevBin = 0;
82 fNextBin = 0;
83 fNextRow = 0;
84
85 fStartPadParams = 0;
86 fEndPadParams = 0;
f644512a 87 fLUTr = 0;
a8ffd46b 88 fLUTforwardZ = 0;
a8ffd46b 89 fLUTbackwardZ = 0;
a8ffd46b 90
0bd0c1ef 91}
92
4aa41877 93AliHLTHoughTransformerRow::~AliHLTHoughTransformerRow()
0bd0c1ef 94{
917e711b 95 //Destructor
a8ffd46b 96 if(fLastTransformer) return;
7e4eb453 97 DeleteHistograms();
0bd0c1ef 98#ifdef do_mc
a8ffd46b 99 if(fTrackID)
0bd0c1ef 100 {
101 for(Int_t i=0; i<GetNEtaSegments(); i++)
102 {
a8ffd46b 103 if(!fTrackID[i]) continue;
104 delete fTrackID[i];
0bd0c1ef 105 }
a8ffd46b 106 delete [] fTrackID;
107 fTrackID = 0;
0bd0c1ef 108 }
109#endif
110
a8ffd46b 111 if(fGapCount)
0bd0c1ef 112 {
113 for(Int_t i=0; i<GetNEtaSegments(); i++)
114 {
a8ffd46b 115 if(!fGapCount[i]) continue;
116 delete [] fGapCount[i];
0bd0c1ef 117 }
a8ffd46b 118 delete [] fGapCount;
119 fGapCount = 0;
0bd0c1ef 120 }
a8ffd46b 121 if(fCurrentRowCount)
0bd0c1ef 122 {
123 for(Int_t i=0; i<GetNEtaSegments(); i++)
124 {
a8ffd46b 125 if(fCurrentRowCount[i])
126 delete [] fCurrentRowCount[i];
0bd0c1ef 127 }
a8ffd46b 128 delete [] fCurrentRowCount;
129 fCurrentRowCount = 0;
130 }
131 if(fTrackNRows)
132 {
133 delete [] fTrackNRows;
134 fTrackNRows = 0;
135 }
136 if(fTrackFirstRow)
137 {
138 delete [] fTrackFirstRow;
139 fTrackFirstRow = 0;
140 }
141 if(fTrackLastRow)
142 {
143 delete [] fTrackLastRow;
144 fTrackLastRow = 0;
145 }
146 if(fInitialGapCount)
147 {
148 delete [] fInitialGapCount;
149 fInitialGapCount = 0;
0bd0c1ef 150 }
a8ffd46b 151 if(fPrevBin)
0bd0c1ef 152 {
153 for(Int_t i=0; i<GetNEtaSegments(); i++)
154 {
a8ffd46b 155 if(!fPrevBin[i]) continue;
156 delete [] fPrevBin[i];
157 }
158 delete [] fPrevBin;
159 fPrevBin = 0;
160 }
161 if(fNextBin)
162 {
163 for(Int_t i=0; i<GetNEtaSegments(); i++)
164 {
165 if(!fNextBin[i]) continue;
166 delete [] fNextBin[i];
167 }
168 delete [] fNextBin;
169 fNextBin = 0;
170 }
171 if(fNextRow)
172 {
173 for(Int_t i=0; i<GetNEtaSegments(); i++)
174 {
175 if(!fNextRow[i]) continue;
176 delete [] fNextRow[i];
177 }
178 delete [] fNextRow;
179 fNextRow = 0;
180 }
181 if(fStartPadParams)
182 {
4aa41877 183 for(Int_t i = AliHLTTransform::GetFirstRow(0); i<=AliHLTTransform::GetLastRow(5); i++)
a8ffd46b 184 {
185 if(!fStartPadParams[i]) continue;
186 delete [] fStartPadParams[i];
187 }
188 delete [] fStartPadParams;
189 fStartPadParams = 0;
190 }
191 if(fEndPadParams)
192 {
4aa41877 193 for(Int_t i = AliHLTTransform::GetFirstRow(0); i<=AliHLTTransform::GetLastRow(5); i++)
a8ffd46b 194 {
195 if(!fEndPadParams[i]) continue;
196 delete [] fEndPadParams[i];
197 }
198 delete [] fEndPadParams;
199 fEndPadParams = 0;
200 }
f644512a 201 if(fLUTr)
a8ffd46b 202 {
4aa41877 203 for(Int_t i = AliHLTTransform::GetFirstRow(0); i<=AliHLTTransform::GetLastRow(5); i++)
a8ffd46b 204 {
f644512a 205 if(!fLUTr[i]) continue;
206 delete [] fLUTr[i];
0bd0c1ef 207 }
f644512a 208 delete [] fLUTr;
209 fLUTr = 0;
0bd0c1ef 210 }
a8ffd46b 211 if(fLUTforwardZ)
de3c3890 212 {
a8ffd46b 213 delete[] fLUTforwardZ;
214 fLUTforwardZ=0;
de3c3890 215 }
a8ffd46b 216 if(fLUTbackwardZ)
de3c3890 217 {
a8ffd46b 218 delete[] fLUTbackwardZ;
219 fLUTbackwardZ=0;
220 }
0bd0c1ef 221}
222
4aa41877 223void AliHLTHoughTransformerRow::DeleteHistograms()
0bd0c1ef 224{
917e711b 225 // Clean up
0bd0c1ef 226 if(!fParamSpace)
227 return;
228 for(Int_t i=0; i<GetNEtaSegments(); i++)
229 {
230 if(!fParamSpace[i]) continue;
231 delete fParamSpace[i];
232 }
233 delete [] fParamSpace;
234}
235
4aa41877 236struct AliHLTTrackLength {
ce4c28d5 237 // Structure is used for temporarely storage of the LUT
238 // which contains the track lengths associated to each hough
239 // space bin
240 Bool_t fIsFilled; // Is bin already filled?
241 UInt_t fFirstRow; // First TPC row crossed by the track
242 UInt_t fLastRow; // Last TPC row crossed by the track
243 Float_t fTrackPt; // Pt of the track
a8ffd46b 244};
245
4aa41877 246void AliHLTHoughTransformerRow::CreateHistograms(Int_t nxbin,Float_t xmin,Float_t xmax,
0bd0c1ef 247 Int_t nybin,Float_t ymin,Float_t ymax)
248{
917e711b 249 //Create the histograms (parameter space)
250 //nxbin = #bins in X
251 //nybin = #bins in Y
252 //xmin xmax ymin ymax = histogram limits in X and Y
7e4eb453 253 if(fLastTransformer) {
4aa41877 254 SetTransformerArrays((AliHLTHoughTransformerRow *)fLastTransformer);
7e4eb453 255 return;
256 }
4aa41877 257 fParamSpace = new AliHLTHistogram*[GetNEtaSegments()];
0bd0c1ef 258
259 Char_t histname[256];
260 for(Int_t i=0; i<GetNEtaSegments(); i++)
261 {
262 sprintf(histname,"paramspace_%d",i);
4aa41877 263 fParamSpace[i] = new AliHLTHistogram(histname,"",nxbin,xmin,xmax,nybin,ymin,ymax);
0bd0c1ef 264 }
a8ffd46b 265#ifdef do_mc
266 {
4aa41877 267 AliHLTHistogram *hist = fParamSpace[0];
a8ffd46b 268 Int_t ncellsx = (hist->GetNbinsX()+3)/2;
269 Int_t ncellsy = (hist->GetNbinsY()+3)/2;
270 Int_t ncells = ncellsx*ncellsy;
271 if(!fTrackID)
272 {
4aa41877 273 LOG(AliHLTLog::kInformational,"AliHLTHoughTransformerRow::CreateHistograms()","")
274 <<"Transformer: Allocating "<<GetNEtaSegments()*ncells*sizeof(AliHLTTrackIndex)<<" bytes to fTrackID"<<ENDLOG;
275 fTrackID = new AliHLTTrackIndex*[GetNEtaSegments()];
a8ffd46b 276 for(Int_t i=0; i<GetNEtaSegments(); i++)
4aa41877 277 fTrackID[i] = new AliHLTTrackIndex[ncells];
a8ffd46b 278 }
279 }
0bd0c1ef 280#endif
4aa41877 281 AliHLTHistogram *hist = fParamSpace[0];
0bd0c1ef 282 Int_t ncells = (hist->GetNbinsX()+2)*(hist->GetNbinsY()+2);
a8ffd46b 283 if(!fGapCount)
0bd0c1ef 284 {
4aa41877 285 LOG(AliHLTLog::kInformational,"AliHLTHoughTransformerRow::CreateHistograms()","")
a8ffd46b 286 <<"Transformer: Allocating "<<GetNEtaSegments()*ncells*sizeof(UChar_t)<<" bytes to fGapCount"<<ENDLOG;
287 fGapCount = new UChar_t*[GetNEtaSegments()];
0bd0c1ef 288 for(Int_t i=0; i<GetNEtaSegments(); i++)
a8ffd46b 289 fGapCount[i] = new UChar_t[ncells];
0bd0c1ef 290 }
a8ffd46b 291 if(!fCurrentRowCount)
0bd0c1ef 292 {
4aa41877 293 LOG(AliHLTLog::kInformational,"AliHLTHoughTransformerRow::CreateHistograms()","")
a8ffd46b 294 <<"Transformer: Allocating "<<GetNEtaSegments()*ncells*sizeof(UChar_t)<<" bytes to fCurrentRowCount"<<ENDLOG;
295 fCurrentRowCount = new UChar_t*[GetNEtaSegments()];
0bd0c1ef 296 for(Int_t i=0; i<GetNEtaSegments(); i++)
a8ffd46b 297 fCurrentRowCount[i] = new UChar_t[ncells];
0bd0c1ef 298 }
a8ffd46b 299 if(!fPrevBin)
0bd0c1ef 300 {
4aa41877 301 LOG(AliHLTLog::kInformational,"AliHLTHoughTransformerRow::CreateHistograms()","")
a8ffd46b 302 <<"Transformer: Allocating "<<GetNEtaSegments()*ncells*sizeof(UChar_t)<<" bytes to fPrevBin"<<ENDLOG;
303 fPrevBin = new UChar_t*[GetNEtaSegments()];
0bd0c1ef 304 for(Int_t i=0; i<GetNEtaSegments(); i++)
a8ffd46b 305 fPrevBin[i] = new UChar_t[ncells];
306 }
307 if(!fNextBin)
308 {
4aa41877 309 LOG(AliHLTLog::kInformational,"AliHLTHoughTransformerRow::CreateHistograms()","")
a8ffd46b 310 <<"Transformer: Allocating "<<GetNEtaSegments()*ncells*sizeof(UChar_t)<<" bytes to fNextBin"<<ENDLOG;
311 fNextBin = new UChar_t*[GetNEtaSegments()];
312 for(Int_t i=0; i<GetNEtaSegments(); i++)
313 fNextBin[i] = new UChar_t[ncells];
314 }
315 Int_t ncellsy = hist->GetNbinsY()+2;
316 if(!fNextRow)
317 {
4aa41877 318 LOG(AliHLTLog::kInformational,"AliHLTHoughTransformerRow::CreateHistograms()","")
a8ffd46b 319 <<"Transformer: Allocating "<<GetNEtaSegments()*ncellsy*sizeof(UChar_t)<<" bytes to fNextRow"<<ENDLOG;
320 fNextRow = new UChar_t*[GetNEtaSegments()];
321 for(Int_t i=0; i<GetNEtaSegments(); i++)
322 fNextRow[i] = new UChar_t[ncellsy];
0bd0c1ef 323 }
324
a8ffd46b 325 if(!fTrackNRows)
de3c3890 326 {
4aa41877 327 LOG(AliHLTLog::kInformational,"AliHLTHoughTransformerRow::CreateHistograms()","")
a8ffd46b 328 <<"Transformer: Allocating "<<ncells*sizeof(UChar_t)<<" bytes to fTrackNRows"<<ENDLOG;
329 fTrackNRows = new UChar_t[ncells];
4aa41877 330 LOG(AliHLTLog::kInformational,"AliHLTHoughTransformerRow::CreateHistograms()","")
a8ffd46b 331 <<"Transformer: Allocating "<<ncells*sizeof(UChar_t)<<" bytes to fTrackFirstRow"<<ENDLOG;
332 fTrackFirstRow = new UChar_t[ncells];
4aa41877 333 LOG(AliHLTLog::kInformational,"AliHLTHoughTransformerRow::CreateHistograms()","")
a8ffd46b 334 <<"Transformer: Allocating "<<ncells*sizeof(UChar_t)<<" bytes to fTrackLastRow"<<ENDLOG;
335 fTrackLastRow = new UChar_t[ncells];
4aa41877 336 LOG(AliHLTLog::kInformational,"AliHLTHoughTransformerRow::CreateHistograms()","")
a8ffd46b 337 <<"Transformer: Allocating "<<ncells*sizeof(UChar_t)<<" bytes to fInitialGapCount"<<ENDLOG;
338 fInitialGapCount = new UChar_t[ncells];
339
de3c3890 340
4aa41877 341 AliHLTHoughTrack track;
de3c3890 342 Int_t xmin = hist->GetFirstXbin();
343 Int_t xmax = hist->GetLastXbin();
a8ffd46b 344 Int_t xmiddle = (hist->GetNbinsX()+1)/2;
de3c3890 345 Int_t ymin = hist->GetFirstYbin();
346 Int_t ymax = hist->GetLastYbin();
347 Int_t nxbins = hist->GetNbinsX()+2;
a8ffd46b 348 Int_t nxgrid = (hist->GetNbinsX()+3)/2+1;
349 Int_t nygrid = hist->GetNbinsY()+3;
350
4aa41877 351 AliHLTTrackLength *tracklength = new AliHLTTrackLength[nxgrid*nygrid];
352 memset(tracklength,0,nxgrid*nygrid*sizeof(AliHLTTrackLength));
a8ffd46b 353
74e77d1b 354 for(Int_t ybin=ymin-1; ybin<=(ymax+1); ybin++)
de3c3890 355 {
a8ffd46b 356 for(Int_t xbin=xmin-1; xbin<=xmiddle; xbin++)
de3c3890 357 {
a8ffd46b 358 fTrackNRows[xbin + ybin*nxbins] = 255;
359 for(Int_t deltay = 0; deltay <= 1; deltay++) {
360 for(Int_t deltax = 0; deltax <= 1; deltax++) {
361
4aa41877 362 AliHLTTrackLength *curtracklength = &tracklength[(xbin + deltax) + (ybin + deltay)*nxgrid];
de3c3890 363 UInt_t maxfirstrow = 0;
364 UInt_t maxlastrow = 0;
a8ffd46b 365 Float_t maxtrackpt = 0;
366 if(curtracklength->fIsFilled) {
367 maxfirstrow = curtracklength->fFirstRow;
368 maxlastrow = curtracklength->fLastRow;
369 maxtrackpt = curtracklength->fTrackPt;
370 }
371 else {
372 Float_t xtrack = hist->GetPreciseBinCenterX((Float_t)xbin+0.5*(Float_t)(2*deltax-1));
373 Float_t ytrack = hist->GetPreciseBinCenterY((Float_t)ybin+0.5*(Float_t)(2*deltay-1));
374
375 Float_t psi = atan((xtrack-ytrack)/(fgBeta1-fgBeta2));
376 Float_t kappa = 2.0*(xtrack*cos(psi)-fgBeta1*sin(psi));
377 track.SetTrackParameters(kappa,psi,1);
7e4eb453 378 maxtrackpt = track.GetPt();
4aa41877 379 if(maxtrackpt < 0.9*0.1*AliHLTTransform::GetSolenoidField())
7e4eb453 380 {
381 maxfirstrow = maxlastrow = 0;
382 curtracklength->fIsFilled = kTRUE;
383 curtracklength->fFirstRow = maxfirstrow;
384 curtracklength->fLastRow = maxlastrow;
385 curtracklength->fTrackPt = maxtrackpt;
386 }
387 else
388 {
389 Bool_t firstrow = kFALSE;
390 UInt_t curfirstrow = 0;
391 UInt_t curlastrow = 0;
a8ffd46b 392
7e4eb453 393 Double_t centerx = track.GetCenterX();
394 Double_t centery = track.GetCenterY();
395 Double_t radius = track.GetRadius();
a8ffd46b 396
4aa41877 397 for(Int_t j=AliHLTTransform::GetFirstRow(0); j<=AliHLTTransform::GetLastRow(5); j++)
7e4eb453 398 {
399 Float_t hit[3];
400 // if(!track.GetCrossingPoint(j,hit)) continue;
4aa41877 401 hit[0] = AliHLTTransform::Row2X(j);
7e4eb453 402 Double_t aa = (hit[0] - centerx)*(hit[0] - centerx);
403 Double_t r2 = radius*radius;
404 if(aa > r2)
405 continue;
406
407 Double_t aa2 = sqrt(r2 - aa);
408 Double_t y1 = centery + aa2;
409 Double_t y2 = centery - aa2;
410 hit[1] = y1;
411 if(fabs(y2) < fabs(y1)) hit[1] = y2;
a8ffd46b 412
7e4eb453 413 hit[2] = 0;
a8ffd46b 414
4aa41877 415 AliHLTTransform::LocHLT2Raw(hit,0,j);
7e4eb453 416 hit[1] += 0.5;
4aa41877 417 if(hit[1]>=0 && hit[1]<AliHLTTransform::GetNPads(j))
7e4eb453 418 {
419 if(!firstrow) {
420 curfirstrow = j;
421 firstrow = kTRUE;
422 }
423 curlastrow = j;
424 }
425 else {
426 if(firstrow) {
427 firstrow = kFALSE;
428 if((curlastrow-curfirstrow) >= (maxlastrow-maxfirstrow)) {
429 maxfirstrow = curfirstrow;
430 maxlastrow = curlastrow;
431 }
432 }
a8ffd46b 433 }
de3c3890 434 }
7e4eb453 435 if((curlastrow-curfirstrow) >= (maxlastrow-maxfirstrow)) {
436 maxfirstrow = curfirstrow;
437 maxlastrow = curlastrow;
de3c3890 438 }
a8ffd46b 439
7e4eb453 440 curtracklength->fIsFilled = kTRUE;
441 curtracklength->fFirstRow = maxfirstrow;
442 curtracklength->fLastRow = maxlastrow;
443 curtracklength->fTrackPt = maxtrackpt;
444 }
de3c3890 445 }
a8ffd46b 446 if((maxlastrow-maxfirstrow) < fTrackNRows[xbin + ybin*nxbins]) {
447 fTrackNRows[xbin + ybin*nxbins] = maxlastrow-maxfirstrow;
448 fInitialGapCount[xbin + ybin*nxbins] = 1;
9c0ddf1e 449 if((maxlastrow-maxfirstrow+1)<=MIN_TRACK_LENGTH)
f644512a 450 fInitialGapCount[xbin + ybin*nxbins] = MAX_N_GAPS+1;
4aa41877 451 if(maxtrackpt < 0.9*0.1*AliHLTTransform::GetSolenoidField())
a8ffd46b 452 fInitialGapCount[xbin + ybin*nxbins] = MAX_N_GAPS;
453 fTrackFirstRow[xbin + ybin*nxbins] = maxfirstrow;
454 fTrackLastRow[xbin + ybin*nxbins] = maxlastrow;
455
456 Int_t xmirror = xmax - xbin + 1;
457 Int_t ymirror = ymax - ybin + 1;
458 fTrackNRows[xmirror + ymirror*nxbins] = fTrackNRows[xbin + ybin*nxbins];
459 fInitialGapCount[xmirror + ymirror*nxbins] = fInitialGapCount[xbin + ybin*nxbins];
460 fTrackFirstRow[xmirror + ymirror*nxbins] = fTrackFirstRow[xbin + ybin*nxbins];
461 fTrackLastRow[xmirror + ymirror*nxbins] = fTrackLastRow[xbin + ybin*nxbins];
de3c3890 462 }
463 }
464 }
a8ffd46b 465 // cout<<" fTrackNRows "<<(Int_t)fInitialGapCount[xbin + ybin*nxbins]<<" "<<xbin<<" "<<ybin<<" "<<(Int_t)fTrackNRows[xbin + ybin*nxbins]<<" "<<(Int_t)fTrackFirstRow[xbin + ybin*nxbins]<<" "<<(Int_t)fTrackLastRow[xbin + ybin*nxbins]<<" "<<endl;
de3c3890 466 }
467 }
a109e73e 468 delete [] tracklength;
de3c3890 469 }
0bd0c1ef 470
a8ffd46b 471 if(!fStartPadParams)
472 {
4aa41877 473 Int_t nrows = AliHLTTransform::GetLastRow(5) - AliHLTTransform::GetFirstRow(0) + 1;
474 LOG(AliHLTLog::kInformational,"AliHLTHoughTransformerRow::CreateHistograms()","")
475 <<"Transformer: Allocating about "<<nrows*100*sizeof(AliHLTPadHoughParams)<<" bytes to fStartPadParams"<<ENDLOG;
476 fStartPadParams = new AliHLTPadHoughParams*[nrows];
477 LOG(AliHLTLog::kInformational,"AliHLTHoughTransformerRow::CreateHistograms()","")
478 <<"Transformer: Allocating about "<<nrows*100*sizeof(AliHLTPadHoughParams)<<" bytes to fEndPadParams"<<ENDLOG;
479 fEndPadParams = new AliHLTPadHoughParams*[nrows];
480 LOG(AliHLTLog::kInformational,"AliHLTHoughTransformerRow::CreateHistograms()","")
f644512a 481 <<"Transformer: Allocating about "<<nrows*100*sizeof(Float_t)<<" bytes to fLUTr"<<ENDLOG;
482 fLUTr = new Float_t*[nrows];
a8ffd46b 483
484 Float_t beta1 = fgBeta1;
485 Float_t beta2 = fgBeta2;
486 Float_t beta1minusbeta2 = fgBeta1 - fgBeta2;
487 Float_t ymin = hist->GetYmin();
488 Float_t histbin = hist->GetBinWidthY();
489 Float_t xmin = hist->GetXmin();
490 Float_t xmax = hist->GetXmax();
491 Float_t xbin = (xmax-xmin)/hist->GetNbinsX();
492 Int_t firstbinx = hist->GetFirstXbin();
493 Int_t lastbinx = hist->GetLastXbin();
494 Int_t nbinx = hist->GetNbinsX()+2;
495 Int_t firstbin = hist->GetFirstYbin();
496 Int_t lastbin = hist->GetLastYbin();
4aa41877 497 for(Int_t i=AliHLTTransform::GetFirstRow(0); i<=AliHLTTransform::GetLastRow(5); i++)
a8ffd46b 498 {
4aa41877 499 Int_t npads = AliHLTTransform::GetNPads(i);
500 Int_t ipatch = AliHLTTransform::GetPatch(i);
501 Double_t padpitch = AliHLTTransform::GetPadPitchWidth(ipatch);
502 Float_t x = AliHLTTransform::Row2X(i);
a8ffd46b 503 Float_t x2 = x*x;
504
4aa41877 505 fStartPadParams[i] = new AliHLTPadHoughParams[npads];
506 fEndPadParams[i] = new AliHLTPadHoughParams[npads];
f644512a 507 fLUTr[i] = new Float_t[npads];
a8ffd46b 508 for(Int_t pad=0; pad<npads; pad++)
509 {
510 Float_t y = (pad-0.5*(npads-1))*padpitch;
f644512a 511 fLUTr[i][pad] = sqrt(x2 + y*y);
a8ffd46b 512 Float_t starty = (pad-0.5*npads)*padpitch;
513 Float_t r1 = x2 + starty*starty;
514 Float_t xoverr1 = x/r1;
515 Float_t startyoverr1 = starty/r1;
516 Float_t endy = (pad-0.5*(npads-2))*padpitch;
517 Float_t r2 = x2 + endy*endy;
518 Float_t xoverr2 = x/r2;
519 Float_t endyoverr2 = endy/r2;
520 Float_t a1 = beta1minusbeta2/(xoverr1-beta2);
521 Float_t b1 = (xoverr1-beta1)/(xoverr1-beta2);
522 Float_t a2 = beta1minusbeta2/(xoverr2-beta2);
523 Float_t b2 = (xoverr2-beta1)/(xoverr2-beta2);
524
525 Float_t alpha1 = (a1*startyoverr1+b1*ymin-xmin)/xbin;
526 Float_t deltaalpha1 = b1*histbin/xbin;
527 if(b1<0)
528 alpha1 += deltaalpha1;
529 Float_t alpha2 = (a2*endyoverr2+b2*ymin-xmin)/xbin;
530 Float_t deltaalpha2 = b2*histbin/xbin;
531 if(b2>=0)
532 alpha2 += deltaalpha2;
533
534 fStartPadParams[i][pad].fAlpha = alpha1;
535 fStartPadParams[i][pad].fDeltaAlpha = deltaalpha1;
536 fEndPadParams[i][pad].fAlpha = alpha2;
537 fEndPadParams[i][pad].fDeltaAlpha = deltaalpha2;
538
539 //Find the first and last bin rows to be filled
540 Bool_t binfound1 = kFALSE;
541 Bool_t binfound2 = kFALSE;
542 Int_t firstbin1 = lastbin;
543 Int_t firstbin2 = lastbin;
544 Int_t lastbin1 = firstbin;
545 Int_t lastbin2 = firstbin;
546 for(Int_t b=firstbin; b<=lastbin; b++, alpha1 += deltaalpha1, alpha2 += deltaalpha2)
547 {
548 Int_t binx1 = 1 + (Int_t)alpha1;
549 if(binx1<=lastbinx) {
550 UChar_t initialgapcount;
551 if(binx1>=firstbinx)
552 initialgapcount = fInitialGapCount[binx1 + b*nbinx];
553 else
554 initialgapcount = fInitialGapCount[firstbinx + b*nbinx];
555 if(initialgapcount != MAX_N_GAPS) {
556 if(!binfound1) {
557 firstbin1 = b;
558 binfound1 = kTRUE;
559 }
560 lastbin1 = b;
561 }
562 }
563 Int_t binx2 = 1 + (Int_t)alpha2;
564 if(binx2>=firstbin) {
565 UChar_t initialgapcount;
566 if(binx2<=lastbinx)
567 initialgapcount = fInitialGapCount[binx2 + b*nbinx];
568 else
569 initialgapcount = fInitialGapCount[lastbinx + b*nbinx];
570 if(initialgapcount != MAX_N_GAPS) {
571 if(!binfound2) {
572 firstbin2 = b;
573 binfound2 = kTRUE;
574 }
575 lastbin2 = b;
576 }
577 }
578 }
579 fStartPadParams[i][pad].fFirstBin=firstbin1;
580 fStartPadParams[i][pad].fLastBin=lastbin1;
581 fEndPadParams[i][pad].fFirstBin=firstbin2;
582 fEndPadParams[i][pad].fLastBin=lastbin2;
583 }
584 }
585 }
586
0bd0c1ef 587 //create lookup table for z of the digits
a8ffd46b 588 if(!fLUTforwardZ)
589 {
4aa41877 590 Int_t ntimebins = AliHLTTransform::GetNTimeBins();
591 LOG(AliHLTLog::kInformational,"AliHLTHoughTransformerRow::CreateHistograms()","")
a8ffd46b 592 <<"Transformer: Allocating "<<ntimebins*sizeof(Float_t)<<" bytes to fLUTforwardZ"<<ENDLOG;
593 fLUTforwardZ = new Float_t[ntimebins];
4aa41877 594 LOG(AliHLTLog::kInformational,"AliHLTHoughTransformerRow::CreateHistograms()","")
a8ffd46b 595 <<"Transformer: Allocating "<<ntimebins*sizeof(Float_t)<<" bytes to fLUTbackwardZ"<<ENDLOG;
596 fLUTbackwardZ = new Float_t[ntimebins];
a8ffd46b 597 for(Int_t i=0; i<ntimebins; i++){
598 Float_t z;
4aa41877 599 z=AliHLTTransform::GetZFast(0,i,GetZVertex());
a8ffd46b 600 fLUTforwardZ[i]=z;
4aa41877 601 z=AliHLTTransform::GetZFast(18,i,GetZVertex());
a8ffd46b 602 fLUTbackwardZ[i]=z;
a8ffd46b 603 }
604 }
0bd0c1ef 605}
606
4aa41877 607void AliHLTHoughTransformerRow::Reset()
0bd0c1ef 608{
917e711b 609 //Reset all the histograms. Should be done when processing new slice
7e4eb453 610 if(fLastTransformer) return;
0bd0c1ef 611 if(!fParamSpace)
612 {
4aa41877 613 LOG(AliHLTLog::kWarning,"AliHLTHoughTransformer::Reset","Histograms")
0bd0c1ef 614 <<"No histograms to reset"<<ENDLOG;
615 return;
616 }
617
618 for(Int_t i=0; i<GetNEtaSegments(); i++)
619 fParamSpace[i]->Reset();
620
621#ifdef do_mc
a8ffd46b 622 {
4aa41877 623 AliHLTHistogram *hist = fParamSpace[0];
a8ffd46b 624 Int_t ncellsx = (hist->GetNbinsX()+3)/2;
625 Int_t ncellsy = (hist->GetNbinsY()+3)/2;
626 Int_t ncells = ncellsx*ncellsy;
627 for(Int_t i=0; i<GetNEtaSegments(); i++)
4aa41877 628 memset(fTrackID[i],0,ncells*sizeof(AliHLTTrackIndex));
a8ffd46b 629 }
0bd0c1ef 630#endif
4aa41877 631 AliHLTHistogram *hist = fParamSpace[0];
0bd0c1ef 632 Int_t ncells = (hist->GetNbinsX()+2)*(hist->GetNbinsY()+2);
633 for(Int_t i=0; i<GetNEtaSegments(); i++)
634 {
a8ffd46b 635 memcpy(fGapCount[i],fInitialGapCount,ncells*sizeof(UChar_t));
636 memcpy(fCurrentRowCount[i],fTrackFirstRow,ncells*sizeof(UChar_t));
0bd0c1ef 637 }
638}
639
4aa41877 640Int_t AliHLTHoughTransformerRow::GetEtaIndex(Double_t eta) const
0bd0c1ef 641{
642 //Return the histogram index of the corresponding eta.
643
644 Double_t etaslice = (GetEtaMax() - GetEtaMin())/GetNEtaSegments();
645 Double_t index = (eta-GetEtaMin())/etaslice;
646 return (Int_t)index;
647}
648
4aa41877 649inline AliHLTHistogram *AliHLTHoughTransformerRow::GetHistogram(Int_t etaindex)
0bd0c1ef 650{
917e711b 651 // Return a pointer to the histogram which contains etaindex eta slice
652 if(!fParamSpace || etaindex >= GetNEtaSegments() || etaindex < 0)
0bd0c1ef 653 return 0;
917e711b 654 if(!fParamSpace[etaindex])
0bd0c1ef 655 return 0;
917e711b 656 return fParamSpace[etaindex];
0bd0c1ef 657}
658
4aa41877 659Double_t AliHLTHoughTransformerRow::GetEta(Int_t etaindex,Int_t /*slice*/) const
0bd0c1ef 660{
917e711b 661 // Return eta calculated in the middle of the eta slice
662 Double_t etaslice = (GetEtaMax()-GetEtaMin())/GetNEtaSegments();
0bd0c1ef 663 Double_t eta=0;
917e711b 664 eta=(Double_t)((etaindex+0.5)*etaslice);
0bd0c1ef 665 return eta;
666}
667
4aa41877 668void AliHLTHoughTransformerRow::TransformCircle()
a8ffd46b 669{
ce4c28d5 670 // This method contains the hough transformation
671 // Depending on the option selected, it reads as an input
672 // either the preloaded array with digits or directly raw
673 // data stream (option fast_raw)
a8ffd46b 674 if(GetDataPointer())
675 TransformCircleFromDigitArray();
676 else if(fTPCRawStream)
677 TransformCircleFromRawStream();
678}
4aa41877 679void AliHLTHoughTransformerRow::TransformCircleFromDigitArray()
0bd0c1ef 680{
917e711b 681 //Do the Hough Transform
de3c3890 682
f644512a 683 //Load the parameters used by the fast calculation of eta
684 Double_t etaparam1 = GetEtaCalcParam1();
685 Double_t etaparam2 = GetEtaCalcParam2();
686 Double_t etaparam3 = GetEtaCalcParam3();
687
917e711b 688 Int_t netasegments = GetNEtaSegments();
689 Double_t etamin = GetEtaMin();
690 Double_t etaslice = (GetEtaMax() - etamin)/netasegments;
0bd0c1ef 691
917e711b 692 Int_t lowerthreshold = GetLowerThreshold();
0bd0c1ef 693
694 //Assumes that all the etaslice histos are the same!
4aa41877 695 AliHLTHistogram *h = fParamSpace[0];
a8ffd46b 696 Int_t firstbiny = h->GetFirstYbin();
697 Int_t firstbinx = h->GetFirstXbin();
698 Int_t lastbinx = h->GetLastXbin();
0bd0c1ef 699 Int_t nbinx = h->GetNbinsX()+2;
700
917e711b 701 UChar_t lastpad;
3cb5c013 702 Int_t lastetaindex=-1;
4aa41877 703 AliHLTEtaRow *etaclust = new AliHLTEtaRow[netasegments];
0bd0c1ef 704
4aa41877 705 AliHLTDigitRowData *tempPt = GetDataPointer();
0bd0c1ef 706 if(!tempPt)
707 {
4aa41877 708 LOG(AliHLTLog::kError,"AliHLTHoughTransformer::TransformCircle","Data")
0bd0c1ef 709 <<"No input data "<<ENDLOG;
710 return;
711 }
712
713 Int_t ipatch = GetPatch();
a8ffd46b 714 Int_t ilastpatch = GetLastPatch();
0bd0c1ef 715 Int_t islice = GetSlice();
a8ffd46b 716 Float_t *lutz;
f644512a 717 if(islice < 18)
a8ffd46b 718 lutz = fLUTforwardZ;
f644512a 719 else
a8ffd46b 720 lutz = fLUTbackwardZ;
0bd0c1ef 721
0bd0c1ef 722 //Loop over the padrows:
4aa41877 723 for(UChar_t i=AliHLTTransform::GetFirstRow(ipatch); i<=AliHLTTransform::GetLastRow(ipatch); i++)
0bd0c1ef 724 {
917e711b 725 lastpad = 255;
0bd0c1ef 726 //Flush eta clusters array
4aa41877 727 memset(etaclust,0,netasegments*sizeof(AliHLTEtaRow));
0bd0c1ef 728
f644512a 729 Float_t radius=0;
0bd0c1ef 730
0bd0c1ef 731 //Get the data on this padrow:
4aa41877 732 AliHLTDigitData *digPt = tempPt->fDigitData;
0bd0c1ef 733 if((Int_t)i != (Int_t)tempPt->fRow)
734 {
4aa41877 735 LOG(AliHLTLog::kError,"AliHLTHoughTransformerRow::TransformCircle","Data")
736 <<"AliHLTHoughTransform::TransformCircle : Mismatching padrow numbering "<<(Int_t)i<<" "<<(Int_t)tempPt->fRow<<ENDLOG;
0bd0c1ef 737 continue;
738 }
739 // cout<<" Starting row "<<i<<endl;
740 //Loop over the data on this padrow:
741 for(UInt_t j=0; j<tempPt->fNDigit; j++)
742 {
743 UShort_t charge = digPt[j].fCharge;
917e711b 744 if((Int_t)charge <= lowerthreshold)
0bd0c1ef 745 continue;
746 UChar_t pad = digPt[j].fPad;
747 UShort_t time = digPt[j].fTime;
de3c3890 748
917e711b 749 if(pad != lastpad)
0bd0c1ef 750 {
f644512a 751 radius = fLUTr[(Int_t)i][(Int_t)pad];
917e711b 752 lastetaindex = -1;
0bd0c1ef 753 }
754
a8ffd46b 755 Float_t z = lutz[(Int_t)time];
f644512a 756 Double_t radiuscorr = radius*(1.+etaparam3*radius*radius);
757 Double_t zovr = z/radiuscorr;
758 Double_t eta = (etaparam1-etaparam2*fabs(zovr))*zovr;
0bd0c1ef 759 //Get the corresponding index, which determines which histogram to fill:
917e711b 760 Int_t etaindex = (Int_t)((eta-etamin)/etaslice);
0bd0c1ef 761
762#ifndef do_mc
917e711b 763 if(etaindex == lastetaindex) continue;
0bd0c1ef 764#endif
917e711b 765 // cout<<" Digit at patch "<<ipatch<<" row "<<i<<" pad "<<(Int_t)pad<<" time "<<time<<" etaslice "<<etaindex<<endl;
0bd0c1ef 766
917e711b 767 if(etaindex < 0 || etaindex >= netasegments)
0bd0c1ef 768 continue;
769
917e711b 770 if(!etaclust[etaindex].fIsFound)
0bd0c1ef 771 {
917e711b 772 etaclust[etaindex].fStartPad = pad;
773 etaclust[etaindex].fEndPad = pad;
917e711b 774 etaclust[etaindex].fIsFound = 1;
0bd0c1ef 775#ifdef do_mc
8cc2b1a5 776 FillClusterMCLabels(digPt[j],&etaclust[etaindex]);
0bd0c1ef 777#endif
778 continue;
779 }
780 else
781 {
917e711b 782 if(pad <= (etaclust[etaindex].fEndPad + 1))
0bd0c1ef 783 {
917e711b 784 etaclust[etaindex].fEndPad = pad;
0bd0c1ef 785#ifdef do_mc
8cc2b1a5 786 FillClusterMCLabels(digPt[j],&etaclust[etaindex]);
0bd0c1ef 787#endif
788 }
789 else
790 {
a8ffd46b 791 FillCluster(i,etaindex,etaclust,ilastpatch,firstbinx,lastbinx,nbinx,firstbiny);
0bd0c1ef 792
917e711b 793 etaclust[etaindex].fStartPad = pad;
794 etaclust[etaindex].fEndPad = pad;
0bd0c1ef 795
796#ifdef do_mc
a8ffd46b 797 memset(etaclust[etaindex].fMcLabels,0,MaxTrack);
8cc2b1a5 798 FillClusterMCLabels(digPt[j],&etaclust[etaindex]);
0bd0c1ef 799#endif
800 }
801 }
917e711b 802 lastpad = pad;
803 lastetaindex = etaindex;
0bd0c1ef 804 }
805 //Write remaining clusters
917e711b 806 for(Int_t etaindex = 0;etaindex < netasegments;etaindex++)
0bd0c1ef 807 {
808 //Check for empty row
917e711b 809 if((etaclust[etaindex].fStartPad == 0) && (etaclust[etaindex].fEndPad == 0)) continue;
0bd0c1ef 810
a8ffd46b 811 FillCluster(i,etaindex,etaclust,ilastpatch,firstbinx,lastbinx,nbinx,firstbiny);
0bd0c1ef 812 }
813
814 //Move the data pointer to the next padrow:
4aa41877 815 AliHLTMemHandler::UpdateRowPointer(tempPt);
0bd0c1ef 816 }
817
917e711b 818 delete [] etaclust;
0bd0c1ef 819}
820
4aa41877 821void AliHLTHoughTransformerRow::TransformCircleFromRawStream()
0bd0c1ef 822{
a8ffd46b 823 //Do the Hough Transform
824
f644512a 825 //Load the parameters used by the fast calculation of eta
826 Double_t etaparam1 = GetEtaCalcParam1();
827 Double_t etaparam2 = GetEtaCalcParam2();
828 Double_t etaparam3 = GetEtaCalcParam3();
829
a8ffd46b 830 Int_t netasegments = GetNEtaSegments();
831 Double_t etamin = GetEtaMin();
832 Double_t etaslice = (GetEtaMax() - etamin)/netasegments;
833
834 Int_t lowerthreshold = GetLowerThreshold();
835
836 //Assumes that all the etaslice histos are the same!
4aa41877 837 AliHLTHistogram *h = fParamSpace[0];
a8ffd46b 838 Int_t firstbiny = h->GetFirstYbin();
839 Int_t firstbinx = h->GetFirstXbin();
840 Int_t lastbinx = h->GetLastXbin();
841 Int_t nbinx = h->GetNbinsX()+2;
842
3cb5c013 843 Int_t lastetaindex = -1;
4aa41877 844 AliHLTEtaRow *etaclust = new AliHLTEtaRow[netasegments];
a8ffd46b 845
846 if(!fTPCRawStream)
0bd0c1ef 847 {
4aa41877 848 LOG(AliHLTLog::kError,"AliHLTHoughTransformer::TransformCircle","Data")
a8ffd46b 849 <<"No input data "<<ENDLOG;
850 return;
851 }
852
853 Int_t ipatch = GetPatch();
4aa41877 854 UChar_t rowmin = AliHLTTransform::GetFirstRowOnDDL(ipatch);
855 UChar_t rowmax = AliHLTTransform::GetLastRowOnDDL(ipatch);
856 // Int_t ntimebins = AliHLTTransform::GetNTimeBins();
a8ffd46b 857 Int_t ilastpatch = GetLastPatch();
858 Int_t islice = GetSlice();
859 Float_t *lutz;
f644512a 860 if(islice < 18)
a8ffd46b 861 lutz = fLUTforwardZ;
f644512a 862 else
a8ffd46b 863 lutz = fLUTbackwardZ;
a8ffd46b 864
865 //Flush eta clusters array
4aa41877 866 memset(etaclust,0,netasegments*sizeof(AliHLTEtaRow));
a8ffd46b 867
868 UChar_t i=0;
869 Int_t npads=0;
f644512a 870 Float_t radius=0;
a8ffd46b 871 UChar_t pad=0;
872
873 //Loop over the rawdata:
874 while (fTPCRawStream->Next()) {
875
876 if(fTPCRawStream->IsNewSector() || fTPCRawStream->IsNewRow()) {
877
878 //Write remaining clusters
879 for(Int_t etaindex = 0;etaindex < netasegments;etaindex++)
880 {
881 //Check for empty row
882 if((etaclust[etaindex].fStartPad == 0) && (etaclust[etaindex].fEndPad == 0)) continue;
883
884 FillCluster(i,etaindex,etaclust,ilastpatch,firstbinx,lastbinx,nbinx,firstbiny);
885 }
886
887 Int_t sector=fTPCRawStream->GetSector();
888 Int_t row=fTPCRawStream->GetRow();
889 Int_t slice,srow;
4aa41877 890 AliHLTTransform::Sector2Slice(slice,srow,sector,row);
a8ffd46b 891 if(slice!=islice){
4aa41877 892 LOG(AliHLTLog::kError,"AliHLTDDLDataFileHandler::DDLDigits2Memory","Slice")
893 <<AliHLTLog::kDec<<"Found slice "<<slice<<", expected "<<islice<<ENDLOG;
a8ffd46b 894 continue;
895 }
896
897 i=(UChar_t)srow;
4aa41877 898 npads = AliHLTTransform::GetNPads(srow)-1;
a8ffd46b 899
900 //Flush eta clusters array
4aa41877 901 memset(etaclust,0,netasegments*sizeof(AliHLTEtaRow));
a8ffd46b 902
f644512a 903 radius=0;
a8ffd46b 904
905 }
906
907 if((i<rowmin)||(i>rowmax))continue;
908
909 // cout<<" Starting row "<<(UInt_t)i<<endl;
910 //Loop over the data on this padrow:
911 if(fTPCRawStream->IsNewRow() || fTPCRawStream->IsNewPad()) {
912 pad=fTPCRawStream->GetPad();
913 /*
914 if((pad<0)||(pad>=(npads+1))){
4aa41877 915 LOG(AliHLTLog::kError,"AliHLTDDLDataFileHandler::DDLDigits2Memory","Pad")
916 <<AliHLTLog::kDec<<"Pad value out of bounds "<<pad<<" "
a8ffd46b 917 <<npads+1<<ENDLOG;
918 continue;
919 }
920 */
f644512a 921 radius = fLUTr[(Int_t)i][(Int_t)pad];
a8ffd46b 922 lastetaindex = -1;
923 }
924
925 UShort_t time=fTPCRawStream->GetTime();
926 /*
927 if((time<0)||(time>=ntimebins)){
4aa41877 928 LOG(AliHLTLog::kError,"AliHLTDDLDataFileHandler::DDLDigits2Memory","Time")
929 <<AliHLTLog::kDec<<"Time out of bounds "<<time<<" "
930 <<AliHLTTransform::GetNTimeBins()<<ENDLOG;
a8ffd46b 931 continue;
932 }
933 */
934
935 if(fTPCRawStream->GetSignal() <= lowerthreshold)
936 continue;
937
a8ffd46b 938 Float_t z = lutz[(Int_t)time];
f644512a 939 Double_t radiuscorr = radius*(1.+etaparam3*radius*radius);
940 Double_t zovr = z/radiuscorr;
941 Double_t eta = (etaparam1-etaparam2*fabs(zovr))*zovr;
a8ffd46b 942 //Get the corresponding index, which determines which histogram to fill:
943 Int_t etaindex = (Int_t)((eta-etamin)/etaslice);
944
945#ifndef do_mc
946 if(etaindex == lastetaindex) continue;
947#endif
948 // cout<<" Digit at patch "<<ipatch<<" row "<<i<<" pad "<<(Int_t)pad<<" time "<<time<<" etaslice "<<etaindex<<endl;
949
950 if(etaindex < 0 || etaindex >= netasegments)
951 continue;
952
953 if(!etaclust[etaindex].fIsFound)
954 {
955 etaclust[etaindex].fStartPad = pad;
956 etaclust[etaindex].fEndPad = pad;
957 etaclust[etaindex].fIsFound = 1;
958 continue;
959 }
960 else
961 {
962 if(pad <= (etaclust[etaindex].fEndPad + 1))
963 {
964 etaclust[etaindex].fEndPad = pad;
965 }
966 else
967 {
968 FillCluster(i,etaindex,etaclust,ilastpatch,firstbinx,lastbinx,nbinx,firstbiny);
969
970 etaclust[etaindex].fStartPad = pad;
971 etaclust[etaindex].fEndPad = pad;
972
973 }
974 }
975 lastetaindex = etaindex;
976 }
977
978 //Write remaining clusters
979 for(Int_t etaindex = 0;etaindex < netasegments;etaindex++)
980 {
981 //Check for empty row
982 if((etaclust[etaindex].fStartPad == 0) && (etaclust[etaindex].fEndPad == 0)) continue;
983
984 FillCluster(i,etaindex,etaclust,ilastpatch,firstbinx,lastbinx,nbinx,firstbiny);
0bd0c1ef 985 }
986
a8ffd46b 987 delete [] etaclust;
988}
989
a8ffd46b 990#ifndef do_mc
4aa41877 991Int_t AliHLTHoughTransformerRow::GetTrackID(Int_t /*etaindex*/,Double_t /*alpha1*/,Double_t /*alpha2*/) const
3cb5c013 992{
993 // Does nothing if do_mc undefined
4aa41877 994 LOG(AliHLTLog::kWarning,"AliHLTHoughTransformerRow::GetTrackID","Data")
a8ffd46b 995 <<"Flag switched off"<<ENDLOG;
996 return -1;
997#else
4aa41877 998Int_t AliHLTHoughTransformerRow::GetTrackID(Int_t etaindex,Double_t alpha1,Double_t alpha2) const
3cb5c013 999{
1000 // Returns the MC label for a given peak found in the Hough space
917e711b 1001 if(etaindex < 0 || etaindex > GetNEtaSegments())
0bd0c1ef 1002 {
4aa41877 1003 LOG(AliHLTLog::kWarning,"AliHLTHoughTransformerRow::GetTrackID","Data")
917e711b 1004 <<"Wrong etaindex "<<etaindex<<ENDLOG;
0bd0c1ef 1005 return -1;
1006 }
4aa41877 1007 AliHLTHistogram *hist = fParamSpace[etaindex];
e81040b7 1008 Int_t bin = hist->FindLabelBin(alpha1,alpha2);
e81ffe36 1009 if(bin==-1) {
4aa41877 1010 LOG(AliHLTLog::kWarning,"AliHLTHoughTransformerRow::GetTrackID()","")
e81040b7 1011 <<"Track candidate outside Hough space boundaries: "<<alpha1<<" "<<alpha2<<ENDLOG;
e81ffe36 1012 return -1;
1013 }
0bd0c1ef 1014 Int_t label=-1;
1015 Int_t max=0;
1016 for(UInt_t i=0; i<(MaxTrack-1); i++)
1017 {
a8ffd46b 1018 Int_t nhits=fTrackID[etaindex][bin].fNHits[i];
0bd0c1ef 1019 if(nhits == 0) break;
1020 if(nhits > max)
1021 {
1022 max = nhits;
a8ffd46b 1023 label = fTrackID[etaindex][bin].fLabel[i];
0bd0c1ef 1024 }
1025 }
1026 Int_t label2=-1;
1027 Int_t max2=0;
1028 for(UInt_t i=0; i<(MaxTrack-1); i++)
1029 {
a8ffd46b 1030 Int_t nhits=fTrackID[etaindex][bin].fNHits[i];
0bd0c1ef 1031 if(nhits == 0) break;
1032 if(nhits > max2)
1033 {
a8ffd46b 1034 if(fTrackID[etaindex][bin].fLabel[i]!=label) {
0bd0c1ef 1035 max2 = nhits;
a8ffd46b 1036 label2 = fTrackID[etaindex][bin].fLabel[i];
0bd0c1ef 1037 }
1038 }
1039 }
917e711b 1040 if(max2 !=0 ) {
4aa41877 1041 LOG(AliHLTLog::kDebug,"AliHLTHoughTransformerRow::GetTrackID()","")
a8ffd46b 1042 <<" TrackID"<<" label "<<label<<" max "<<max<<" label2 "<<label2<<" max2 "<<max2<<" "<<(Float_t)max2/(Float_t)max<<" "<<fTrackID[etaindex][bin].fLabel[MaxTrack-1]<<" "<<(Int_t)fTrackID[etaindex][bin].fNHits[MaxTrack-1]<<ENDLOG;
917e711b 1043 }
0bd0c1ef 1044 return label;
1045#endif
0bd0c1ef 1046}
1047
4aa41877 1048Int_t AliHLTHoughTransformerRow::GetTrackLength(Double_t alpha1,Double_t alpha2,Int_t *rows) const
f644512a 1049{
1050 // Returns the track length for a given peak found in the Hough space
1051
4aa41877 1052 AliHLTHistogram *hist = fParamSpace[0];
f644512a 1053 Int_t bin = hist->FindBin(alpha1,alpha2);
1054 if(bin==-1) {
4aa41877 1055 LOG(AliHLTLog::kWarning,"AliHLTHoughTransformerRow::GetTrackLength()","")
f644512a 1056 <<"Track candidate outside Hough space boundaries: "<<alpha1<<" "<<alpha2<<ENDLOG;
1057 return -1;
1058 }
1059 rows[0] = fTrackFirstRow[bin];
1060 rows[1] = fTrackLastRow[bin];
1061
1062 return 0;
1063}
1064
4aa41877 1065inline void AliHLTHoughTransformerRow::FillClusterRow(UChar_t i,Int_t binx1,Int_t binx2,UChar_t *ngaps2,UChar_t *currentrow2,UChar_t *lastrow2
a8ffd46b 1066#ifdef do_mc
4aa41877 1067 ,AliHLTEtaRow etaclust,AliHLTTrackIndex *trackid
a8ffd46b 1068#endif
1069 )
1070{
ce4c28d5 1071 // The method is a part of the fast hough transform.
1072 // It fills one row of the hough space.
1073 // It is called by FillCluster() method inside the
1074 // loop over alpha2 bins
1075
a8ffd46b 1076 for(Int_t bin=binx1;bin<=binx2;bin++)
1077 {
1078 if(ngaps2[bin] < MAX_N_GAPS) {
1079 if(i < lastrow2[bin] && i > currentrow2[bin]) {
1080 ngaps2[bin] += (i-currentrow2[bin]-1);
1081 currentrow2[bin]=i;
1082 }
1083#ifdef do_mc
1084 if(i < lastrow2[bin] && i >= currentrow2[bin]) {
1085 for(UInt_t t=0;t<(MaxTrack-1); t++)
1086 {
1087 Int_t label = etaclust.fMcLabels[t];
1088 if(label == 0) break;
1089 UInt_t c;
1090 Int_t tempbin2 = (Int_t)(bin/2);
1091 for(c=0; c<(MaxTrack-1); c++)
1092 if(trackid[tempbin2].fLabel[c] == label || trackid[tempbin2].fNHits[c] == 0)
1093 break;
1094 trackid[tempbin2].fLabel[c] = label;
1095 if(trackid[tempbin2].fCurrentRow[c] != i) {
1096 trackid[tempbin2].fNHits[c]++;
1097 trackid[tempbin2].fCurrentRow[c] = i;
1098 }
1099 }
1100 }
1101#endif
1102 }
1103 }
1104
1105}
1106
4aa41877 1107inline void AliHLTHoughTransformerRow::FillCluster(UChar_t i,Int_t etaindex,AliHLTEtaRow *etaclust,Int_t ilastpatch,Int_t firstbinx,Int_t lastbinx,Int_t nbinx,Int_t firstbiny)
a8ffd46b 1108{
ce4c28d5 1109 // The method is a part of the fast hough transform.
1110 // It fills a TPC cluster into the hough space.
1111
a8ffd46b 1112 UChar_t *ngaps = fGapCount[etaindex];
1113 UChar_t *currentrow = fCurrentRowCount[etaindex];
1114 UChar_t *lastrow = fTrackLastRow;
1115 UChar_t *prevbin = fPrevBin[etaindex];
1116 UChar_t *nextbin = fNextBin[etaindex];
1117 UChar_t *nextrow = fNextRow[etaindex];
1118#ifdef do_mc
4aa41877 1119 AliHLTTrackIndex *trackid = fTrackID[etaindex];
a8ffd46b 1120#endif
1121
1122 //Do the transformation:
4aa41877 1123 AliHLTPadHoughParams *startparams = &fStartPadParams[(Int_t)i][etaclust[etaindex].fStartPad];
1124 AliHLTPadHoughParams *endparams = &fEndPadParams[(Int_t)i][etaclust[etaindex].fEndPad];
a8ffd46b 1125
1126 Float_t alpha1 = startparams->fAlpha;
1127 Float_t deltaalpha1 = startparams->fDeltaAlpha;
1128 Float_t alpha2 = endparams->fAlpha;
1129 Float_t deltaalpha2 = endparams->fDeltaAlpha;
1130
1131 Int_t firstbin1 = startparams->fFirstBin;
1132 Int_t firstbin2 = endparams->fFirstBin;
1133 Int_t firstbin = firstbin1;
1134 if(firstbin>firstbin2) firstbin = firstbin2;
1135
1136 Int_t lastbin1 = startparams->fLastBin;
1137 Int_t lastbin2 = endparams->fLastBin;
1138 Int_t lastbin = lastbin1;
1139 if(lastbin<lastbin2) lastbin = lastbin2;
1140
1141 alpha1 += (firstbin-firstbiny)*deltaalpha1;
1142 alpha2 += (firstbin-firstbiny)*deltaalpha2;
1143
1144 //Fill the histogram along the alpha2 range
1145 if(ilastpatch == -1) {
1146 for(Int_t b=firstbin; b<=lastbin; b++, alpha1 += deltaalpha1, alpha2 += deltaalpha2)
1147 {
1148 Int_t binx1 = 1 + (Int_t)alpha1;
1149 if(binx1>lastbinx) continue;
1150 if(binx1<firstbinx) binx1 = firstbinx;
1151 Int_t binx2 = 1 + (Int_t)alpha2;
1152 if(binx2<firstbinx) continue;
1153 if(binx2>lastbinx) binx2 = lastbinx;
1154#ifdef do_mc
1155 if(binx2<binx1) {
4aa41877 1156 LOG(AliHLTLog::kWarning,"AliHLTHoughTransformerRow::TransformCircle()","")
a8ffd46b 1157 <<"Wrong filling "<<binx1<<" "<<binx2<<" "<<i<<" "<<etaclust[etaindex].fStartPad<<" "<<etaclust[etaindex].fEndPad<<ENDLOG;
1158 }
1159#endif
1160 Int_t tempbin = b*nbinx;
1161 UChar_t *ngaps2 = ngaps + tempbin;
1162 UChar_t *currentrow2 = currentrow + tempbin;
1163 UChar_t *lastrow2 = lastrow + tempbin;
1164#ifdef do_mc
1165 Int_t tempbin2 = ((Int_t)(b/2))*((Int_t)((nbinx+1)/2));
4aa41877 1166 AliHLTTrackIndex *trackid2 = trackid + tempbin2;
a8ffd46b 1167#endif
1168 FillClusterRow(i,binx1,binx2,ngaps2,currentrow2,lastrow2
1169#ifdef do_mc
1170 ,etaclust[etaindex],trackid2
1171#endif
1172 );
1173 }
1174 }
1175 else {
1176 for(Int_t b=firstbin; b<=lastbin; b++, alpha1 += deltaalpha1, alpha2 += deltaalpha2)
1177 {
1178 Int_t binx1 = 1 + (Int_t)alpha1;
1179 if(binx1>lastbinx) continue;
1180 if(binx1<firstbinx) binx1 = firstbinx;
1181 Int_t binx2 = 1 + (Int_t)alpha2;
1182 if(binx2<firstbinx) continue;
1183 if(binx2>lastbinx) binx2 = lastbinx;
1184#ifdef do_mc
1185 if(binx2<binx1) {
4aa41877 1186 LOG(AliHLTLog::kWarning,"AliHLTHoughTransformerRow::TransformCircle()","")
a8ffd46b 1187 <<"Wrong filling "<<binx1<<" "<<binx2<<" "<<i<<" "<<etaclust[etaindex].fStartPad<<" "<<etaclust[etaindex].fEndPad<<ENDLOG;
1188 }
1189#endif
fbd52b13 1190 if(nextrow[b] > b) {
a8ffd46b 1191 Int_t deltab = (nextrow[b] - b - 1);
1192 b += deltab;
1193 alpha1 += deltaalpha1*deltab;
1194 alpha2 += deltaalpha2*deltab;
1195 continue;
1196 }
1197 Int_t tempbin = b*nbinx;
1198 binx1 = (UInt_t)nextbin[binx1 + tempbin];
1199 binx2 = (UInt_t)prevbin[binx2 + tempbin];
1200 if(binx2<binx1) continue;
1201 UChar_t *ngaps2 = ngaps + tempbin;
1202 UChar_t *currentrow2 = currentrow + tempbin;
1203 UChar_t *lastrow2 = lastrow + tempbin;
1204#ifdef do_mc
1205 Int_t tempbin2 = ((Int_t)(b/2))*((Int_t)((nbinx+1)/2));
4aa41877 1206 AliHLTTrackIndex *trackid2 = trackid + tempbin2;
a8ffd46b 1207#endif
1208 FillClusterRow(i,binx1,binx2,ngaps2,currentrow2,lastrow2
1209#ifdef do_mc
1210 ,etaclust[etaindex],trackid2
1211#endif
1212 );
1213 }
1214 }
1215}
1216
1217#ifdef do_mc
4aa41877 1218inline void AliHLTHoughTransformerRow::FillClusterMCLabels(AliHLTDigitData digpt,AliHLTEtaRow *etaclust)
a8ffd46b 1219{
ce4c28d5 1220 // The method is a part of the fast hough transform.
1221 // It fills the MC labels of a TPC cluster into a
1222 // special hough space array.
a8ffd46b 1223 for(Int_t t=0; t<3; t++)
1224 {
1225 Int_t label = digpt.fTrackID[t];
1226 if(label < 0) break;
1227 UInt_t c;
1228 for(c=0; c<(MaxTrack-1); c++)
8cc2b1a5 1229 if(etaclust->fMcLabels[c] == label || etaclust->fMcLabels[c] == 0)
a8ffd46b 1230 break;
1231
8cc2b1a5 1232 etaclust->fMcLabels[c] = label;
a8ffd46b 1233 }
de3c3890 1234}
a8ffd46b 1235#endif
ce4c28d5 1236
4aa41877 1237void AliHLTHoughTransformerRow::SetTransformerArrays(AliHLTHoughTransformerRow *tr)
ce4c28d5 1238{
1239 // In case of sequential filling of the hough space, the method is used to
1240 // transmit the pointers to the hough arrays from one transformer to the
1241 // following one.
1242
1243 fGapCount = tr->fGapCount;
1244 fCurrentRowCount = tr->fCurrentRowCount;
1245#ifdef do_mc
1246 fTrackID = tr->fTrackID;
1247#endif
1248 fTrackNRows = tr->fTrackNRows;
1249 fTrackFirstRow = tr->fTrackFirstRow;
1250 fTrackLastRow = tr->fTrackLastRow;
1251 fInitialGapCount = tr->fInitialGapCount;
1252
1253 fPrevBin = tr->fPrevBin;
1254 fNextBin = tr->fNextBin;
1255 fNextRow = tr->fNextRow;
1256
1257 fStartPadParams = tr->fStartPadParams;
1258 fEndPadParams = tr->fEndPadParams;
f644512a 1259 fLUTr = tr->fLUTr;
ce4c28d5 1260 fLUTforwardZ = tr->fLUTforwardZ;
ce4c28d5 1261 fLUTbackwardZ = tr->fLUTbackwardZ;
ce4c28d5 1262
7e4eb453 1263 fParamSpace = tr->fParamSpace;
1264
ce4c28d5 1265 return;
1266}