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