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