59cc2cbe98411b2e4adc12139310d933694ed770
[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 "AliL3HoughMaxFinder.h"
20 #include "AliL3Benchmark.h"
21 #ifdef use_aliroot
22 #include "AliL3FileHandler.h"
23 #else
24 #include "AliL3MemHandler.h"
25 #endif
26 #include "AliL3DataHandler.h"
27 #include "AliL3DigitData.h"
28 #include "AliL3HoughEval.h"
29 #include "AliL3Transform.h"
30 #include "AliL3TrackArray.h"
31 #include "AliL3HoughTrack.h"
32
33 #if GCCVERSION == 3
34 using namespace std;
35 #endif
36
37 /** /class AliL3Hough
38 //<pre>
39 //_____________________________________________________________
40 // AliL3Hough
41 //
42 // Interface class for the Hough transform
43 //
44 // Example how to use:
45 //
46 // AliL3Hough *hough = new AliL3Hough(path,kTRUE,NumberOfEtaSegments);
47 // hough->ReadData(slice);
48 // hough->Transform();
49 // hough->FindTrackCandidates();
50 // 
51 // AliL3TrackArray *tracks = hough->GetTracks(patch);
52 //</pre>
53 */
54
55 ClassImp(AliL3Hough)
56
57 AliL3Hough::AliL3Hough()
58 {
59   //Constructor
60   
61   fBinary        = kFALSE;
62   fAddHistograms = kFALSE;
63   fDoIterative   = kFALSE; 
64   fWriteDigits   = kFALSE;
65   fUse8bits      = kFALSE;
66
67   fMemHandler       = 0;
68   fHoughTransformer = 0;
69   fEval             = 0;
70   fPeakFinder       = 0;
71   fTracks           = 0;
72   fMerger           = 0;
73   fInterMerger      = 0;
74   fGlobalMerger     = 0;
75   fBenchmark        = 0;
76   
77   fNEtaSegments     = 0;
78   fNPatches         = 0;
79   fVersion          = 0;
80   fCurrentSlice     = 0;
81   
82   SetTransformerParams();
83   SetThreshold();
84   SetNSaveIterations();
85   SetPeakThreshold();
86 }
87
88 AliL3Hough::AliL3Hough(Char_t *path,Bool_t binary,Int_t n_eta_segments,Bool_t bit8,Int_t tv)
89 {
90   //Default ctor.
91
92   fBinary = binary;
93   strcpy(fPath,path);
94   fNEtaSegments  = n_eta_segments;
95   fAddHistograms = kFALSE;
96   fDoIterative   = kFALSE; 
97   fWriteDigits   = kFALSE;
98   fUse8bits      = bit8;
99   fVersion       = tv;
100 }
101
102 AliL3Hough::~AliL3Hough()
103 {
104   //dtor
105
106   CleanUp();
107   if(fMerger)
108     delete fMerger;
109   if(fInterMerger)
110     delete fInterMerger;
111   if(fPeakFinder)
112     delete fPeakFinder;
113   if(fGlobalMerger)
114     delete fGlobalMerger;
115   if(fBenchmark)
116     delete fBenchmark;
117 }
118
119 void AliL3Hough::CleanUp()
120 {
121   //Cleanup memory
122   
123   for(Int_t i=0; i<fNPatches; i++)
124     {
125       if(fTracks[i]) delete fTracks[i];
126       if(fEval[i]) delete fEval[i];
127       if(fHoughTransformer[i]) delete fHoughTransformer[i];
128       if(fMemHandler[i]) delete fMemHandler[i];
129     }
130   
131   if(fTracks) delete [] fTracks;
132   if(fEval) delete [] fEval;
133   if(fHoughTransformer) delete [] fHoughTransformer;
134   if(fMemHandler) delete [] fMemHandler;
135   
136 }
137
138 void AliL3Hough::Init(Char_t *path,Bool_t binary,Int_t n_eta_segments,Bool_t bit8,Int_t tv)
139 {
140   fBinary = binary;
141   strcpy(fPath,path);
142   fNEtaSegments = n_eta_segments;
143   fWriteDigits  = kFALSE;
144   fUse8bits     = bit8;
145   fVersion      = tv;
146
147   Init(); //do the rest
148 }
149
150 void AliL3Hough::Init(Bool_t doit, Bool_t addhists)
151 {
152   fDoIterative   = doit; 
153   fAddHistograms = addhists;
154
155   AliL3Transform::Init(fPath,!fBinary);
156   fNPatches = AliL3Transform::GetNPatches();
157   
158   fHoughTransformer = new AliL3HoughBaseTransformer*[fNPatches];
159   fMemHandler = new AliL3MemHandler*[fNPatches];
160
161   fTracks = new AliL3TrackArray*[fNPatches];
162   fEval = new AliL3HoughEval*[fNPatches];
163
164   for(Int_t i=0; i<fNPatches; i++)
165     {
166       switch (fVersion){ //choose Transformer
167       case 1: 
168         fHoughTransformer[i] = new AliL3HoughTransformerLUT(0,i,fNEtaSegments);
169         break;
170       case 2:
171         fHoughTransformer[i] = new AliL3HoughClusterTransformer(0,i,fNEtaSegments);
172         break;
173       case 3:
174         fHoughTransformer[i] = new AliL3HoughTransformerVhdl(0,i,fNEtaSegments,fNSaveIterations);
175         break;
176       default:
177         fHoughTransformer[i] = new AliL3HoughTransformer(0,i,fNEtaSegments,kFALSE,kFALSE);
178       }
179
180       fHoughTransformer[i]->CreateHistograms(fNBinX,fLowPt,fNBinY,-fPhi,fPhi);
181       fHoughTransformer[i]->SetLowerThreshold(fThreshold);
182  
183       LOG(AliL3Log::kInformational,"AliL3Hough::Init","Version")
184         <<"Initializing Hough transformer version "<<fVersion<<ENDLOG;
185       
186       fEval[i] = new AliL3HoughEval();
187       fTracks[i] = new AliL3TrackArray("AliL3HoughTrack");
188       if(fUse8bits)
189         fMemHandler[i] = new AliL3DataHandler();
190       else
191 #ifdef use_aliroot
192         {
193           fMemHandler[i] = new AliL3FileHandler();
194           if(!fBinary)
195             {
196               Char_t filename[1024];
197               sprintf(filename,"%s/digitfile.root",fPath);
198               fMemHandler[i]->SetAliInput(filename);
199             }
200         }
201 #else
202       fMemHandler[i] = new AliL3MemHandler();
203 #endif
204     }
205
206   fPeakFinder = new AliL3HoughMaxFinder("KappaPhi",1000);
207   fMerger = new AliL3HoughMerger(fNPatches);
208   fInterMerger = new AliL3HoughIntMerger();
209   fGlobalMerger = 0;
210   fBenchmark = new AliL3Benchmark();
211
212 }
213
214 void AliL3Hough::DoBench(Char_t *name)
215 {
216   fBenchmark->Analyze(name);
217 }
218
219 void AliL3Hough::Process(Int_t minslice,Int_t maxslice)
220 {
221   //Process all slices [minslice,maxslice].
222   fGlobalMerger = new AliL3HoughGlobalMerger(minslice,maxslice);
223   
224   for(Int_t i=minslice; i<=maxslice; i++)
225     {
226       ReadData(i);
227       Transform();
228       if(fAddHistograms)
229         AddAllHistograms();
230       FindTrackCandidates();
231       //Evaluate();
232       //fGlobalMerger->FillTracks(fTracks[0],i);
233     }
234 }
235
236 void AliL3Hough::ReadData(Int_t slice,Int_t eventnr)
237 {
238   //Read data from files, binary or root.
239   
240   fCurrentSlice = slice;
241   for(Int_t i=0; i<fNPatches; i++)
242     {
243       fMemHandler[i]->Free();
244       UInt_t ndigits=0;
245       AliL3DigitRowData *digits =0;
246       Char_t name[256];
247       fMemHandler[i]->Init(slice,i);
248       if(fBinary)//take input data from binary files
249         {
250           if(fUse8bits)
251             sprintf(name,"%sdigits_c8_%d_%d.raw",fPath,slice,i);
252           else
253             sprintf(name,"%sdigits_%d_%d.raw",fPath,slice,i);
254
255           fMemHandler[i]->SetBinaryInput(name);
256           digits = (AliL3DigitRowData *)fMemHandler[i]->CompBinary2Memory(ndigits);
257           fMemHandler[i]->CloseBinaryInput();
258         }
259       else //read data from root file
260         {
261 #ifdef use_aliroot
262           digits=(AliL3DigitRowData *)fMemHandler[i]->AliAltroDigits2Memory(ndigits,eventnr);
263           fMemHandler[i]->FreeDigitsTree();
264 #else
265           cerr<<"You cannot read from rootfile now"<<endl;
266 #endif
267         }
268
269       //set input data and init transformer
270       fHoughTransformer[i]->SetInputData(ndigits,digits);
271       fHoughTransformer[i]->Init(slice,i,fNEtaSegments);
272     }
273 }
274
275 void AliL3Hough::Transform(Int_t *row_range)
276 {
277   //Transform all data given to the transformer within the given slice
278   //(after ReadData(slice))
279   
280   Double_t initTime,cpuTime;
281   initTime = GetCpuTime();
282   for(Int_t i=0; i<fNPatches; i++)
283     {
284       fHoughTransformer[i]->Reset();//Reset the histograms
285       fBenchmark->Start("Hough Transform");
286       if(!row_range)
287         fHoughTransformer[i]->TransformCircle();
288       else
289         fHoughTransformer[i]->TransformCircleC(row_range,1);
290       fBenchmark->Stop("Hough Transform");
291     }
292   cpuTime = GetCpuTime() - initTime;
293   LOG(AliL3Log::kInformational,"AliL3Hough::Transform()","Timing")
294     <<"Transform done in average per patch of "<<cpuTime*1000/fNPatches<<" ms"<<ENDLOG;
295 }
296
297 void AliL3Hough::MergePatches()
298 {
299   if(fAddHistograms) //Nothing to merge here
300     return;
301   fMerger->MergePatches(kTRUE);
302 }
303
304 void AliL3Hough::MergeInternally()
305 {
306   if(fAddHistograms)
307     fInterMerger->FillTracks(fTracks[0]);
308   else
309     fInterMerger->FillTracks(fMerger->GetOutTracks());
310   
311   fInterMerger->MMerge();
312 }
313
314 void AliL3Hough::ProcessSliceIter()
315 {
316   //Process current slice (after ReadData(slice)) iteratively.
317   
318   if(!fAddHistograms)
319     {
320       for(Int_t i=0; i<fNPatches; i++)
321         {
322           ProcessPatchIter(i);
323           fMerger->FillTracks(fTracks[i],i); //Copy tracks to merger
324         }
325     }
326   else
327     {
328       for(Int_t i=0; i<10; i++)
329         {
330           Transform();
331           AddAllHistograms();
332           InitEvaluate();
333           AliL3HoughBaseTransformer *tr = fHoughTransformer[0];
334           for(Int_t j=0; j<fNEtaSegments; j++)
335             {
336               AliL3Histogram *hist = tr->GetHistogram(j);
337               if(hist->GetNEntries()==0) continue;
338               fPeakFinder->Reset();
339               fPeakFinder->SetHistogram(hist);
340               fPeakFinder->FindAbsMaxima();
341               AliL3HoughTrack *track = (AliL3HoughTrack*)fTracks[0]->NextTrack();
342               track->SetTrackParameters(fPeakFinder->GetXPeak(0),fPeakFinder->GetYPeak(0),fPeakFinder->GetWeight(0));
343               track->SetEtaIndex(j);
344               track->SetEta(tr->GetEta(j,fCurrentSlice));
345               for(Int_t k=0; k<fNPatches; k++)
346                 {
347                   fEval[i]->SetNumOfPadsToLook(2);
348                   fEval[i]->SetNumOfRowsToMiss(2);
349                   fEval[i]->RemoveFoundTracks();
350                   /*
351                   Int_t nrows=0;
352                   if(!fEval[i]->LookInsideRoad(track,nrows))
353                     {
354                       fTracks[0]->Remove(fTracks[0]->GetNTracks()-1);
355                       fTracks[0]->Compress();
356                     }
357                   */
358                 }
359             }
360           
361         }
362       
363     }
364 }
365
366 void AliL3Hough::ProcessPatchIter(Int_t patch)
367 {
368   //Process patch in a iterative way. 
369   //transform + peakfinding + evaluation + transform +...
370
371   Int_t num_of_tries = 5;
372   AliL3HoughBaseTransformer *tr = fHoughTransformer[patch];
373   AliL3TrackArray *tracks = fTracks[patch];
374   tracks->Reset();
375   AliL3HoughEval *ev = fEval[patch];
376   ev->InitTransformer(tr);
377   //ev->RemoveFoundTracks();
378   ev->SetNumOfRowsToMiss(3);
379   ev->SetNumOfPadsToLook(2);
380   AliL3Histogram *hist;
381   for(Int_t t=0; t<num_of_tries; t++)
382     {
383       tr->Reset();
384       tr->TransformCircle();
385       for(Int_t i=0; i<fNEtaSegments; i++)
386         {
387           hist = tr->GetHistogram(i);
388           if(hist->GetNEntries()==0) continue;
389           fPeakFinder->Reset();
390           fPeakFinder->SetHistogram(hist);
391           fPeakFinder->FindAbsMaxima();
392           //fPeakFinder->FindPeak1();
393           AliL3HoughTrack *track = (AliL3HoughTrack*)tracks->NextTrack();
394           track->SetTrackParameters(fPeakFinder->GetXPeak(0),fPeakFinder->GetYPeak(0),fPeakFinder->GetWeight(0));
395           track->SetEtaIndex(i);
396           track->SetEta(tr->GetEta(i,fCurrentSlice));
397           /*
398           Int_t nrows=0;
399           if(!ev->LookInsideRoad(track,nrows))
400             {   
401               tracks->Remove(tracks->GetNTracks()-1);
402               tracks->Compress();
403             }
404           */
405         }
406     }
407   fTracks[0]->QSort();
408   LOG(AliL3Log::kInformational,"AliL3Hough::ProcessPatch","NTracks")
409     <<AliL3Log::kDec<<"Found "<<tracks->GetNTracks()<<" tracks in patch "<<patch<<ENDLOG;
410 }
411
412 void AliL3Hough::AddAllHistograms()
413 {
414   //Add the histograms within one etaslice.
415   //Resulting histogram are in patch=0.
416
417   Double_t initTime,cpuTime;
418   initTime = GetCpuTime();
419   fBenchmark->Start("Add Histograms");
420   for(Int_t i=0; i<fNEtaSegments; i++)
421     {
422       AliL3Histogram *hist0 = fHoughTransformer[0]->GetHistogram(i);
423       for(Int_t j=1; j<fNPatches; j++)
424         {
425           AliL3Histogram *hist = fHoughTransformer[j]->GetHistogram(i);
426           hist0->Add(hist);
427         }
428     }
429   fBenchmark->Stop("Add Histograms");
430   fAddHistograms = kTRUE;
431   cpuTime = GetCpuTime() - initTime;
432   LOG(AliL3Log::kInformational,"AliL3Hough::AddAllHistograms()","Timing")
433     <<"Adding histograms in "<<cpuTime*1000<<" ms"<<ENDLOG;
434 }
435
436 void AliL3Hough::FindTrackCandidates()
437 {
438   //Look for peaks in histograms, and find the track candidates
439   
440   Int_t n_patches;
441   if(fAddHistograms)
442     n_patches = 1; //Histograms have been added.
443   else
444     n_patches = fNPatches;
445   
446   Double_t initTime,cpuTime;
447   initTime = GetCpuTime();
448   fBenchmark->Start("Find Maxima");
449   for(Int_t i=0; i<n_patches; i++)
450     {
451       AliL3HoughBaseTransformer *tr = fHoughTransformer[i];
452       fTracks[i]->Reset();
453       
454       for(Int_t j=0; j<fNEtaSegments; j++)
455         {
456           AliL3Histogram *hist = tr->GetHistogram(j);
457           if(hist->GetNEntries()==0) continue;
458           fPeakFinder->Reset();
459           fPeakFinder->SetHistogram(hist);
460           //fPeakFinder->FindPeak1(2,1);
461           fPeakFinder->FindMaxima(0,0); //Simple maxima finder
462           //fPeakFinder->FindAbsMaxima();
463           for(Int_t k=0; k<fPeakFinder->GetEntries(); k++)
464             {
465               if(fPeakFinder->GetWeight(k) == 0 || fPeakFinder->GetWeight(k) < fPeakThreshold) continue;
466               AliL3HoughTrack *track = (AliL3HoughTrack*)fTracks[i]->NextTrack();
467               track->SetTrackParameters(fPeakFinder->GetXPeak(k),fPeakFinder->GetYPeak(k),fPeakFinder->GetWeight(k));
468               track->SetEtaIndex(j);
469               track->SetEta(tr->GetEta(j,fCurrentSlice));
470               track->SetRowRange(AliL3Transform::GetFirstRow(0),AliL3Transform::GetLastRow(5));
471             }
472         }
473       fTracks[i]->QSort();
474     }
475   fBenchmark->Stop("Find Maxima");
476   cpuTime = GetCpuTime() - initTime;
477   LOG(AliL3Log::kInformational,"AliL3Hough::FindTrackCandidates()","Timing")
478     <<"Maxima finding done in "<<cpuTime*1000<<" ms"<<ENDLOG;
479 }
480
481 void AliL3Hough::InitEvaluate()
482 {
483   //Pass the transformer objects to the AliL3HoughEval objects:
484   //This will provide the evaluation objects with all the necessary
485   //data and parameters it needs.
486   
487   for(Int_t i=0; i<fNPatches; i++) 
488     fEval[i]->InitTransformer(fHoughTransformer[i]);
489 }
490
491 Int_t AliL3Hough::Evaluate(Int_t road_width,Int_t nrowstomiss)
492 {
493   //Evaluate the tracks, by looking along the road in the raw data.
494   //If track does not cross all padrows - rows2miss, it is removed from the arrray.
495   //If histograms were not added, the check is done locally in patch,
496   //meaning that nrowstomiss is the number of padrows the road can miss with respect
497   //to the number of rows in the patch.
498   //If the histograms were added, the comparison is done globally in the _slice_, 
499   //meaing that nrowstomiss is the number of padrows the road can miss with
500   //respect to the total number of padrows in the slice.
501   //
502   //Return value = number of tracks which were removed (only in case of fAddHistograms)
503   
504   if(!fTracks[0])
505     {
506       LOG(AliL3Log::kError,"AliL3Hough::Evaluate","Track Array")
507         <<"No tracks to work with..."<<ENDLOG;
508       return 0;
509     }
510   
511   Int_t removed_tracks=0;
512   AliL3TrackArray *tracks=0;
513
514   if(fAddHistograms)
515     {
516       tracks = fTracks[0];
517       for(Int_t i=0; i<tracks->GetNTracks(); i++)
518         {
519           AliL3Track *track = tracks->GetCheckedTrack(i);
520           if(!track) continue;
521           track->SetNHits(0);
522         }
523     }
524   
525   for(Int_t i=0; i<fNPatches; i++)
526     EvaluatePatch(i,road_width,nrowstomiss);
527   
528   //Here we check the tracks globally; 
529   //how many good rows (padrows with signal) 
530   //did it cross in the slice
531   if(fAddHistograms) 
532     {
533       for(Int_t j=0; j<tracks->GetNTracks(); j++)
534         {
535           AliL3HoughTrack *track = (AliL3HoughTrack*)tracks->GetCheckedTrack(j);
536           
537           if(track->GetNHits() < AliL3Transform::GetNRows() - nrowstomiss)
538             {
539               tracks->Remove(j);
540               removed_tracks++;
541             }
542         }
543       tracks->Compress();
544       tracks->QSort();
545     }
546     
547   return removed_tracks;
548 }
549
550 void AliL3Hough::EvaluatePatch(Int_t i,Int_t road_width,Int_t nrowstomiss)
551 {
552   //Evaluate patch i.
553   
554   fEval[i]->InitTransformer(fHoughTransformer[i]);
555   fEval[i]->SetNumOfPadsToLook(road_width);
556   fEval[i]->SetNumOfRowsToMiss(nrowstomiss);
557   //fEval[i]->RemoveFoundTracks();
558   
559   AliL3TrackArray *tracks=0;
560   
561   if(!fAddHistograms)
562     tracks = fTracks[i];
563   else
564     tracks = fTracks[0];
565   
566   Int_t nrows=0;
567   for(Int_t j=0; j<tracks->GetNTracks(); j++)
568     {
569       AliL3HoughTrack *track = (AliL3HoughTrack*)tracks->GetCheckedTrack(j);
570       if(!track)
571         {
572           LOG(AliL3Log::kWarning,"AliL3Hough::EvaluatePatch","Track array")
573             <<"Track object missing!"<<ENDLOG;
574           continue;
575         } 
576       nrows=0;
577       Int_t rowrange[2] = {AliL3Transform::GetFirstRow(i),AliL3Transform::GetLastRow(i)};
578       Bool_t result = fEval[i]->LookInsideRoad(track,nrows,rowrange);
579       if(fAddHistograms)
580         {
581           Int_t pre=track->GetNHits();
582           track->SetNHits(pre+nrows);
583         }
584       else//the track crossed too few good padrows (padrows with signal) in the patch, so remove it
585         {
586           if(result == kFALSE)
587             tracks->Remove(j);
588         }
589     }
590   
591   tracks->Compress();
592
593 }
594
595 void AliL3Hough::MergeEtaSlices()
596 {
597   //Merge tracks found in neighbouring eta slices.
598   //Removes the track with the lower weight.
599   
600   fBenchmark->Start("Merge Eta-slices");
601   AliL3TrackArray *tracks = fTracks[0];
602   if(!tracks)
603     {
604       cerr<<"AliL3Hough::MergeEtaSlices : No tracks "<<endl;
605       return;
606     }
607   for(Int_t j=0; j<tracks->GetNTracks(); j++)
608     {
609       AliL3HoughTrack *track1 = (AliL3HoughTrack*)tracks->GetCheckedTrack(j);
610       if(!track1) continue;
611       for(Int_t k=j+1; k<tracks->GetNTracks(); k++)
612         {
613           AliL3HoughTrack *track2 = (AliL3HoughTrack*)tracks->GetCheckedTrack(k);
614           if(!track2) continue;
615           if(abs(track1->GetEtaIndex() - track2->GetEtaIndex()) != 1) continue;
616           if(fabs(track1->GetKappa()-track2->GetKappa()) < 0.006 && 
617              fabs(track1->GetPsi()- track2->GetPsi()) < 0.1)
618             {
619               //cout<<"Merging track in slices "<<track1->GetEtaIndex()<<" "<<track2->GetEtaIndex()<<endl;
620               if(track1->GetWeight() > track2->GetWeight())
621                 tracks->Remove(k);
622               else
623                 tracks->Remove(j);
624             }
625         }
626     }
627   fBenchmark->Stop("Merge Eta-slices");
628   tracks->Compress();
629 }
630
631 void AliL3Hough::WriteTracks(Int_t slice,Char_t *path)
632 {
633   
634   AliL3MemHandler mem;
635   Char_t fname[100];
636   if(fAddHistograms)
637     {
638       sprintf(fname,"%s/tracks_ho_%d.raw",path,slice);
639       mem.SetBinaryOutput(fname);
640       mem.TrackArray2Binary(fTracks[0]);
641       mem.CloseBinaryOutput();
642     }
643   else 
644     {
645       for(Int_t i=0; i<fNPatches; i++)
646         {
647           sprintf(fname,"%s/tracks_ho_%d_%d.raw",path,slice,i);
648           mem.SetBinaryOutput(fname);
649           mem.TrackArray2Binary(fTracks[i]);
650           mem.CloseBinaryOutput();
651         }
652     }
653 }
654
655 void AliL3Hough::WriteDigits(Char_t *outfile)
656 {
657 #ifdef use_aliroot  
658   //Write the current data to a new rootfile.
659
660   for(Int_t i=0; i<fNPatches; i++)
661     {
662       AliL3DigitRowData *tempPt = (AliL3DigitRowData*)fHoughTransformer[i]->GetDataPointer();
663       fMemHandler[i]->AliDigits2RootFile(tempPt,outfile);
664     }
665 #else
666   cerr<<"AliL3Hough::WriteDigits : You need to compile with AliROOT!"<<endl;
667   return;
668 #endif  
669 }
670
671 Double_t AliL3Hough::GetCpuTime()
672 {
673   //Return the Cputime in seconds.
674  struct timeval tv;
675  gettimeofday( &tv, NULL );
676  return tv.tv_sec+(((Double_t)tv.tv_usec)/1000000.);
677  //return (Double_t)(clock()) / CLOCKS_PER_SEC;
678 }
679