Forgot to remove #include in previous checkin
[u/mrichter/AliRoot.git] / HLT / hough / AliL3Hough.cxx
1 //$Id$
2
3 // Author: Anders Vestbo <mailto:vestbo@fi.uib.no>
4 //*-- Copyright &copy ASV 
5
6
7 #include <string.h>
8
9 #include "AliL3HoughMerger.h"
10 #include "AliL3HoughIntMerger.h"
11 #include "AliL3HoughGlobalMerger.h"
12 #include "AliL3Logging.h"
13 #include "AliL3Histogram.h"
14 #include "AliL3Hough.h"
15 #include "AliL3HoughTransformer.h"
16 #include "AliL3HoughMaxFinder.h"
17 #ifdef use_aliroot
18 #include "AliL3FileHandler.h"
19 #else
20 #include "AliL3MemHandler.h"
21 #endif
22 #include "AliL3DigitData.h"
23 #include "AliL3HoughEval.h"
24 #include "AliL3Transform.h"
25 #include "AliL3TrackArray.h"
26 #include "AliL3HoughTrack.h"
27
28
29 //_____________________________________________________________
30 // AliL3Hough
31 //
32 // Interface class for the Hough transform
33 //
34 // Example how to use:
35 //
36 // AliL3Hough *hough = new AliL3Hough(path,kTRUE,NumberOfEtaSegments);
37 // hough->ReadData(slice);
38 // hough->Transform();
39 // hough->FindTrackCandidates();
40 // 
41 // AliL3TrackArray *tracks = hough->GetTracks(patch);
42
43 ClassImp(AliL3Hough)
44
45 AliL3Hough::AliL3Hough()
46 {
47   //Constructor
48   
49   fBinary = kFALSE;
50   fNEtaSegments = 0;
51   fAddHistograms = kFALSE;
52   fDoIterative = kFALSE; 
53   fWriteDigits=kFALSE;
54   fNPatches=0;
55   fMemHandler = 0;
56   fHoughTransformer = 0;
57   fEval = 0;
58   fPeakFinder = 0;
59   fTracks = 0;
60   fMerger = 0;
61   fInterMerger = 0;
62   fGlobalMerger = 0;
63 }
64
65
66 AliL3Hough::AliL3Hough(Char_t *path,Bool_t binary,Int_t n_eta_segments)
67 {
68   //Default ctor.
69
70   fBinary = binary;
71   strcpy(fPath,path);
72   fNEtaSegments = n_eta_segments;
73   fAddHistograms = kFALSE;
74   fDoIterative = kFALSE; 
75   fWriteDigits = kFALSE;
76   Init();
77 }
78
79
80 AliL3Hough::~AliL3Hough()
81 {
82   //dtor
83
84   CleanUp();
85   if(fMerger)
86     delete fMerger;
87   if(fInterMerger)
88     delete fInterMerger;
89   if(fPeakFinder)
90     delete fPeakFinder;
91   if(fGlobalMerger)
92     delete fGlobalMerger;
93 }
94
95 void AliL3Hough::CleanUp()
96 {
97   //Cleanup memory
98   
99   for(Int_t i=0; i<fNPatches; i++)
100     {
101       if(fTracks[i]) delete fTracks[i];
102       if(fEval[i]) delete fEval[i];
103       if(fHoughTransformer[i]) delete fHoughTransformer[i];
104       if(fMemHandler[i]) delete fMemHandler[i];
105     }
106   
107   /*    
108         if(fTracks) delete [] fTracks;
109         if(fEval) delete [] fEval;
110         if(fHoughTransformer) delete [] fHoughTransformer;
111         if(fMemHandler) delete [] fMemHandler;
112   */
113 }
114
115 void AliL3Hough::Init()
116 {
117   AliL3Transform::Init(fPath);
118   fPeakThreshold = 0;
119   fNPatches = AliL3Transform::GetNPatches();
120   fHoughTransformer = new AliL3HoughBaseTransformer*[fNPatches];
121 #ifdef use_aliroot
122   fMemHandler = new AliL3FileHandler*[fNPatches];
123 #else
124   fMemHandler = new AliL3MemHandler*[fNPatches];
125 #endif
126   fTracks = new AliL3TrackArray*[fNPatches];
127   fEval = new AliL3HoughEval*[fNPatches];
128   for(Int_t i=0; i<fNPatches; i++)
129     {
130       fHoughTransformer[i] = new AliL3HoughTransformer(1,i,fNEtaSegments);
131       //fHoughTransformer[i]->CreateHistograms(64,-0.003,0.003,64,-0.26,0.26);
132       fHoughTransformer[i]->CreateHistograms(64,0.1,64,-30,30);
133       fHoughTransformer[i]->SetThreshold(3);
134       fEval[i] = new AliL3HoughEval();
135       fTracks[i] = new AliL3TrackArray("AliL3HoughTrack");
136 #ifdef use_aliroot
137       fMemHandler[i] = new AliL3FileHandler();
138       if(!fBinary)
139         fMemHandler[i]->SetAliInput(fPath);
140 #else
141       fMemHandler[i] = new AliL3MemHandler();
142 #endif
143       
144     }
145   fPeakFinder = new AliL3HoughMaxFinder("KappaPhi");
146   fMerger = new AliL3HoughMerger(fNPatches);
147   fInterMerger = new AliL3HoughIntMerger();
148   fGlobalMerger = 0;
149 }
150
151 void AliL3Hough::Process(Int_t minslice,Int_t maxslice)
152 {
153   //Process all slices [minslice,maxslice].
154   fGlobalMerger = new AliL3HoughGlobalMerger(minslice,maxslice);
155   
156   for(Int_t i=minslice; i<=maxslice; i++)
157     {
158       ReadData(i);
159       Transform();
160       if(fAddHistograms)
161         AddAllHistograms();
162       FindTrackCandidates();
163       Evaluate();
164       fGlobalMerger->FillTracks(fTracks[0],i);
165     }
166   
167   
168 }
169
170 void AliL3Hough::ReadData(Int_t slice)
171 {
172   //Read data from files, binary or root.
173
174   for(Int_t i=0; i<fNPatches; i++)
175     {
176       fMemHandler[i]->Free();
177       UInt_t ndigits=0;
178       AliL3DigitRowData *digits =0;
179       Char_t name[256];
180       if(fBinary)//take input data from binary files
181         {
182           sprintf(name,"%sdigits_%d_%d.raw",fPath,slice,i);
183           fMemHandler[i]->SetBinaryInput(name);
184           digits = (AliL3DigitRowData *)fMemHandler[i]->CompBinary2Memory(ndigits);
185           fMemHandler[i]->CloseBinaryInput();
186         }
187       else //read data from root file
188         {
189 #ifdef use_aliroot
190           const Int_t rows[2] = {AliL3Transform::GetFirstRow(i),AliL3Transform::GetLastRow(i)};
191           fMemHandler[i]->Init(slice,i,rows);
192           digits=(AliL3DigitRowData *)fMemHandler[i]->AliDigits2Memory(ndigits); 
193 #else
194           cerr<<"You cannot read from rootfile now"<<endl;
195 #endif
196         }
197       fHoughTransformer[i]->SetInputData(ndigits,digits);
198     }
199 }
200
201 void AliL3Hough::Transform(Int_t row_range)
202 {
203   //Transform all data given to the transformer within the given slice
204   //(after ReadData(slice))
205
206   for(Int_t i=0; i<fNPatches; i++)
207     {
208       fHoughTransformer[i]->Reset();//Reset the histograms
209       if(row_range < 0)
210         fHoughTransformer[i]->TransformCircle();
211       else
212         fHoughTransformer[i]->TransformCircleC(row_range);
213     }
214 }
215
216 void AliL3Hough::MergePatches()
217 {
218   if(fAddHistograms) //Nothing to merge here
219     return;
220   fMerger->MergePatches(kTRUE);
221 }
222
223 void AliL3Hough::MergeInternally()
224 {
225   if(fAddHistograms)
226     fInterMerger->FillTracks(fTracks[0]);
227   else
228     fInterMerger->FillTracks(fMerger->GetOutTracks());
229   
230   fInterMerger->MMerge();
231 }
232
233 void AliL3Hough::ProcessSliceIter()
234 {
235   //Process current slice (after ReadData(slice)) iteratively.
236   
237   for(Int_t i=0; i<fNPatches; i++)
238     {
239       ProcessPatchIter(i);
240       fMerger->FillTracks(fTracks[i],i); //Copy tracks to merger
241     }
242   
243 }
244
245 void AliL3Hough::ProcessPatchIter(Int_t patch)
246 {
247   //Process patch in a iterative way. 
248   //transform + peakfinding + evaluation + transform +...
249
250   Int_t num_of_tries = 10;
251   AliL3HoughBaseTransformer *tr = fHoughTransformer[patch];
252   AliL3TrackArray *tracks = fTracks[patch];
253   tracks->Reset();
254   AliL3HoughEval *ev = fEval[patch];
255   ev->InitTransformer(tr);
256   ev->RemoveFoundTracks();
257   ev->SetNumOfRowsToMiss(2);
258   ev->SetNumOfPadsToLook(2);
259   AliL3Histogram *hist;
260   for(Int_t t=0; t<num_of_tries; t++)
261     {
262       tr->Reset();
263       tr->TransformCircle();
264       for(Int_t i=0; i<fNEtaSegments; i++)
265         {
266           hist = tr->GetHistogram(i);
267           if(hist->GetNEntries()==0) continue;
268           fPeakFinder->SetHistogram(hist);
269           //Int_t n=1;
270           Float_t x,y;
271           //fPeakFinder->FindAbsMaxima(*x,*y);
272           fPeakFinder->FindPeak(3,0.95,5,x,y);
273           AliL3HoughTrack *track = (AliL3HoughTrack*)tracks->NextTrack();
274           track->SetTrackParameters(x,y,1);
275           if(!ev->LookInsideRoad(track,i))
276             {   
277               tracks->Remove(tracks->GetNTracks()-1);
278               tracks->Compress();
279             }
280         }
281     }
282   LOG(AliL3Log::kInformational,"AliL3Hough::ProcessPatch","NTracks")
283     <<AliL3Log::kDec<<"Found "<<tracks->GetNTracks()<<" tracks in patch "<<patch<<ENDLOG;
284 }
285
286
287 void AliL3Hough::AddAllHistograms()
288 {
289   //Add the histograms within one etaslice.
290   //Resulting histogram are in patch=0.
291
292   for(Int_t i=0; i<fNEtaSegments; i++)
293     {
294       AliL3Histogram *hist0 = fHoughTransformer[0]->GetHistogram(i);
295       for(Int_t j=1; j<fNPatches; j++)
296         {
297           AliL3Histogram *hist = fHoughTransformer[j]->GetHistogram(i);
298           hist0->Add(hist);
299         }
300     }
301   fAddHistograms = kTRUE;
302 }
303
304 void AliL3Hough::FindTrackCandidates()
305 {
306   //Look for peaks in histograms, and find the track candidates
307   
308   Int_t n_patches;
309   if(fAddHistograms)
310     n_patches = 1; //Histograms has been added.
311   else
312     n_patches = fNPatches;
313
314   
315   for(Int_t i=0; i<n_patches; i++)
316     {
317       AliL3HoughBaseTransformer *tr = fHoughTransformer[i];
318       Double_t eta_slice = (tr->GetEtaMax()-tr->GetEtaMin()/tr->GetNEtaSegments());
319       fTracks[i]->Reset();
320       for(Int_t j=0; j<fNEtaSegments; j++)
321         {
322           AliL3Histogram *hist = tr->GetHistogram(j);
323           if(hist->GetNEntries()==0) continue;
324           fPeakFinder->SetHistogram(hist);
325           fPeakFinder->SetThreshold(fPeakThreshold);
326           Int_t n=20;
327           Float_t x[n];
328           Float_t y[n];
329           Int_t weight[n];
330           //fPeakFinder->FindPeak1(x,y,weight,n,2,1);
331           fPeakFinder->FindMaxima(x,y,weight,n);
332           for(Int_t k=0; k<n; k++)
333             {
334               if(weight[k] == 0) continue;
335               
336               AliL3HoughTrack *track = (AliL3HoughTrack*)fTracks[i]->NextTrack();
337               track->SetTrackParameters(x[k],y[k],weight[k]);
338               track->SetEtaIndex(j);
339               track->SetEta((Double_t)(j*eta_slice));
340               track->SetRowRange(AliL3Transform::GetFirstRow(0),AliL3Transform::GetLastRow(5));
341             }
342         }
343       fTracks[i]->QSort();
344     }
345 }
346
347 void AliL3Hough::Evaluate(Int_t road_width)
348 {
349   //Evaluate the tracks, by looking along the road in the raw data.
350   
351
352   if(!fTracks[0])
353     {
354       LOG(AliL3Log::kError,"AliL3Hough::Evaluate","Track Array")
355         <<"No tracks to work with..."<<ENDLOG;
356       return;
357     }
358   
359   printf("Number of tracks before evaluation %d\n",fTracks[0]->GetNTracks());
360   AliL3TrackArray *tracks;
361   for(Int_t i=0; i<fNPatches; i++)
362     {
363       fEval[i]->InitTransformer(fHoughTransformer[i]);
364       continue;
365       fEval[i]->SetNumOfRowsToMiss(2);
366       fEval[i]->SetNumOfPadsToLook(road_width);
367       if(fAddHistograms)
368         tracks = fTracks[0];
369       else
370         tracks = fTracks[i];
371       for(Int_t j=0; j<tracks->GetNTracks(); j++)
372         {
373           AliL3HoughTrack *track = (AliL3HoughTrack*)tracks->GetCheckedTrack(j);
374           if(!track)
375             {
376               LOG(AliL3Log::kWarning,"AliL3Hough::Evaluate","Track array")
377                 <<"Track object missing!"<<ENDLOG;
378               continue;
379             }
380           
381           if(!fEval[i]->LookInsideRoad(track,track->GetEtaIndex()))
382             tracks->Remove(j);
383           if(fAddHistograms)
384             track->SetRowRange(AliL3Transform::GetFirstRow(0),AliL3Transform::GetLastRow(5));//All rows included
385         }
386       tracks->Compress();
387       tracks->QSort(); //Sort the tracks according to weight
388       
389       if(!fAddHistograms)
390         fMerger->FillTracks(tracks,i); //Copy tracks to the track merger
391     }
392   
393 }
394
395 void AliL3Hough::EvaluateWithEta()
396 {
397   if(!fTracks[0])
398     {
399       printf("AliL3Hough::EvaluateWithEta: NO TRACKS\n");
400       return;
401     }
402   printf("Number of tracks before evaluation %d\n",fTracks[0]->GetNTracks());
403  
404   for(Int_t i=0; i<fNPatches; i++)
405     {
406       fEval[i]->InitTransformer(fHoughTransformer[i]);
407       fEval[i]->FindEta(fTracks[0]);
408     }
409   fMerger->FillTracks(fTracks[0],0);
410 }
411
412 void AliL3Hough::WriteTracks(Char_t *path)
413 {
414   AliL3MemHandler *mem = new AliL3MemHandler();
415   Char_t fname[100];
416   if(fAddHistograms)
417     {
418       sprintf(fname,"%s/tracks.raw",path);
419       mem->SetBinaryOutput(fname);
420       mem->TrackArray2Binary(fTracks[0]);
421       mem->CloseBinaryOutput();
422     }
423   else 
424     {
425       for(Int_t i=0; i<fNPatches; i++)
426         {
427           sprintf(fname,"%s/tracks_%d.raw",path,i);
428           mem->SetBinaryOutput(fname);
429           mem->TrackArray2Binary(fTracks[i]);
430           mem->CloseBinaryOutput();
431         }
432     }
433   delete mem;
434   
435 }
436 #ifdef use_aliroot
437 void AliL3Hough::WriteDigits(Char_t *outfile)
438 {
439   //Write the current data to a new rootfile.
440
441   for(Int_t i=0; i<fNPatches; i++)
442     {
443       AliL3DigitRowData *tempPt = (AliL3DigitRowData*)fHoughTransformer[i]->GetDataPointer();
444       fMemHandler[i]->AliDigits2RootFile(tempPt,outfile);
445     }
446   
447 }
448 #endif