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