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