New version of SPD raw-data reconstruction. The format now correponds to the actual...
[u/mrichter/AliRoot.git] / HLT / src / AliHLTGlobalMerger.cxx
1 // @(#) $Id$
2
3 // Author: Uli Frankenfeld <mailto:franken@fi.uib.no>
4 //*-- Copyright &copy ALICE HLT Group 
5
6 #include "AliHLTStandardIncludes.h"
7
8 #include "AliHLTLogging.h"
9 #include "AliHLTGlobalMerger.h"
10 #include "AliHLTTrack.h"
11 #include "AliHLTTransform.h"
12 #include "AliHLTTrackArray.h"
13
14 /** \class AliHLTGlobalMerger
15 <pre>
16 //_____________________________________________________________
17 // AliHLTGlobalMerger
18 //
19 // The L3 Slice merger
20 //
21 </pre>
22 */
23
24 #if __GNUC__ >= 3
25 using namespace std;
26 #endif
27
28 ClassImp(AliHLTGlobalMerger)
29
30 AliHLTGlobalMerger::AliHLTGlobalMerger()
31 {
32   //Default constructor. Use Setup to specify and setup the necessary parameters and arrays.
33   Is2Global(kTRUE);
34   SetParameter(0,0,0,0,0);
35   fNSlices=0;
36   fFirst=0;
37   fLast=0;
38 }
39
40
41 AliHLTGlobalMerger::~AliHLTGlobalMerger()
42 {
43   //Destructor
44 }
45
46 void AliHLTGlobalMerger::Setup(Int_t first,Int_t last)
47 {
48   //Used to setup the arrays and everything
49   
50   fNSlices = last-first+1;
51   fFirst = first;
52   fLast = last;
53   InitMerger(last-first+1);
54 }
55
56 void AliHLTGlobalMerger::InitSlice(Int_t slice)
57 {
58   // 
59   // Select Sector The following FillTracks call will 
60   // fill this Sector
61   //
62   fSlice = slice;
63   fCurrentTracks = fSlice - fFirst; 
64 }
65
66 Double_t AliHLTGlobalMerger::CheckTracks(AliHLTTrack *innertrack,AliHLTTrack *outertrack,Int_t slice)
67 {
68   //Compare the tracks by propagating the outermost track to the last and first point plane
69   //of the innermost track. This plane is defined by the padrow plane where these points
70   //are.
71   
72   if(innertrack->GetCharge()!=outertrack->GetCharge()) return -1;
73   
74   Float_t angle = 0;//perpendicular to the padrowplane (in local system)
75   AliHLTTransform::Local2GlobalAngle(&angle,slice);
76   Double_t dx[2],dy[2],dz[2];
77   Double_t diff =-1;
78   AliHLTTrack *tracks[2];
79   tracks[0] = innertrack;
80   tracks[1] = outertrack;
81   SortGlobalTracks(tracks,2);
82   innertrack = tracks[0]; 
83   outertrack = tracks[1];
84   
85   Float_t point[3];
86   
87   point[0]=innertrack->GetLastPointX();
88   point[1]=innertrack->GetLastPointY();
89   point[2]=innertrack->GetLastPointZ();
90   AliHLTTransform::Global2LocHLT(point,slice);
91   
92   outertrack->CalculateReferencePoint(angle,point[0]);//local x = global distance to padrowplane
93   if(!outertrack->IsPoint()) return diff;
94   dx[0] = fabs(outertrack->GetPointX()-innertrack->GetLastPointX());
95   dy[0] = fabs(outertrack->GetPointY()-innertrack->GetLastPointY());
96   dz[0] = fabs(outertrack->GetPointZ()-innertrack->GetLastPointZ());
97   
98   point[0]=innertrack->GetFirstPointX();
99   point[1]=innertrack->GetFirstPointY();
100   point[2]=innertrack->GetFirstPointZ();
101   AliHLTTransform::Global2LocHLT(point,slice);
102   
103   outertrack->CalculateReferencePoint(angle,point[0]);//local x = global distance to padrowplane
104   if(!outertrack->IsPoint()) return diff;
105   dx[1] = fabs(outertrack->GetPointX()-innertrack->GetFirstPointX());
106   dy[1] = fabs(outertrack->GetPointY()-innertrack->GetFirstPointY());
107   dz[1] = fabs(outertrack->GetPointZ()-innertrack->GetFirstPointZ());
108   
109   diff=0;//This was a tough bug to find....
110   for(Int_t i=0; i<2; i++)
111     diff += sqrt(dx[i]*dx[i] + dy[i]*dy[i] + dz[i]*dz[i]);
112   return diff;
113 }
114
115 void AliHLTGlobalMerger::SlowMerge(Char_t *path)
116 {
117   //Tuning of parameters. This matches _all_ tracks between two neighbouring
118   //slices, and merges the ones which are closest in space. The difference
119   //is written to a ntuppel, which can be used as a input for the SetParameters
120   //when using normal Merge function.
121   
122   
123   void* ntuple=GetNtuple();
124   AliHLTTrack *track[2];
125   AliHLTTrackArray *tout = GetOutTracks();
126   if(fNSlices<2)
127     {
128       LOG(AliHLTLog::kWarning,"AliHLTGlobalMerger::SlowMerge","Slice Number")
129         <<"Need more than one Slice!"<<ENDLOG;
130       return;
131     }
132   
133   for(Int_t i=0; i<fNSlices; i++)
134     {
135       //if(fNSlices!=18 && i+1 == fNSlices) continue; 
136       Int_t slice = fFirst + i;
137       AliHLTTrackArray *ttt0=GetInTracks(i);
138       Int_t slice2 = i+1;
139       //if(slice2==fNSlices) slice2 =0; 
140       
141       //Make sure slices are on the same side of the TPC
142       if(slice2 == 18) slice2=0;
143       else if(slice2 == 36) slice2=18;
144       AliHLTTrackArray *ttt1=GetInTracks(slice2);
145       //10 degrees -> the border of the slices in local coordinates
146       Float_t angle = AliHLTTransform::Pi()/18; 
147       AliHLTTransform::Local2GlobalAngle(&angle,slice);
148       
149       //In the two following cases, the angle is 2*pi, so set it back to 0 in order for
150       //the calculation of crossing points to be correct.
151       if(slice==17 || slice==35) 
152         angle=0;
153       if(i==0)
154         ttt0->QSort();
155       ttt1->QSort();
156       for(Int_t s0=0;s0<ttt0->GetNTracks();s0++)
157         {
158           AliHLTTrack *track0=ttt0->GetCheckedTrack(s0);
159           if(!track0) continue;
160           track0->CalculateHelix();
161           track0->CalculateEdgePoint(angle);
162           //      if(track0->IsPoint()) AddTrack(tout,track0);
163         }
164       for(Int_t s1=0;s1<ttt1->GetNTracks();s1++)
165         {
166           AliHLTTrack *track1=ttt1->GetCheckedTrack(s1);
167           if(!track1) continue;
168           track1->CalculateHelix();
169           track1->CalculateEdgePoint(angle);
170           //      if(track1->IsPoint())  AddTrack(tout,track1); 
171         }
172       Bool_t merge = kTRUE;
173       while(merge)
174         {
175           Int_t min0=-1,min1=-1;
176           Double_t min=5;
177           Int_t n0=ttt0->GetNTracks(),n1=ttt1->GetNTracks();
178           for(Int_t s0=0;s0<n0;s0++)
179             {
180               AliHLTTrack *track0=ttt0->GetCheckedTrack(s0);
181               if(!track0) continue;
182               if(!track0->IsPoint()) continue;
183               for(Int_t s1=0;s1<n1;s1++)
184                 {
185                   AliHLTTrack *track1=ttt1->GetCheckedTrack(s1);
186                   if(!track1) continue;
187                   if(!track1->IsPoint()) continue;
188                   
189                   //Double_t diff = TrackDiff(track0,track1,angle);
190                   Double_t diff = CheckTracks(track0,track1,slice);
191                   //PrintDiff(track0,track1);
192                   if(diff>=0&&diff<min)
193                     {
194                       min=diff;
195                       min0=s0;
196                       min1=s1;
197                     }
198                 }
199             }
200           if(min0>=0&&min1>=0)
201             {
202               AliHLTTrack *track0=ttt0->GetTrack(min0);
203               AliHLTTrack *track1=ttt1->GetTrack(min1);
204               track[0] = track0;
205               track[1] = track1;
206               SortGlobalTracks(track,2);
207               track1->CalculateEdgePoint((angle+AliHLTTransform::Pi()/9));
208               if(track1->IsPoint())//Check if the track will cross the boundary of yet another slice.
209                 MultiMerge(ttt1,track,2);
210               else
211                 MultiMerge(tout,track,2); 
212               track0->CalculateReferencePoint(angle);
213               track1->CalculateReferencePoint(angle);
214               //PrintDiff(track0,track1);
215               FillNtuple(ntuple,track0,track1);
216               ttt0->Remove(min0);
217               ttt1->Remove(min1);
218               
219             }
220           else merge = kFALSE;
221         }
222       ttt0->Compress();
223       ttt1->Compress();
224       LOG(AliHLTLog::kInformational,"AliHLTGlobalMerger::SlowMerge","Result")
225         <<AliHLTLog::kDec<<"Merged Tracks: "<<tout->GetNTracks()<<" at:"
226         <<angle<<ENDLOG;
227     }
228   Char_t fname[1024];
229   sprintf(fname,"%s/merge_parameters.root",path);
230   WriteNtuple(fname,ntuple);
231 }
232
233 void AliHLTGlobalMerger::Merge()
234 {
235   //Normal merging procedure. Matches tracks which are within limits
236   //set by SetParameters. Parameters can be tuned by SlowMerge.
237   
238   AliHLTTrack *track[2];
239   AliHLTTrackArray *tout = GetOutTracks();
240   if(fNSlices<2)
241     {
242       LOG(AliHLTLog::kWarning,"AliHLTGlobalMerger::Merge","Slice Number")
243         <<"Need more than one Slice!"<<ENDLOG;
244       return;
245     }
246   for(Int_t i=0; i<fNSlices; i++)
247     {
248       //if(fNSlices!=18 && i+1 == fNSlices) continue; 
249       Int_t slice = fFirst + i;
250       AliHLTTrackArray *ttt0=GetInTracks(i);
251       Int_t slice2 = i+1;
252       //if(slice2==fNSlices) slice2 =0;
253       
254       //Make sure slices are on the same side of the TPC
255       if(slice2 == 18) slice2=0;
256       else if(slice2 == 36) slice2=18;
257       AliHLTTrackArray *ttt1=GetInTracks(slice2);
258       //10 degrees -> the border of the slices in local coordinates
259       Float_t angle = AliHLTTransform::Pi()/18; 
260       AliHLTTransform::Local2GlobalAngle(&angle,slice);
261       
262       //In the two following cases, the angle is 2*pi, so set it back to 0 in order for
263       //the calculation of crossing points to be correct.
264       if(slice==17 || slice==35)
265         angle=0;
266       if(i==0)
267         ttt0->QSort();
268       ttt1->QSort();
269       Bool_t *ismatched0  = new Bool_t[ttt0->GetNTracks()];
270       Bool_t *ismatched1  = new Bool_t[ttt1->GetNTracks()];
271       Int_t n0=0,n1=0;
272       for(Int_t s0=0;s0<ttt0->GetNTracks();s0++)
273         {
274           ismatched0[s0]=kFALSE;
275           AliHLTTrack *track0=ttt0->GetCheckedTrack(s0);
276           if(!track0) continue;
277           track0->CalculateHelix();
278           track0->CalculateEdgePoint(angle);
279           if(track0->IsPoint()) 
280             {
281               n0++;
282               track0->CalculateReferencePoint(angle);
283             }
284         }
285       for(Int_t s1=0;s1<ttt1->GetNTracks();s1++)
286         {
287           ismatched1[s1]=kFALSE;
288           AliHLTTrack *track1=ttt1->GetCheckedTrack(s1);
289           if(!track1) continue;
290           track1->CalculateHelix();
291           track1->CalculateEdgePoint(angle);
292           if(track1->IsPoint()) 
293             {
294               n1++;
295               track1->CalculateReferencePoint(angle);
296             }
297         }
298       for(Int_t s0=0;s0<ttt0->GetNTracks();s0++)
299         {
300           if(ismatched0[s0]) continue;
301           AliHLTTrack *track0=ttt0->GetCheckedTrack(s0);
302           if(!track0) continue;
303           if(!track0->IsPoint()) continue;
304           for(Int_t s1=0;s1<ttt1->GetNTracks();s1++)
305             {
306               if(ismatched1[s1]) continue;
307               AliHLTTrack *track1=ttt1->GetCheckedTrack(s1);
308               if(!track1) continue;
309               if(!track1->IsPoint()) continue;
310               if(IsRTrack(track0,track1))
311                 {
312                   track[0] = track0;
313                   track[1] = track1;
314                   SortGlobalTracks(track,2);
315                   Double_t r0 = pow(track[0]->GetLastPointX(),2)+
316                     pow(track[0]->GetLastPointY(),2);
317                   Double_t r1 = pow(track[1]->GetFirstPointX(),2)+
318                     pow(track[1]->GetFirstPointY(),2);
319                   if(r0<r1)
320                     {
321                       MultiMerge(tout,track,2); 
322                       ismatched0[s0]=kTRUE;
323                       ismatched1[s1]=kTRUE;
324                       ttt0->Remove(s0);
325                       ttt1->Remove(s1);
326                       break;
327                       /*
328                         The track is merged, so we will _not_ look for more matches.
329                         Because there could easily be more matches, if a track is being
330                         split within the sector....
331                       */
332                     }
333                 }
334             }
335         }
336       LOG(AliHLTLog::kInformational,"AliHLTGlobalMerger::Merge","Result")
337         <<AliHLTLog::kDec<<"slice0: "<<n0<<" slice1: "<<n1
338         <<" Merged Tracks: "<<tout->GetNTracks()<<ENDLOG;
339       delete [] ismatched0;
340       delete [] ismatched1;
341     }
342 }