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