]> git.uio.no Git - u/mrichter/AliRoot.git/blame_incremental - HLT/hough/AliL3Hough.cxx
TGenat3 instantiation eliminated. It should be done in advance
[u/mrichter/AliRoot.git] / HLT / hough / AliL3Hough.cxx
... / ...
CommitLineData
1// @(#) $Id$
2
3// Author: Anders Vestbo <mailto:vestbo@fi.uib.no>
4//*-- Copyright &copy 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
39using 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
61ClassImp(AliL3Hough)
62
63AliL3Hough::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
111AliL3Hough::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
148AliL3Hough::~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
178void 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
204void 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
215void 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
243void 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
354void 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/*
397void 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*/
418void 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
441void 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
475void 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
488void 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
501void 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
514void AliL3Hough::DoBench(Char_t *name)
515{
516 fBenchmark->Analyze(name);
517}
518
519void 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
540void 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
590void 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
630void 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
638void 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
649void 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
701void 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
747void 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
771void 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
824void 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
967void 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
988void 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
1066void 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
1119void 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
1129Int_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
1188void 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
1233void 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
1269void 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
1284void 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
1310Int_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
1337void 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
1353Double_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
1361void *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
1380void 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
1397Int_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}