]>
Commit | Line | Data |
---|---|---|
1 | // @(#) $Id$ | |
2 | ||
3 | // Author: Anders Vestbo <mailto:vestbo@fi.uib.no> | |
4 | //*-- Copyright © ALICE HLT Group | |
5 | ||
6 | #include "AliL3StandardIncludes.h" | |
7 | #include <sys/time.h> | |
8 | ||
9 | #include "AliL3Logging.h" | |
10 | #include "AliL3HoughMerger.h" | |
11 | #include "AliL3HoughIntMerger.h" | |
12 | #include "AliL3HoughGlobalMerger.h" | |
13 | #include "AliL3Histogram.h" | |
14 | #include "AliL3Hough.h" | |
15 | #include "AliL3HoughTransformer.h" | |
16 | #include "AliL3HoughClusterTransformer.h" | |
17 | #include "AliL3HoughTransformerLUT.h" | |
18 | #include "AliL3HoughTransformerVhdl.h" | |
19 | #include "AliL3HoughTransformerRow.h" | |
20 | #include "AliL3HoughMaxFinder.h" | |
21 | #include "AliL3Benchmark.h" | |
22 | #ifdef use_aliroot | |
23 | #include "AliL3FileHandler.h" | |
24 | #else | |
25 | #include "AliL3MemHandler.h" | |
26 | #endif | |
27 | #include "AliL3DataHandler.h" | |
28 | #include "AliL3DigitData.h" | |
29 | #include "AliL3HoughEval.h" | |
30 | #include "AliL3Transform.h" | |
31 | #include "AliL3TrackArray.h" | |
32 | #include "AliL3HoughTrack.h" | |
33 | #include "AliL3DDLDataFileHandler.h" | |
34 | #include "AliL3HoughKalmanTrack.h" | |
35 | ||
36 | #include "TThread.h" | |
37 | ||
38 | #if __GNUC__ == 3 | |
39 | using namespace std; | |
40 | #endif | |
41 | ||
42 | /** /class AliL3Hough | |
43 | //<pre> | |
44 | //_____________________________________________________________ | |
45 | // AliL3Hough | |
46 | // | |
47 | // Interface class for the Hough transform | |
48 | // | |
49 | // Example how to use: | |
50 | // | |
51 | // AliL3Hough *hough = new AliL3Hough(path,kTRUE,NumberOfEtaSegments); | |
52 | // hough->ReadData(slice); | |
53 | // hough->Transform(); | |
54 | // hough->FindTrackCandidates(); | |
55 | // | |
56 | // AliL3TrackArray *tracks = hough->GetTracks(patch); | |
57 | // | |
58 | //</pre> | |
59 | */ | |
60 | ||
61 | ClassImp(AliL3Hough) | |
62 | ||
63 | AliL3Hough::AliL3Hough() | |
64 | { | |
65 | //Constructor | |
66 | ||
67 | fBinary = kFALSE; | |
68 | fAddHistograms = kFALSE; | |
69 | fDoIterative = kFALSE; | |
70 | fWriteDigits = kFALSE; | |
71 | fUse8bits = kFALSE; | |
72 | ||
73 | fMemHandler = 0; | |
74 | fHoughTransformer = 0; | |
75 | fEval = 0; | |
76 | fPeakFinder = 0; | |
77 | fTracks = 0; | |
78 | fGlobalTracks = 0; | |
79 | fMerger = 0; | |
80 | fInterMerger = 0; | |
81 | fGlobalMerger = 0; | |
82 | fBenchmark = 0; | |
83 | ||
84 | fNEtaSegments = 0; | |
85 | fNPatches = 0; | |
86 | fLastPatch =-1; | |
87 | fVersion = 0; | |
88 | fCurrentSlice = 0; | |
89 | fEvent = 0; | |
90 | ||
91 | fKappaSpread = 6; | |
92 | fPeakRatio = 0.5; | |
93 | fInputFile = 0; | |
94 | fInputPtr = 0; | |
95 | fRawEvent = 0; | |
96 | ||
97 | SetTransformerParams(); | |
98 | SetThreshold(); | |
99 | SetNSaveIterations(); | |
100 | SetPeakThreshold(); | |
101 | #ifdef use_aliroot | |
102 | //just be sure that index is empty for new event | |
103 | AliL3FileHandler::CleanStaticIndex(); | |
104 | #ifdef use_newio | |
105 | fRunLoader = 0; | |
106 | #endif | |
107 | #endif | |
108 | fThread = 0; | |
109 | } | |
110 | ||
111 | AliL3Hough::AliL3Hough(Char_t *path,Bool_t binary,Int_t netasegments,Bool_t bit8,Int_t tv,Char_t *infile,Char_t *ptr) | |
112 | { | |
113 | //Normal constructor | |
114 | fBinary = binary; | |
115 | strcpy(fPath,path); | |
116 | fNEtaSegments = netasegments; | |
117 | fAddHistograms = kFALSE; | |
118 | fDoIterative = kFALSE; | |
119 | fWriteDigits = kFALSE; | |
120 | fUse8bits = bit8; | |
121 | fVersion = tv; | |
122 | fKappaSpread=6; | |
123 | fPeakRatio=0.5; | |
124 | if(!fBinary) { | |
125 | if(infile) { | |
126 | fInputFile = infile; | |
127 | fInputPtr = 0; | |
128 | } | |
129 | else { | |
130 | fInputFile = 0; | |
131 | fInputPtr = ptr; | |
132 | } | |
133 | } | |
134 | else { | |
135 | fInputFile = 0; | |
136 | fInputPtr = 0; | |
137 | } | |
138 | #ifdef use_aliroot | |
139 | //just be sure that index is empty for new event | |
140 | AliL3FileHandler::CleanStaticIndex(); | |
141 | #ifdef use_newio | |
142 | fRunLoader = 0; | |
143 | #endif | |
144 | #endif | |
145 | fThread = 0; | |
146 | } | |
147 | ||
148 | AliL3Hough::~AliL3Hough() | |
149 | { | |
150 | //dtor | |
151 | ||
152 | CleanUp(); | |
153 | if(fMerger) | |
154 | delete fMerger; | |
155 | //cout << "Cleaned class merger " << endl; | |
156 | if(fInterMerger) | |
157 | delete fInterMerger; | |
158 | //cout << "Cleaned class inter " << endl; | |
159 | if(fPeakFinder) | |
160 | delete fPeakFinder; | |
161 | //cout << "Cleaned class peak " << endl; | |
162 | if(fGlobalMerger) | |
163 | delete fGlobalMerger; | |
164 | //cout << "Cleaned class global " << endl; | |
165 | if(fBenchmark) | |
166 | delete fBenchmark; | |
167 | //cout << "Cleaned class bench " << endl; | |
168 | if(fGlobalTracks) | |
169 | delete fGlobalTracks; | |
170 | //cout << "Cleaned class globaltracks " << endl; | |
171 | if(fThread) { | |
172 | // fThread->Delete(); | |
173 | delete fThread; | |
174 | fThread = 0; | |
175 | } | |
176 | } | |
177 | ||
178 | void AliL3Hough::CleanUp() | |
179 | { | |
180 | //Cleanup memory | |
181 | ||
182 | for(Int_t i=0; i<fNPatches; i++) | |
183 | { | |
184 | if(fTracks[i]) delete fTracks[i]; | |
185 | //cout << "Cleaned tracks " << i << endl; | |
186 | if(fEval[i]) delete fEval[i]; | |
187 | //cout << "Cleaned eval " << i << endl; | |
188 | if(fHoughTransformer[i]) delete fHoughTransformer[i]; | |
189 | //cout << "Cleaned traf " << i << endl; | |
190 | if(fMemHandler[i]) delete fMemHandler[i]; | |
191 | //cout << "Cleaned mem " << i << endl; | |
192 | } | |
193 | ||
194 | if(fTracks) delete [] fTracks; | |
195 | //cout << "Cleaned class tracks " << endl; | |
196 | if(fEval) delete [] fEval; | |
197 | //cout << "Cleaned class eval " << endl; | |
198 | if(fHoughTransformer) delete [] fHoughTransformer; | |
199 | //cout << "Cleaned cleass trafo " << endl; | |
200 | if(fMemHandler) delete [] fMemHandler; | |
201 | //cout << "Cleaned class mem " << endl; | |
202 | } | |
203 | ||
204 | void AliL3Hough::Init(Int_t netasegments,Int_t tv,AliRawEvent *rawevent,Float_t zvertex) | |
205 | { | |
206 | //Normal constructor | |
207 | fNEtaSegments = netasegments; | |
208 | fVersion = tv; | |
209 | fRawEvent = rawevent; | |
210 | fZVertex = zvertex; | |
211 | ||
212 | Init(); | |
213 | } | |
214 | ||
215 | void AliL3Hough::Init(Char_t *path,Bool_t binary,Int_t netasegments,Bool_t bit8,Int_t tv,Char_t *infile,Char_t *ptr,Float_t zvertex) | |
216 | { | |
217 | //Normal init of the AliL3Hough | |
218 | fBinary = binary; | |
219 | strcpy(fPath,path); | |
220 | fNEtaSegments = netasegments; | |
221 | fWriteDigits = kFALSE; | |
222 | fUse8bits = bit8; | |
223 | fVersion = tv; | |
224 | if(!fBinary) { | |
225 | if(infile) { | |
226 | fInputFile = infile; | |
227 | fInputPtr = 0; | |
228 | } | |
229 | else { | |
230 | fInputFile = 0; | |
231 | fInputPtr = ptr; | |
232 | } | |
233 | } | |
234 | else { | |
235 | fInputFile = 0; | |
236 | fInputPtr = 0; | |
237 | } | |
238 | fZVertex = zvertex; | |
239 | ||
240 | Init(); //do the rest | |
241 | } | |
242 | ||
243 | void AliL3Hough::Init(Bool_t doit, Bool_t addhists) | |
244 | { | |
245 | // Init | |
246 | fDoIterative = doit; | |
247 | fAddHistograms = addhists; | |
248 | ||
249 | fNPatches = AliL3Transform::GetNPatches(); | |
250 | fHoughTransformer = new AliL3HoughBaseTransformer*[fNPatches]; | |
251 | fMemHandler = new AliL3MemHandler*[fNPatches]; | |
252 | ||
253 | fTracks = new AliL3TrackArray*[fNPatches]; | |
254 | fEval = new AliL3HoughEval*[fNPatches]; | |
255 | ||
256 | fGlobalTracks = new AliL3TrackArray("AliL3HoughTrack"); | |
257 | ||
258 | AliL3HoughBaseTransformer *lasttransformer = 0; | |
259 | ||
260 | for(Int_t i=0; i<fNPatches; i++) | |
261 | { | |
262 | switch (fVersion){ //choose Transformer | |
263 | case 1: | |
264 | fHoughTransformer[i] = new AliL3HoughTransformerLUT(0,i,fNEtaSegments); | |
265 | break; | |
266 | case 2: | |
267 | fHoughTransformer[i] = new AliL3HoughClusterTransformer(0,i,fNEtaSegments); | |
268 | break; | |
269 | case 3: | |
270 | fHoughTransformer[i] = new AliL3HoughTransformerVhdl(0,i,fNEtaSegments,fNSaveIterations); | |
271 | break; | |
272 | case 4: | |
273 | fHoughTransformer[i] = new AliL3HoughTransformerRow(0,i,fNEtaSegments,kFALSE,fZVertex); | |
274 | break; | |
275 | default: | |
276 | fHoughTransformer[i] = new AliL3HoughTransformer(0,i,fNEtaSegments,kFALSE,kFALSE); | |
277 | } | |
278 | ||
279 | fHoughTransformer[i]->SetLastTransformer(lasttransformer); | |
280 | lasttransformer = fHoughTransformer[i]; | |
281 | // fHoughTransformer[i]->CreateHistograms(fNBinX[i],fLowPt[i],fNBinY[i],-fPhi[i],fPhi[i]); | |
282 | fHoughTransformer[i]->CreateHistograms(fNBinX[i],-fLowPt[i],fLowPt[i],fNBinY[i],-fPhi[i],fPhi[i]); | |
283 | //fHoughTransformer[i]->CreateHistograms(fLowPt[i],fUpperPt[i],fPtRes[i],fNBinY[i],fPhi[i]); | |
284 | ||
285 | fHoughTransformer[i]->SetLowerThreshold(fThreshold[i]); | |
286 | fHoughTransformer[i]->SetUpperThreshold(100); | |
287 | ||
288 | LOG(AliL3Log::kInformational,"AliL3Hough::Init","Version") | |
289 | <<"Initializing Hough transformer version "<<fVersion<<ENDLOG; | |
290 | ||
291 | fEval[i] = new AliL3HoughEval(); | |
292 | fTracks[i] = new AliL3TrackArray("AliL3HoughTrack"); | |
293 | if(fUse8bits) | |
294 | fMemHandler[i] = new AliL3DataHandler(); | |
295 | else | |
296 | #ifdef use_aliroot | |
297 | { | |
298 | if(!fRawEvent) { | |
299 | if(!fInputFile) { | |
300 | if(!fInputPtr) { | |
301 | /* In case of reading digits file */ | |
302 | fMemHandler[i] = new AliL3FileHandler(kTRUE); //use static index | |
303 | if(!fBinary) { | |
304 | #if use_newio | |
305 | if(!fRunLoader) { | |
306 | #endif | |
307 | Char_t filename[1024]; | |
308 | sprintf(filename,"%s/digitfile.root",fPath); | |
309 | fMemHandler[i]->SetAliInput(filename); | |
310 | #if use_newio | |
311 | } | |
312 | else { | |
313 | fMemHandler[i]->SetAliInput(fRunLoader); | |
314 | } | |
315 | #endif | |
316 | } | |
317 | } | |
318 | else { | |
319 | /* In case of reading from DATE */ | |
320 | fMemHandler[i] = new AliL3DDLDataFileHandler(); | |
321 | fMemHandler[i]->SetReaderInput(fInputPtr,-1); | |
322 | } | |
323 | } | |
324 | else { | |
325 | /* In case of reading rawdata from ROOT file */ | |
326 | fMemHandler[i] = new AliL3DDLDataFileHandler(); | |
327 | fMemHandler[i]->SetReaderInput(fInputFile); | |
328 | } | |
329 | } | |
330 | else { | |
331 | /* In case of reading rawdata using AliRawEvent */ | |
332 | fMemHandler[i] = new AliL3DDLDataFileHandler(); | |
333 | fMemHandler[i]->SetReaderInput(fRawEvent); | |
334 | } | |
335 | } | |
336 | #else | |
337 | fMemHandler[i] = new AliL3MemHandler(); | |
338 | #endif | |
339 | } | |
340 | ||
341 | fPeakFinder = new AliL3HoughMaxFinder("KappaPhi",50000); | |
342 | if(fVersion!=4) { | |
343 | fMerger = new AliL3HoughMerger(fNPatches); | |
344 | fInterMerger = new AliL3HoughIntMerger(); | |
345 | } | |
346 | else { | |
347 | fMerger = 0; | |
348 | fInterMerger = 0; | |
349 | } | |
350 | fGlobalMerger = 0; | |
351 | fBenchmark = new AliL3Benchmark(); | |
352 | } | |
353 | ||
354 | void AliL3Hough::SetTransformerParams(Float_t ptres,Float_t ptmin,Float_t ptmax,Int_t ny,Int_t patch) | |
355 | { | |
356 | // Setup the parameters for the Hough Transformer | |
357 | // This includes the bin size and limits for | |
358 | // the parameter space histograms | |
359 | ||
360 | Int_t mrow; | |
361 | Float_t psi=0; | |
362 | if(patch==-1) | |
363 | mrow = 80; | |
364 | else | |
365 | mrow = AliL3Transform::GetLastRow(patch); | |
366 | if(ptmin) | |
367 | { | |
368 | Double_t lineradius = sqrt(pow(AliL3Transform::Row2X(mrow),2) + pow(AliL3Transform::GetMaxY(mrow),2)); | |
369 | Double_t kappa = -1*AliL3Transform::GetBField()*AliL3Transform::GetBFact()/ptmin; | |
370 | psi = AliL3Transform::Deg2Rad(10) - asin(lineradius*kappa/2); | |
371 | cout<<"Calculated psi range "<<psi<<" in patch "<<patch<<endl; | |
372 | } | |
373 | ||
374 | if(patch==-1) | |
375 | { | |
376 | Int_t i=0; | |
377 | while(i < 6) | |
378 | { | |
379 | fPtRes[i] = ptres; | |
380 | fLowPt[i] = ptmin; | |
381 | fUpperPt[i] = ptmax; | |
382 | fNBinY[i] = ny; | |
383 | fPhi[i] = psi; | |
384 | fNBinX[i]=0; | |
385 | i++; | |
386 | } | |
387 | return; | |
388 | } | |
389 | ||
390 | fPtRes[patch] = ptres; | |
391 | fLowPt[patch] = ptmin; | |
392 | fUpperPt[patch] = ptmax; | |
393 | fNBinY[patch] = ny; | |
394 | fPhi[patch] = psi; | |
395 | } | |
396 | /* | |
397 | void AliL3Hough::SetTransformerParams(Int_t nx,Int_t ny,Float_t ptmin,Int_t patch) | |
398 | { | |
399 | // Setup the parameters for the Hough Transformer | |
400 | ||
401 | Int_t mrow=80; | |
402 | Double_t lineradius = sqrt(pow(AliL3Transform::Row2X(mrow),2) + pow(AliL3Transform::GetMaxY(mrow),2)); | |
403 | Double_t kappa = -1*AliL3Transform::GetBField()*AliL3Transform::GetBFact()/ptmin; | |
404 | Double_t psi = AliL3Transform::Deg2Rad(10) - asin(lineradius*kappa/2); | |
405 | cout<<"Calculated psi range "<<psi<<" in patch "<<patch<<endl; | |
406 | ||
407 | Int_t i=0; | |
408 | while(i < 6) | |
409 | { | |
410 | fLowPt[i] = ptmin; | |
411 | fNBinY[i] = ny; | |
412 | fNBinX[i] = nx; | |
413 | fPhi[i] = psi; | |
414 | i++; | |
415 | } | |
416 | } | |
417 | */ | |
418 | void AliL3Hough::SetTransformerParams(Int_t nx,Int_t ny,Float_t ptmin,Int_t /*patch*/) | |
419 | { | |
420 | // Setup the parameters for the Hough Transformer | |
421 | ||
422 | Double_t lineradius = 1.0/(AliL3HoughTransformerRow::GetBeta1()*sqrt(1.0+tan(AliL3Transform::Pi()*10/180)*tan(AliL3Transform::Pi()*10/180))); | |
423 | Double_t alpha1 = AliL3HoughTransformerRow::GetBeta1()*tan(AliL3Transform::Pi()*10/180); | |
424 | Double_t kappa = 1*AliL3Transform::GetBField()*AliL3Transform::GetBFact()/(ptmin*0.9); | |
425 | Double_t psi = AliL3Transform::Deg2Rad(10) - asin(lineradius*kappa/2); | |
426 | // cout<<"Calculated psi range "<<psi<<" in patch "<<patch<<endl; | |
427 | Double_t alpha2 = alpha1 - (AliL3HoughTransformerRow::GetBeta1()-AliL3HoughTransformerRow::GetBeta2())*tan(psi); | |
428 | // cout<<"Calculated alphas range "<<alpha1<<" "<<alpha2<<" in patch "<<patch<<endl; | |
429 | ||
430 | Int_t i=0; | |
431 | while(i < 6) | |
432 | { | |
433 | fLowPt[i] = 1.1*alpha1; | |
434 | fNBinY[i] = ny; | |
435 | fNBinX[i] = nx; | |
436 | fPhi[i] = alpha2; | |
437 | i++; | |
438 | } | |
439 | } | |
440 | ||
441 | void AliL3Hough::CalcTransformerParams(Float_t ptmin) | |
442 | { | |
443 | // Setup the parameters for the Row Hough Transformer | |
444 | // Automatically adjusts the number of bins in X and Y in a way | |
445 | // that the size of the hough bin is 2x (in X) and 2.5 (in Y) the | |
446 | // size of the tpc pads | |
447 | ||
448 | Double_t lineradius = 1.0/(AliL3HoughTransformerRow::GetBeta1()*sqrt(1.0+tan(AliL3Transform::Pi()*10/180)*tan(AliL3Transform::Pi()*10/180))); | |
449 | Double_t alpha1 = AliL3HoughTransformerRow::GetBeta1()*tan(AliL3Transform::Pi()*10/180); | |
450 | Double_t kappa = 1*AliL3Transform::GetBField()*AliL3Transform::GetBFact()/(ptmin*0.9); | |
451 | Double_t psi = AliL3Transform::Deg2Rad(10) - asin(lineradius*kappa/2); | |
452 | // cout<<"Calculated psi range "<<psi<<endl; | |
453 | Double_t alpha2 = alpha1 - (AliL3HoughTransformerRow::GetBeta1()-AliL3HoughTransformerRow::GetBeta2())*tan(psi); | |
454 | alpha1 *= 1.1; | |
455 | // cout<<"Calculated alphas range "<<alpha1<<" "<<alpha2<<endl; | |
456 | ||
457 | Double_t sizex = 2.0*AliL3Transform::GetPadPitchWidthLow()*AliL3HoughTransformerRow::GetBeta1()*AliL3HoughTransformerRow::GetBeta1(); | |
458 | Double_t sizey = 2.5*AliL3Transform::GetPadPitchWidthUp()*AliL3HoughTransformerRow::GetBeta2()*AliL3HoughTransformerRow::GetBeta2(); | |
459 | ||
460 | Int_t nx = 2*(Int_t)(alpha1/sizex)+1; | |
461 | Int_t ny = 2*(Int_t)(alpha2/sizey)+1; | |
462 | // cout<<"Calculated number of bins "<<nx<<" "<<ny<<endl; | |
463 | ||
464 | Int_t i=0; | |
465 | while(i < 6) | |
466 | { | |
467 | fLowPt[i] = alpha1; | |
468 | fNBinY[i] = ny; | |
469 | fNBinX[i] = nx; | |
470 | fPhi[i] = alpha2; | |
471 | i++; | |
472 | } | |
473 | } | |
474 | ||
475 | void AliL3Hough::SetTransformerParams(Int_t nx,Int_t ny,Float_t lpt,Float_t phi) | |
476 | { | |
477 | Int_t i=0; | |
478 | while(i < 6) | |
479 | { | |
480 | fLowPt[i] = lpt; | |
481 | fNBinY[i] = ny; | |
482 | fNBinX[i] = nx; | |
483 | fPhi[i] = phi; | |
484 | i++; | |
485 | } | |
486 | } | |
487 | ||
488 | void AliL3Hough::SetThreshold(Int_t t3,Int_t patch) | |
489 | { | |
490 | // Set digits threshold | |
491 | if(patch==-1) | |
492 | { | |
493 | Int_t i=0; | |
494 | while(i < 6) | |
495 | fThreshold[i++]=t3; | |
496 | return; | |
497 | } | |
498 | fThreshold[patch]=t3; | |
499 | } | |
500 | ||
501 | void AliL3Hough::SetPeakThreshold(Int_t threshold,Int_t patch) | |
502 | { | |
503 | // Set Peak Finder threshold | |
504 | if(patch==-1) | |
505 | { | |
506 | Int_t i=0; | |
507 | while(i < 6) | |
508 | fPeakThreshold[i++]=threshold; | |
509 | return; | |
510 | } | |
511 | fPeakThreshold[patch]=threshold; | |
512 | } | |
513 | ||
514 | void AliL3Hough::DoBench(Char_t *name) | |
515 | { | |
516 | fBenchmark->Analyze(name); | |
517 | } | |
518 | ||
519 | void AliL3Hough::Process(Int_t minslice,Int_t maxslice) | |
520 | { | |
521 | //Process all slices [minslice,maxslice]. | |
522 | fGlobalMerger = new AliL3HoughGlobalMerger(minslice,maxslice); | |
523 | ||
524 | for(Int_t i=minslice; i<=maxslice; i++) | |
525 | { | |
526 | ReadData(i); | |
527 | Transform(); | |
528 | if(fAddHistograms) { | |
529 | if(fVersion != 4) | |
530 | AddAllHistograms(); | |
531 | else | |
532 | AddAllHistogramsRows(); | |
533 | } | |
534 | FindTrackCandidates(); | |
535 | //Evaluate(); | |
536 | //fGlobalMerger->FillTracks(fTracks[0],i); | |
537 | } | |
538 | } | |
539 | ||
540 | void AliL3Hough::ReadData(Int_t slice,Int_t eventnr) | |
541 | { | |
542 | //Read data from files, binary or root. | |
543 | ||
544 | #ifdef use_aliroot | |
545 | if(fEvent!=eventnr) //just be sure that index is empty for new event | |
546 | AliL3FileHandler::CleanStaticIndex(); | |
547 | #endif | |
548 | fCurrentSlice = slice; | |
549 | ||
550 | for(Int_t i=0; i<fNPatches; i++) | |
551 | { | |
552 | fMemHandler[i]->Free(); | |
553 | UInt_t ndigits=0; | |
554 | AliL3DigitRowData *digits =0; | |
555 | Char_t name[256]; | |
556 | fMemHandler[i]->Init(slice,i); | |
557 | if(fBinary)//take input data from binary files | |
558 | { | |
559 | if(fUse8bits) | |
560 | sprintf(name,"%s/binaries/digits_c8_%d_%d_%d.raw",fPath,eventnr,slice,i); | |
561 | else | |
562 | sprintf(name,"%s/binaries/digits_%d_%d_%d.raw",fPath,eventnr,slice,i); | |
563 | ||
564 | fMemHandler[i]->SetBinaryInput(name); | |
565 | digits = (AliL3DigitRowData *)fMemHandler[i]->CompBinary2Memory(ndigits); | |
566 | fMemHandler[i]->CloseBinaryInput(); | |
567 | } | |
568 | else //read data from root file | |
569 | { | |
570 | #ifdef use_aliroot | |
571 | if(fEvent!=eventnr) | |
572 | fMemHandler[i]->FreeDigitsTree();//or else the new event is not loaded | |
573 | digits=(AliL3DigitRowData *)fMemHandler[i]->AliAltroDigits2Memory(ndigits,eventnr); | |
574 | #else | |
575 | cerr<<"You cannot read from rootfile now"<<endl; | |
576 | #endif | |
577 | } | |
578 | ||
579 | //Set the pointer to the TPCRawStream in case of fast raw data reading | |
580 | fHoughTransformer[i]->SetTPCRawStream(fMemHandler[i]->GetTPCRawStream()); | |
581 | ||
582 | //set input data and init transformer | |
583 | fHoughTransformer[i]->SetInputData(ndigits,digits); | |
584 | fHoughTransformer[i]->Init(slice,i,fNEtaSegments); | |
585 | } | |
586 | ||
587 | fEvent=eventnr; | |
588 | } | |
589 | ||
590 | void AliL3Hough::Transform(Int_t *rowrange) | |
591 | { | |
592 | //Transform all data given to the transformer within the given slice | |
593 | //(after ReadData(slice)) | |
594 | ||
595 | Double_t initTime,cpuTime; | |
596 | initTime = GetCpuTime(); | |
597 | Int_t patchorder[6] = {5,2,0,1,3,4}; //The order in which patches are processed | |
598 | // Int_t patchorder[6] = {0,1,2,3,4,5}; //The order in which patches are processed | |
599 | // Int_t patchorder[6] = {5,4,3,2,1,0}; //The order in which patches are processed | |
600 | // Int_t patchorder[6] = {5,2,4,3,1,0}; //The order in which patches are processed | |
601 | fLastPatch=-1; | |
602 | for(Int_t i=0; i<fNPatches; i++) | |
603 | { | |
604 | // In case of Row transformer reset the arrays only once | |
605 | if((fVersion != 4) || (i == 0)) { | |
606 | fBenchmark->Start("Hough Reset"); | |
607 | fHoughTransformer[0]->Reset();//Reset the histograms | |
608 | fBenchmark->Stop("Hough Reset"); | |
609 | } | |
610 | fBenchmark->Start("Hough Transform"); | |
611 | PrepareForNextPatch(patchorder[i]); | |
612 | if(!rowrange) { | |
613 | char buf[256]; | |
614 | sprintf(buf,"Patch %d",patchorder[i]); | |
615 | fBenchmark->Start(buf); | |
616 | fHoughTransformer[patchorder[i]]->SetLastPatch(fLastPatch); | |
617 | fHoughTransformer[patchorder[i]]->TransformCircle(); | |
618 | fBenchmark->Stop(buf); | |
619 | } | |
620 | else | |
621 | fHoughTransformer[i]->TransformCircleC(rowrange,1); | |
622 | fBenchmark->Stop("Hough Transform"); | |
623 | fLastPatch=patchorder[i]; | |
624 | } | |
625 | cpuTime = GetCpuTime() - initTime; | |
626 | LOG(AliL3Log::kInformational,"AliL3Hough::Transform()","Timing") | |
627 | <<"Transform done in average per patch of "<<cpuTime*1000/fNPatches<<" ms"<<ENDLOG; | |
628 | } | |
629 | ||
630 | void AliL3Hough::MergePatches() | |
631 | { | |
632 | // Merge patches if they are not summed | |
633 | if(fAddHistograms) //Nothing to merge here | |
634 | return; | |
635 | fMerger->MergePatches(kTRUE); | |
636 | } | |
637 | ||
638 | void AliL3Hough::MergeInternally() | |
639 | { | |
640 | // Merge patches internally | |
641 | if(fAddHistograms) | |
642 | fInterMerger->FillTracks(fTracks[0]); | |
643 | else | |
644 | fInterMerger->FillTracks(fMerger->GetOutTracks()); | |
645 | ||
646 | fInterMerger->MMerge(); | |
647 | } | |
648 | ||
649 | void AliL3Hough::ProcessSliceIter() | |
650 | { | |
651 | //Process current slice (after ReadData(slice)) iteratively. | |
652 | ||
653 | if(!fAddHistograms) | |
654 | { | |
655 | for(Int_t i=0; i<fNPatches; i++) | |
656 | { | |
657 | ProcessPatchIter(i); | |
658 | fMerger->FillTracks(fTracks[i],i); //Copy tracks to merger | |
659 | } | |
660 | } | |
661 | else | |
662 | { | |
663 | for(Int_t i=0; i<10; i++) | |
664 | { | |
665 | Transform(); | |
666 | AddAllHistograms(); | |
667 | InitEvaluate(); | |
668 | AliL3HoughBaseTransformer *tr = fHoughTransformer[0]; | |
669 | for(Int_t j=0; j<fNEtaSegments; j++) | |
670 | { | |
671 | AliL3Histogram *hist = tr->GetHistogram(j); | |
672 | if(hist->GetNEntries()==0) continue; | |
673 | fPeakFinder->Reset(); | |
674 | fPeakFinder->SetHistogram(hist); | |
675 | fPeakFinder->FindAbsMaxima(); | |
676 | AliL3HoughTrack *track = (AliL3HoughTrack*)fTracks[0]->NextTrack(); | |
677 | track->SetTrackParameters(fPeakFinder->GetXPeak(0),fPeakFinder->GetYPeak(0),fPeakFinder->GetWeight(0)); | |
678 | track->SetEtaIndex(j); | |
679 | track->SetEta(tr->GetEta(j,fCurrentSlice)); | |
680 | for(Int_t k=0; k<fNPatches; k++) | |
681 | { | |
682 | fEval[i]->SetNumOfPadsToLook(2); | |
683 | fEval[i]->SetNumOfRowsToMiss(2); | |
684 | fEval[i]->RemoveFoundTracks(); | |
685 | /* | |
686 | Int_t nrows=0; | |
687 | if(!fEval[i]->LookInsideRoad(track,nrows)) | |
688 | { | |
689 | fTracks[0]->Remove(fTracks[0]->GetNTracks()-1); | |
690 | fTracks[0]->Compress(); | |
691 | } | |
692 | */ | |
693 | } | |
694 | } | |
695 | ||
696 | } | |
697 | ||
698 | } | |
699 | } | |
700 | ||
701 | void AliL3Hough::ProcessPatchIter(Int_t patch) | |
702 | { | |
703 | //Process patch in a iterative way. | |
704 | //transform + peakfinding + evaluation + transform +... | |
705 | ||
706 | Int_t numoftries = 5; | |
707 | AliL3HoughBaseTransformer *tr = fHoughTransformer[patch]; | |
708 | AliL3TrackArray *tracks = fTracks[patch]; | |
709 | tracks->Reset(); | |
710 | AliL3HoughEval *ev = fEval[patch]; | |
711 | ev->InitTransformer(tr); | |
712 | //ev->RemoveFoundTracks(); | |
713 | ev->SetNumOfRowsToMiss(3); | |
714 | ev->SetNumOfPadsToLook(2); | |
715 | AliL3Histogram *hist; | |
716 | for(Int_t t=0; t<numoftries; t++) | |
717 | { | |
718 | tr->Reset(); | |
719 | tr->TransformCircle(); | |
720 | for(Int_t i=0; i<fNEtaSegments; i++) | |
721 | { | |
722 | hist = tr->GetHistogram(i); | |
723 | if(hist->GetNEntries()==0) continue; | |
724 | fPeakFinder->Reset(); | |
725 | fPeakFinder->SetHistogram(hist); | |
726 | fPeakFinder->FindAbsMaxima(); | |
727 | //fPeakFinder->FindPeak1(); | |
728 | AliL3HoughTrack *track = (AliL3HoughTrack*)tracks->NextTrack(); | |
729 | track->SetTrackParameters(fPeakFinder->GetXPeak(0),fPeakFinder->GetYPeak(0),fPeakFinder->GetWeight(0)); | |
730 | track->SetEtaIndex(i); | |
731 | track->SetEta(tr->GetEta(i,fCurrentSlice)); | |
732 | /* | |
733 | Int_t nrows=0; | |
734 | if(!ev->LookInsideRoad(track,nrows)) | |
735 | { | |
736 | tracks->Remove(tracks->GetNTracks()-1); | |
737 | tracks->Compress(); | |
738 | } | |
739 | */ | |
740 | } | |
741 | } | |
742 | fTracks[0]->QSort(); | |
743 | LOG(AliL3Log::kInformational,"AliL3Hough::ProcessPatch","NTracks") | |
744 | <<AliL3Log::kDec<<"Found "<<tracks->GetNTracks()<<" tracks in patch "<<patch<<ENDLOG; | |
745 | } | |
746 | ||
747 | void AliL3Hough::AddAllHistograms() | |
748 | { | |
749 | //Add the histograms within one etaslice. | |
750 | //Resulting histogram are in patch=0. | |
751 | ||
752 | Double_t initTime,cpuTime; | |
753 | initTime = GetCpuTime(); | |
754 | fBenchmark->Start("Add Histograms"); | |
755 | for(Int_t i=0; i<fNEtaSegments; i++) | |
756 | { | |
757 | AliL3Histogram *hist0 = fHoughTransformer[0]->GetHistogram(i); | |
758 | for(Int_t j=1; j<fNPatches; j++) | |
759 | { | |
760 | AliL3Histogram *hist = fHoughTransformer[j]->GetHistogram(i); | |
761 | hist0->Add(hist); | |
762 | } | |
763 | } | |
764 | fBenchmark->Stop("Add Histograms"); | |
765 | fAddHistograms = kTRUE; | |
766 | cpuTime = GetCpuTime() - initTime; | |
767 | LOG(AliL3Log::kInformational,"AliL3Hough::AddAllHistograms()","Timing") | |
768 | <<"Adding histograms in "<<cpuTime*1000<<" ms"<<ENDLOG; | |
769 | } | |
770 | ||
771 | void AliL3Hough::AddAllHistogramsRows() | |
772 | { | |
773 | //Add the histograms within one etaslice. | |
774 | //Resulting histogram are in patch=0. | |
775 | ||
776 | Double_t initTime,cpuTime; | |
777 | initTime = GetCpuTime(); | |
778 | fBenchmark->Start("Add HistogramsRows"); | |
779 | ||
780 | UChar_t lastpatchlastrow = AliL3Transform::GetLastRowOnDDL(fLastPatch)+1; | |
781 | ||
782 | UChar_t *tracklastrow = ((AliL3HoughTransformerRow *)fHoughTransformer[0])->GetTrackLastRow(); | |
783 | ||
784 | for(Int_t i=0; i<fNEtaSegments; i++) | |
785 | { | |
786 | UChar_t *gapcount = ((AliL3HoughTransformerRow *)fHoughTransformer[0])->GetGapCount(i); | |
787 | UChar_t *currentrowcount = ((AliL3HoughTransformerRow *)fHoughTransformer[0])->GetCurrentRowCount(i); | |
788 | ||
789 | AliL3Histogram *hist = fHoughTransformer[0]->GetHistogram(i); | |
790 | Int_t xmin = hist->GetFirstXbin(); | |
791 | Int_t xmax = hist->GetLastXbin(); | |
792 | Int_t ymin = hist->GetFirstYbin(); | |
793 | Int_t ymax = hist->GetLastYbin(); | |
794 | Int_t nxbins = hist->GetNbinsX()+2; | |
795 | ||
796 | for(Int_t ybin=ymin; ybin<=ymax; ybin++) | |
797 | { | |
798 | for(Int_t xbin=xmin; xbin<=xmax; xbin++) | |
799 | { | |
800 | Int_t bin = xbin + ybin*nxbins; //Int_t bin = hist->GetBin(xbin,ybin); | |
801 | if(gapcount[bin] < MAX_N_GAPS) { | |
802 | if(tracklastrow[bin] > lastpatchlastrow) { | |
803 | if(lastpatchlastrow > currentrowcount[bin]) | |
804 | gapcount[bin] += (lastpatchlastrow-currentrowcount[bin]-1); | |
805 | } | |
806 | else { | |
807 | if(tracklastrow[bin] > currentrowcount[bin]) | |
808 | gapcount[bin] += (tracklastrow[bin]-currentrowcount[bin]-1); | |
809 | } | |
810 | if(gapcount[bin] < MAX_N_GAPS) | |
811 | hist->AddBinContent(bin,(159-gapcount[bin])); | |
812 | } | |
813 | } | |
814 | } | |
815 | } | |
816 | ||
817 | fBenchmark->Stop("Add HistogramsRows"); | |
818 | fAddHistograms = kTRUE; | |
819 | cpuTime = GetCpuTime() - initTime; | |
820 | LOG(AliL3Log::kInformational,"AliL3Hough::AddAllHistogramsRows()","Timing") | |
821 | <<"Adding histograms in "<<cpuTime*1000<<" ms"<<ENDLOG; | |
822 | } | |
823 | ||
824 | void AliL3Hough::PrepareForNextPatch(Int_t nextpatch) | |
825 | { | |
826 | // Prepare the parameter space for the processing of | |
827 | // the next read patch. According to the already | |
828 | // accumulated number of gaps in parameter space | |
829 | // bins, the routine updates the dynamic | |
830 | // pointers used in order to jump rapidly during the | |
831 | // filling of the parameter space. | |
832 | ||
833 | char buf[256]; | |
834 | sprintf(buf,"Prepare For Patch %d",nextpatch); | |
835 | fBenchmark->Start(buf); | |
836 | ||
837 | UChar_t lastpatchlastrow; | |
838 | if(fLastPatch == -1) | |
839 | lastpatchlastrow = 0; | |
840 | else | |
841 | lastpatchlastrow = AliL3Transform::GetLastRowOnDDL(fLastPatch)+1; | |
842 | UChar_t nextpatchfirstrow; | |
843 | if(nextpatch==0) | |
844 | nextpatchfirstrow = 0; | |
845 | else | |
846 | nextpatchfirstrow = AliL3Transform::GetFirstRowOnDDL(nextpatch)-1; | |
847 | ||
848 | UChar_t *trackfirstrow = ((AliL3HoughTransformerRow *)fHoughTransformer[0])->GetTrackFirstRow(); | |
849 | UChar_t *tracklastrow = ((AliL3HoughTransformerRow *)fHoughTransformer[0])->GetTrackLastRow(); | |
850 | ||
851 | for(Int_t i=0; i<fNEtaSegments; i++) | |
852 | { | |
853 | UChar_t *gapcount = ((AliL3HoughTransformerRow *)fHoughTransformer[0])->GetGapCount(i); | |
854 | UChar_t *currentrowcount = ((AliL3HoughTransformerRow *)fHoughTransformer[0])->GetCurrentRowCount(i); | |
855 | UChar_t *prevbin = ((AliL3HoughTransformerRow *)fHoughTransformer[0])->GetPrevBin(i); | |
856 | UChar_t *nextbin = ((AliL3HoughTransformerRow *)fHoughTransformer[0])->GetNextBin(i); | |
857 | UChar_t *nextrow = ((AliL3HoughTransformerRow *)fHoughTransformer[0])->GetNextRow(i); | |
858 | ||
859 | AliL3Histogram *hist = fHoughTransformer[0]->GetHistogram(i); | |
860 | Int_t xmin = hist->GetFirstXbin(); | |
861 | Int_t xmax = hist->GetLastXbin(); | |
862 | Int_t ymin = hist->GetFirstYbin(); | |
863 | Int_t ymax = hist->GetLastYbin(); | |
864 | Int_t nxbins = hist->GetNbinsX()+2; | |
865 | ||
866 | if(fLastPatch != -1) { | |
867 | UChar_t lastyvalue = 0; | |
868 | Int_t endybin = ymin - 1; | |
869 | for(Int_t ybin=nextrow[ymin]; ybin<=ymax; ybin = nextrow[++ybin]) | |
870 | { | |
871 | UChar_t lastxvalue = 0; | |
872 | UChar_t maxvalue = 0; | |
873 | Int_t endxbin = xmin - 1; | |
874 | for(Int_t xbin=xmin; xbin<=xmax; xbin++) | |
875 | { | |
876 | Int_t bin = xbin + ybin*nxbins; | |
877 | UChar_t value = 0; | |
878 | if(gapcount[bin] < MAX_N_GAPS) { | |
879 | if(tracklastrow[bin] > lastpatchlastrow) { | |
880 | if(lastpatchlastrow > currentrowcount[bin]) | |
881 | gapcount[bin] += (lastpatchlastrow-currentrowcount[bin]-1); | |
882 | } | |
883 | else { | |
884 | if(tracklastrow[bin] > currentrowcount[bin]) | |
885 | gapcount[bin] += (tracklastrow[bin]-currentrowcount[bin]-1); | |
886 | } | |
887 | if(gapcount[bin] < MAX_N_GAPS) { | |
888 | value = 1; | |
889 | maxvalue = 1; | |
890 | if(trackfirstrow[bin] < nextpatchfirstrow) | |
891 | currentrowcount[bin] = nextpatchfirstrow; | |
892 | else | |
893 | currentrowcount[bin] = trackfirstrow[bin]; | |
894 | } | |
895 | } | |
896 | if(value > 0) | |
897 | { | |
898 | nextbin[xbin + ybin*nxbins] = (UChar_t)xbin; | |
899 | prevbin[xbin + ybin*nxbins] = (UChar_t)xbin; | |
900 | if(value > lastxvalue) | |
901 | { | |
902 | UChar_t *tempnextbin = nextbin + endxbin + 1 + ybin*nxbins; | |
903 | memset(tempnextbin,(UChar_t)xbin,xbin-endxbin-1); | |
904 | } | |
905 | endxbin = xbin; | |
906 | } | |
907 | else | |
908 | { | |
909 | prevbin[xbin + ybin*nxbins] = (UChar_t)endxbin; | |
910 | } | |
911 | lastxvalue = value; | |
912 | } | |
913 | UChar_t *tempnextbin = nextbin + endxbin + 1 + ybin*nxbins; | |
914 | memset(tempnextbin,(UChar_t)(xmax+1),xmax-endxbin); | |
915 | if(maxvalue > 0) | |
916 | { | |
917 | nextrow[ybin] = (UChar_t)ybin; | |
918 | if(maxvalue > lastyvalue) | |
919 | { | |
920 | UChar_t *tempnextrow = nextrow + endybin + 1; | |
921 | memset(tempnextrow,(UChar_t)ybin,ybin-endybin-1); | |
922 | } | |
923 | endybin = ybin; | |
924 | } | |
925 | lastyvalue = maxvalue; | |
926 | } | |
927 | UChar_t *tempnextrow = nextrow + endybin + 1; | |
928 | memset(tempnextrow,(UChar_t)(ymax+1),ymax-endybin+1); | |
929 | } | |
930 | else { | |
931 | UChar_t lastyvalue = 0; | |
932 | Int_t endybin = ymin - 1; | |
933 | for(Int_t ybin=ymin; ybin<=ymax; ybin++) | |
934 | { | |
935 | UChar_t maxvalue = 0; | |
936 | for(Int_t xbin=xmin; xbin<=xmax; xbin++) | |
937 | { | |
938 | Int_t bin = xbin + ybin*nxbins; | |
939 | if(gapcount[bin] < MAX_N_GAPS) { | |
940 | maxvalue = 1; | |
941 | if(trackfirstrow[bin] < nextpatchfirstrow) | |
942 | currentrowcount[bin] = nextpatchfirstrow; | |
943 | else | |
944 | currentrowcount[bin] = trackfirstrow[bin]; | |
945 | } | |
946 | } | |
947 | if(maxvalue > 0) | |
948 | { | |
949 | nextrow[ybin] = (UChar_t)ybin; | |
950 | if(maxvalue > lastyvalue) | |
951 | { | |
952 | UChar_t *tempnextrow = nextrow + endybin + 1; | |
953 | memset(tempnextrow,(UChar_t)ybin,ybin-endybin-1); | |
954 | } | |
955 | endybin = ybin; | |
956 | } | |
957 | lastyvalue = maxvalue; | |
958 | } | |
959 | UChar_t *tempnextrow = nextrow + endybin + 1; | |
960 | memset(tempnextrow,(UChar_t)(ymax+1),ymax-endybin+1); | |
961 | } | |
962 | } | |
963 | ||
964 | fBenchmark->Stop(buf); | |
965 | } | |
966 | ||
967 | void AliL3Hough::AddTracks() | |
968 | { | |
969 | // Add current slice slice tracks to the global list of found tracks | |
970 | if(!fTracks[0]) | |
971 | { | |
972 | cerr<<"AliL3Hough::AddTracks : No tracks"<<endl; | |
973 | return; | |
974 | } | |
975 | AliL3TrackArray *tracks = fTracks[0]; | |
976 | for(Int_t i=0; i<tracks->GetNTracks(); i++) | |
977 | { | |
978 | AliL3Track *track = tracks->GetCheckedTrack(i); | |
979 | if(!track) continue; | |
980 | if(track->GetNHits()!=1) cerr<<"NHITS "<<track->GetNHits()<<endl; | |
981 | UInt_t *ids = track->GetHitNumbers(); | |
982 | ids[0] = (fCurrentSlice&0x7f)<<25; | |
983 | } | |
984 | ||
985 | fGlobalTracks->AddTracks(fTracks[0],0,fCurrentSlice); | |
986 | } | |
987 | ||
988 | void AliL3Hough::FindTrackCandidatesRow() | |
989 | { | |
990 | // Find AliL3HoughTransformerRow track candidates | |
991 | if(fVersion != 4) { | |
992 | LOG(AliL3Log::kError,"AliL3Hough::FindTrackCandidatesRow()","") | |
993 | <<"Incompatible Peak Finder version!"<<ENDLOG; | |
994 | return; | |
995 | } | |
996 | ||
997 | //Look for peaks in histograms, and find the track candidates | |
998 | Int_t npatches; | |
999 | if(fAddHistograms) | |
1000 | npatches = 1; //Histograms have been added. | |
1001 | else | |
1002 | npatches = fNPatches; | |
1003 | ||
1004 | Double_t initTime,cpuTime; | |
1005 | initTime = GetCpuTime(); | |
1006 | fBenchmark->Start("Find Maxima"); | |
1007 | for(Int_t i=0; i<npatches; i++) | |
1008 | { | |
1009 | AliL3HoughBaseTransformer *tr = fHoughTransformer[i]; | |
1010 | AliL3Histogram *h = tr->GetHistogram(0); | |
1011 | Float_t deltax = h->GetBinWidthX()*AliL3HoughTransformerRow::GetDAlpha(); | |
1012 | Float_t deltay = h->GetBinWidthY()*AliL3HoughTransformerRow::GetDAlpha(); | |
1013 | Float_t deltaeta = (tr->GetEtaMax()-tr->GetEtaMin())/tr->GetNEtaSegments()*AliL3HoughTransformerRow::GetDEta(); | |
1014 | Float_t zvertex = tr->GetZVertex(); | |
1015 | fTracks[i]->Reset(); | |
1016 | fPeakFinder->Reset(); | |
1017 | ||
1018 | for(Int_t j=0; j<fNEtaSegments; j++) | |
1019 | { | |
1020 | AliL3Histogram *hist = tr->GetHistogram(j); | |
1021 | if(hist->GetNEntries()==0) continue; | |
1022 | fPeakFinder->SetHistogram(hist); | |
1023 | fPeakFinder->SetEtaSlice(j); | |
1024 | fPeakFinder->SetTrackLUTs(((AliL3HoughTransformerRow *)tr)->GetTrackNRows(),((AliL3HoughTransformerRow *)tr)->GetTrackFirstRow(),((AliL3HoughTransformerRow *)tr)->GetTrackLastRow(),((AliL3HoughTransformerRow *)tr)->GetNextRow(j)); | |
1025 | #ifdef do_mc | |
1026 | LOG(AliL3Log::kInformational,"AliL3Hough::FindTrackCandidates()","") | |
1027 | <<"Starting "<<j<<" etaslice"<<ENDLOG; | |
1028 | #endif | |
1029 | fPeakFinder->SetThreshold(fPeakThreshold[i]); | |
1030 | fPeakFinder->FindAdaptedRowPeaks(1,0,0);//Maxima finder for HoughTransformerRow | |
1031 | } | |
1032 | ||
1033 | for(Int_t k=0; k<fPeakFinder->GetEntries(); k++) | |
1034 | { | |
1035 | // if(fPeakFinder->GetWeight(k) < 0) continue; | |
1036 | AliL3HoughTrack *track = (AliL3HoughTrack*)fTracks[i]->NextTrack(); | |
1037 | Double_t starteta = tr->GetEta(fPeakFinder->GetStartEta(k),fCurrentSlice); | |
1038 | Double_t endeta = tr->GetEta(fPeakFinder->GetEndEta(k),fCurrentSlice); | |
1039 | Double_t eta = (starteta+endeta)/2.0; | |
1040 | track->SetTrackParametersRow(fPeakFinder->GetXPeak(k),fPeakFinder->GetYPeak(k),eta,fPeakFinder->GetWeight(k)); | |
1041 | track->SetPterr(deltax); track->SetPsierr(deltay); track->SetTglerr(deltaeta); | |
1042 | track->SetBinXY(fPeakFinder->GetXPeak(k),fPeakFinder->GetYPeak(k),fPeakFinder->GetXPeakSize(k),fPeakFinder->GetYPeakSize(k)); | |
1043 | track->SetZ0(zvertex); | |
1044 | Int_t etaindex = (fPeakFinder->GetStartEta(k)+fPeakFinder->GetEndEta(k))/2; | |
1045 | track->SetEtaIndex(etaindex); | |
1046 | Int_t rows[2]; | |
1047 | ((AliL3HoughTransformerRow *)tr)->GetTrackLength(fPeakFinder->GetXPeak(k),fPeakFinder->GetYPeak(k),rows); | |
1048 | track->SetRowRange(rows[0],rows[1]); | |
1049 | track->SetSector(fCurrentSlice); | |
1050 | track->SetSlice(fCurrentSlice); | |
1051 | #ifdef do_mc | |
1052 | Int_t label = tr->GetTrackID(etaindex,fPeakFinder->GetXPeak(k),fPeakFinder->GetYPeak(k)); | |
1053 | track->SetMCid(label); | |
1054 | #endif | |
1055 | } | |
1056 | LOG(AliL3Log::kInformational,"AliL3Hough::FindTrackCandidates()","") | |
1057 | <<"Found "<<fTracks[i]->GetNTracks()<<" tracks in slice "<<fCurrentSlice<<ENDLOG; | |
1058 | fTracks[i]->QSort(); | |
1059 | } | |
1060 | fBenchmark->Stop("Find Maxima"); | |
1061 | cpuTime = GetCpuTime() - initTime; | |
1062 | LOG(AliL3Log::kInformational,"AliL3Hough::FindTrackCandidates()","Timing") | |
1063 | <<"Maxima finding done in "<<cpuTime*1000<<" ms"<<ENDLOG; | |
1064 | } | |
1065 | ||
1066 | void AliL3Hough::FindTrackCandidates() | |
1067 | { | |
1068 | // Find AliL3HoughTransformer track candidates | |
1069 | if(fVersion == 4) { | |
1070 | LOG(AliL3Log::kError,"AliL3Hough::FindTrackCandidatesRow()","") | |
1071 | <<"Incompatible Peak Finder version!"<<ENDLOG; | |
1072 | return; | |
1073 | } | |
1074 | ||
1075 | Int_t npatches; | |
1076 | if(fAddHistograms) | |
1077 | npatches = 1; //Histograms have been added. | |
1078 | else | |
1079 | npatches = fNPatches; | |
1080 | ||
1081 | Double_t initTime,cpuTime; | |
1082 | initTime = GetCpuTime(); | |
1083 | fBenchmark->Start("Find Maxima"); | |
1084 | for(Int_t i=0; i<npatches; i++) | |
1085 | { | |
1086 | AliL3HoughBaseTransformer *tr = fHoughTransformer[i]; | |
1087 | fTracks[i]->Reset(); | |
1088 | ||
1089 | for(Int_t j=0; j<fNEtaSegments; j++) | |
1090 | { | |
1091 | AliL3Histogram *hist = tr->GetHistogram(j); | |
1092 | if(hist->GetNEntries()==0) continue; | |
1093 | fPeakFinder->Reset(); | |
1094 | fPeakFinder->SetHistogram(hist); | |
1095 | #ifdef do_mc | |
1096 | cout<<"Starting "<<j<<" etaslice"<<endl; | |
1097 | #endif | |
1098 | fPeakFinder->SetThreshold(fPeakThreshold[i]); | |
1099 | fPeakFinder->FindAdaptedPeaks(fKappaSpread,fPeakRatio); | |
1100 | ||
1101 | for(Int_t k=0; k<fPeakFinder->GetEntries(); k++) | |
1102 | { | |
1103 | AliL3HoughTrack *track = (AliL3HoughTrack*)fTracks[i]->NextTrack(); | |
1104 | track->SetTrackParameters(fPeakFinder->GetXPeak(k),fPeakFinder->GetYPeak(k),fPeakFinder->GetWeight(k)); | |
1105 | track->SetEtaIndex(j); | |
1106 | track->SetEta(tr->GetEta(j,fCurrentSlice)); | |
1107 | track->SetRowRange(AliL3Transform::GetFirstRow(0),AliL3Transform::GetLastRow(5)); | |
1108 | } | |
1109 | } | |
1110 | cout<<"Found "<<fTracks[i]->GetNTracks()<<" tracks in patch "<<i<<endl; | |
1111 | fTracks[i]->QSort(); | |
1112 | } | |
1113 | fBenchmark->Stop("Find Maxima"); | |
1114 | cpuTime = GetCpuTime() - initTime; | |
1115 | LOG(AliL3Log::kInformational,"AliL3Hough::FindTrackCandidates()","Timing") | |
1116 | <<"Maxima finding done in "<<cpuTime*1000<<" ms"<<ENDLOG; | |
1117 | } | |
1118 | ||
1119 | void AliL3Hough::InitEvaluate() | |
1120 | { | |
1121 | //Pass the transformer objects to the AliL3HoughEval objects: | |
1122 | //This will provide the evaluation objects with all the necessary | |
1123 | //data and parameters it needs. | |
1124 | ||
1125 | for(Int_t i=0; i<fNPatches; i++) | |
1126 | fEval[i]->InitTransformer(fHoughTransformer[i]); | |
1127 | } | |
1128 | ||
1129 | Int_t AliL3Hough::Evaluate(Int_t roadwidth,Int_t nrowstomiss) | |
1130 | { | |
1131 | //Evaluate the tracks, by looking along the road in the raw data. | |
1132 | //If track does not cross all padrows - rows2miss, it is removed from the arrray. | |
1133 | //If histograms were not added, the check is done locally in patch, | |
1134 | //meaning that nrowstomiss is the number of padrows the road can miss with respect | |
1135 | //to the number of rows in the patch. | |
1136 | //If the histograms were added, the comparison is done globally in the _slice_, | |
1137 | //meaing that nrowstomiss is the number of padrows the road can miss with | |
1138 | //respect to the total number of padrows in the slice. | |
1139 | // | |
1140 | //Return value = number of tracks which were removed (only in case of fAddHistograms) | |
1141 | ||
1142 | if(!fTracks[0]) | |
1143 | { | |
1144 | LOG(AliL3Log::kError,"AliL3Hough::Evaluate","Track Array") | |
1145 | <<"No tracks to work with..."<<ENDLOG; | |
1146 | return 0; | |
1147 | } | |
1148 | ||
1149 | Int_t removedtracks=0; | |
1150 | AliL3TrackArray *tracks=0; | |
1151 | ||
1152 | if(fAddHistograms) | |
1153 | { | |
1154 | tracks = fTracks[0]; | |
1155 | for(Int_t i=0; i<tracks->GetNTracks(); i++) | |
1156 | { | |
1157 | AliL3Track *track = tracks->GetCheckedTrack(i); | |
1158 | if(!track) continue; | |
1159 | track->SetNHits(0); | |
1160 | } | |
1161 | } | |
1162 | ||
1163 | for(Int_t i=0; i<fNPatches; i++) | |
1164 | EvaluatePatch(i,roadwidth,nrowstomiss); | |
1165 | ||
1166 | //Here we check the tracks globally; | |
1167 | //how many good rows (padrows with signal) | |
1168 | //did it cross in the slice | |
1169 | if(fAddHistograms) | |
1170 | { | |
1171 | for(Int_t j=0; j<tracks->GetNTracks(); j++) | |
1172 | { | |
1173 | AliL3HoughTrack *track = (AliL3HoughTrack*)tracks->GetCheckedTrack(j); | |
1174 | ||
1175 | if(track->GetNHits() < AliL3Transform::GetNRows() - nrowstomiss) | |
1176 | { | |
1177 | tracks->Remove(j); | |
1178 | removedtracks++; | |
1179 | } | |
1180 | } | |
1181 | tracks->Compress(); | |
1182 | tracks->QSort(); | |
1183 | } | |
1184 | ||
1185 | return removedtracks; | |
1186 | } | |
1187 | ||
1188 | void AliL3Hough::EvaluatePatch(Int_t i,Int_t roadwidth,Int_t nrowstomiss) | |
1189 | { | |
1190 | //Evaluate patch i. | |
1191 | ||
1192 | fEval[i]->InitTransformer(fHoughTransformer[i]); | |
1193 | fEval[i]->SetNumOfPadsToLook(roadwidth); | |
1194 | fEval[i]->SetNumOfRowsToMiss(nrowstomiss); | |
1195 | //fEval[i]->RemoveFoundTracks(); | |
1196 | ||
1197 | AliL3TrackArray *tracks=0; | |
1198 | ||
1199 | if(!fAddHistograms) | |
1200 | tracks = fTracks[i]; | |
1201 | else | |
1202 | tracks = fTracks[0]; | |
1203 | ||
1204 | Int_t nrows=0; | |
1205 | for(Int_t j=0; j<tracks->GetNTracks(); j++) | |
1206 | { | |
1207 | AliL3HoughTrack *track = (AliL3HoughTrack*)tracks->GetCheckedTrack(j); | |
1208 | if(!track) | |
1209 | { | |
1210 | LOG(AliL3Log::kWarning,"AliL3Hough::EvaluatePatch","Track array") | |
1211 | <<"Track object missing!"<<ENDLOG; | |
1212 | continue; | |
1213 | } | |
1214 | nrows=0; | |
1215 | Int_t rowrange[2] = {AliL3Transform::GetFirstRow(i),AliL3Transform::GetLastRow(i)}; | |
1216 | Bool_t result = fEval[i]->LookInsideRoad(track,nrows,rowrange); | |
1217 | if(fAddHistograms) | |
1218 | { | |
1219 | Int_t pre=track->GetNHits(); | |
1220 | track->SetNHits(pre+nrows); | |
1221 | } | |
1222 | else//the track crossed too few good padrows (padrows with signal) in the patch, so remove it | |
1223 | { | |
1224 | if(result == kFALSE) | |
1225 | tracks->Remove(j); | |
1226 | } | |
1227 | } | |
1228 | ||
1229 | tracks->Compress(); | |
1230 | ||
1231 | } | |
1232 | ||
1233 | void AliL3Hough::MergeEtaSlices() | |
1234 | { | |
1235 | //Merge tracks found in neighbouring eta slices. | |
1236 | //Removes the track with the lower weight. | |
1237 | ||
1238 | fBenchmark->Start("Merge Eta-slices"); | |
1239 | AliL3TrackArray *tracks = fTracks[0]; | |
1240 | if(!tracks) | |
1241 | { | |
1242 | cerr<<"AliL3Hough::MergeEtaSlices : No tracks "<<endl; | |
1243 | return; | |
1244 | } | |
1245 | for(Int_t j=0; j<tracks->GetNTracks(); j++) | |
1246 | { | |
1247 | AliL3HoughTrack *track1 = (AliL3HoughTrack*)tracks->GetCheckedTrack(j); | |
1248 | if(!track1) continue; | |
1249 | for(Int_t k=j+1; k<tracks->GetNTracks(); k++) | |
1250 | { | |
1251 | AliL3HoughTrack *track2 = (AliL3HoughTrack*)tracks->GetCheckedTrack(k); | |
1252 | if(!track2) continue; | |
1253 | if(abs(track1->GetEtaIndex() - track2->GetEtaIndex()) != 1) continue; | |
1254 | if(fabs(track1->GetKappa()-track2->GetKappa()) < 0.006 && | |
1255 | fabs(track1->GetPsi()- track2->GetPsi()) < 0.1) | |
1256 | { | |
1257 | //cout<<"Merging track in slices "<<track1->GetEtaIndex()<<" "<<track2->GetEtaIndex()<<endl; | |
1258 | if(track1->GetWeight() > track2->GetWeight()) | |
1259 | tracks->Remove(k); | |
1260 | else | |
1261 | tracks->Remove(j); | |
1262 | } | |
1263 | } | |
1264 | } | |
1265 | fBenchmark->Stop("Merge Eta-slices"); | |
1266 | tracks->Compress(); | |
1267 | } | |
1268 | ||
1269 | void AliL3Hough::WriteTracks(Char_t *path) | |
1270 | { | |
1271 | // Write found tracks into file | |
1272 | //cout<<"AliL3Hough::WriteTracks : Sorting the tracsk"<<endl; | |
1273 | //fGlobalTracks->QSort(); | |
1274 | ||
1275 | Char_t filename[1024]; | |
1276 | sprintf(filename,"%s/tracks_%d.raw",path,fEvent); | |
1277 | AliL3MemHandler mem; | |
1278 | mem.SetBinaryOutput(filename); | |
1279 | mem.TrackArray2Binary(fGlobalTracks); | |
1280 | mem.CloseBinaryOutput(); | |
1281 | fGlobalTracks->Reset(); | |
1282 | } | |
1283 | ||
1284 | void AliL3Hough::WriteTracks(Int_t slice,Char_t *path) | |
1285 | { | |
1286 | // Write found tracks slice by slice into file | |
1287 | ||
1288 | AliL3MemHandler mem; | |
1289 | Char_t fname[100]; | |
1290 | if(fAddHistograms) | |
1291 | { | |
1292 | sprintf(fname,"%s/tracks_ho_%d_%d.raw",path,fEvent,slice); | |
1293 | mem.SetBinaryOutput(fname); | |
1294 | mem.TrackArray2Binary(fTracks[0]); | |
1295 | mem.CloseBinaryOutput(); | |
1296 | } | |
1297 | else | |
1298 | { | |
1299 | for(Int_t i=0; i<fNPatches; i++) | |
1300 | { | |
1301 | sprintf(fname,"%s/tracks_ho_%d_%d_%d.raw",path,fEvent,slice,i); | |
1302 | mem.SetBinaryOutput(fname); | |
1303 | mem.TrackArray2Binary(fTracks[i]); | |
1304 | mem.CloseBinaryOutput(); | |
1305 | } | |
1306 | } | |
1307 | } | |
1308 | ||
1309 | #ifdef use_aliroot | |
1310 | Int_t AliL3Hough::FillESD(AliESD *esd) | |
1311 | { | |
1312 | // Fill the found hough transform tracks | |
1313 | // into the ESD. The tracks are stored as | |
1314 | // AliESDHLTtrack objects. | |
1315 | ||
1316 | if(!fGlobalTracks) return 0; | |
1317 | Int_t nglobaltracks = 0; | |
1318 | for(Int_t i=0; i<fGlobalTracks->GetNTracks(); i++) | |
1319 | { | |
1320 | AliL3HoughTrack *tpt = (AliL3HoughTrack *)fGlobalTracks->GetCheckedTrack(i); | |
1321 | if(!tpt) continue; | |
1322 | ||
1323 | if(tpt->GetWeight()<0) continue; | |
1324 | AliL3HoughKalmanTrack *tpctrack = new AliL3HoughKalmanTrack(*tpt); | |
1325 | if(!tpctrack) continue; | |
1326 | AliESDtrack *esdtrack2 = new AliESDtrack() ; | |
1327 | esdtrack2->UpdateTrackParams(tpctrack,AliESDtrack::kTPCin); | |
1328 | esd->AddTrack(esdtrack2); | |
1329 | nglobaltracks++; | |
1330 | delete esdtrack2; | |
1331 | delete tpctrack; | |
1332 | } | |
1333 | return nglobaltracks; | |
1334 | } | |
1335 | #endif | |
1336 | ||
1337 | void AliL3Hough::WriteDigits(Char_t *outfile) | |
1338 | { | |
1339 | //Write the current data to a new rootfile. | |
1340 | #ifdef use_aliroot | |
1341 | ||
1342 | for(Int_t i=0; i<fNPatches; i++) | |
1343 | { | |
1344 | AliL3DigitRowData *tempPt = (AliL3DigitRowData*)fHoughTransformer[i]->GetDataPointer(); | |
1345 | fMemHandler[i]->AliDigits2RootFile(tempPt,outfile); | |
1346 | } | |
1347 | #else | |
1348 | cerr<<"AliL3Hough::WriteDigits : You need to compile with AliROOT!"<<endl; | |
1349 | return; | |
1350 | #endif | |
1351 | } | |
1352 | ||
1353 | Double_t AliL3Hough::GetCpuTime() | |
1354 | { | |
1355 | //Return the Cputime in seconds. | |
1356 | struct timeval tv; | |
1357 | gettimeofday( &tv, NULL ); | |
1358 | return tv.tv_sec+(((Double_t)tv.tv_usec)/1000000.); | |
1359 | } | |
1360 | ||
1361 | void *AliL3Hough::ProcessInThread(void *args) | |
1362 | { | |
1363 | // Called in case Hough transform tracking | |
1364 | // is executed in a thread | |
1365 | ||
1366 | AliL3Hough *instance = (AliL3Hough *)args; | |
1367 | Int_t minslice = instance->GetMinSlice(); | |
1368 | Int_t maxslice = instance->GetMaxSlice(); | |
1369 | for(Int_t i=minslice; i<=maxslice; i++) | |
1370 | { | |
1371 | instance->ReadData(i,0); | |
1372 | instance->Transform(); | |
1373 | instance->AddAllHistogramsRows(); | |
1374 | instance->FindTrackCandidatesRow(); | |
1375 | instance->AddTracks(); | |
1376 | } | |
1377 | return (void *)0; | |
1378 | } | |
1379 | ||
1380 | void AliL3Hough::StartProcessInThread(Int_t minslice,Int_t maxslice) | |
1381 | { | |
1382 | // Starts the Hough transform tracking as a | |
1383 | // separate thread. Takes as parameters the | |
1384 | // range of TPC slices (sectors) to be reconstructed | |
1385 | ||
1386 | if(!fThread) { | |
1387 | char buf[255]; | |
1388 | sprintf(buf,"houghtrans_%d_%d",minslice,maxslice); | |
1389 | SetMinMaxSlices(minslice,maxslice); | |
1390 | // fThread = new TThread(buf,(void (*) (void *))&ProcessInThread,(void *)this); | |
1391 | fThread = new TThread(buf,&ProcessInThread,(void *)this); | |
1392 | fThread->Run(); | |
1393 | } | |
1394 | return; | |
1395 | } | |
1396 | ||
1397 | Int_t AliL3Hough::WaitForThreadFinish() | |
1398 | { | |
1399 | // Routine is used in case we run the | |
1400 | // Hough transform tracking in several | |
1401 | // threads and want to sync them before | |
1402 | // writing the results to the ESD | |
1403 | ||
1404 | #if ROOT_VERSION_CODE < 262403 | |
1405 | return TThread::Join(fThread->GetId()); | |
1406 | #else | |
1407 | return fThread->Join(fThread->GetId()); | |
1408 | #endif | |
1409 | } |