#101318: Patch for various problems in AliROOT
[u/mrichter/AliRoot.git] / TRD / AliTRDgtuTMU.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15
16 /* $Id: AliTRDgtuTMU.cxx 28397 2008-09-02 09:33:00Z cblume $ */
17
18 ////////////////////////////////////////////////////////////////////////////
19 //                                                                        //
20 //  Track Matching Unit (TMU) simulation                                  //
21 //                                                                        //
22 //  Author: J. Klein (Jochen.Klein@cern.ch)                               //
23 //                                                                        //
24 ////////////////////////////////////////////////////////////////////////////
25
26 #include "TTree.h"
27 #include "TList.h"
28 #include "TVectorD.h"
29 #include "TMath.h"
30
31 #include "AliESDEvent.h"
32 #include "AliESDTrdTrack.h"
33
34 #include "AliLog.h"
35 #include "AliTRDgeometry.h"
36 #include "AliTRDpadPlane.h"
37
38 #include "AliTRDgtuParam.h"
39 #include "AliTRDgtuTMU.h"
40 #include "AliTRDtrackGTU.h"
41
42 ClassImp(AliTRDgtuTMU)
43
44 AliTRDgtuTMU::AliTRDgtuTMU(Int_t stack, Int_t sector) :
45   TObject(),
46   fTracklets(0x0),
47   fTrackletsPostInput(0x0),
48   fZChannelTracklets(0x0),
49   fTracks(0x0),
50   fGtuParam(0x0),
51   fStack(-1),
52   fSector(-1)
53 {
54   // constructor which initializes the position information of the TMU
55
56   fGtuParam = AliTRDgtuParam::Instance();
57
58   // store tracklets per link
59   fTracklets = new TObjArray*[fGtuParam->GetNLinks()];
60   for (Int_t iLink = 0; iLink < fGtuParam->GetNLinks(); iLink++) {
61     fTracklets[iLink] = new TObjArray();
62   }
63
64   // tracklets after input units per layer
65   fTrackletsPostInput = new TObjArray*[fGtuParam->GetNLayers()];
66   fZChannelTracklets = new TList*[fGtuParam->GetNLayers()];
67   for (Int_t layer = 0;  layer <  fGtuParam->GetNLayers(); layer++) {
68     fZChannelTracklets[layer] = new TList[fGtuParam->GetNZChannels()];
69     fTrackletsPostInput[layer] = new TObjArray();
70   }
71
72   fTracks = new TList*[fGtuParam->GetNZChannels()];
73   for (Int_t zch = 0; zch < fGtuParam->GetNZChannels(); zch++) {
74       fTracks[zch] = new TList[fGtuParam->GetNRefLayers()];
75   }
76
77   if (stack > -1)
78       SetStack(stack);
79   if (sector > -1)
80       SetSector(sector);
81 }
82
83 AliTRDgtuTMU::~AliTRDgtuTMU()
84 {
85   // destructor
86
87   for (Int_t zch = 0; zch < fGtuParam->GetNZChannels(); zch++) {
88     delete [] fTracks[zch];
89   }
90   delete [] fTracks;
91   for (Int_t layer = 0; layer < fGtuParam->GetNLayers(); layer++) {
92     delete [] fZChannelTracklets[layer];
93     //    delete [] fTrackletsPostInput[layer];
94   }
95   delete [] fZChannelTracklets;
96   //  delete [] fTrackletsPostInput;
97
98   for (Int_t iLink = 0; iLink < fGtuParam->GetNLinks(); iLink++) {
99     delete fTracklets[iLink];
100   }
101   delete [] fTracklets;
102 }
103
104 Bool_t AliTRDgtuTMU::SetSector(Int_t sector)
105 {
106   // set the sector
107
108   if (sector > -1 && sector < fGtuParam->GetGeo()->Nsector() ) {
109     fSector = sector;
110     return kTRUE;
111   }
112
113   AliError(Form("Invalid sector given: %i", sector));
114   return kFALSE;
115 }
116
117 Bool_t AliTRDgtuTMU::SetStack(Int_t stack)
118 {
119   // Set the stack (necessary for tracking)
120
121   if (stack > -1 && stack < fGtuParam->GetGeo()->Nstack() ) {
122     fStack = stack;
123     return kTRUE;
124   }
125
126   AliError(Form("Invalid stack given: %i", stack));
127   return kFALSE;
128 }
129
130 Bool_t AliTRDgtuTMU::Reset()
131 {
132   // delete all tracks
133
134   for (Int_t zch = 0; zch < fGtuParam->GetNZChannels(); zch++) {
135     for (Int_t reflayeridx = 0; reflayeridx < fGtuParam->GetNRefLayers(); reflayeridx++) {
136       fTracks[zch][reflayeridx].Clear();
137     }
138   }
139
140   // delete all added tracklets
141   for (Int_t iLink = 0; iLink < fGtuParam->GetNLinks(); iLink++) {
142     fTracklets[iLink]->Clear();
143   }
144   for (Int_t layer = 0; layer < fGtuParam->GetNLinks()/2; layer++) {
145     fTrackletsPostInput[layer]->Clear();
146     for (Int_t zch = 0; zch < fGtuParam->GetNZChannels(); zch++)
147       fZChannelTracklets[layer][zch].Clear();
148   }
149
150   fSector = -1;
151   fStack = -1;
152
153   return kTRUE;
154 }
155
156 Bool_t AliTRDgtuTMU::AddTracklet(AliTRDtrackletGTU *tracklet, Int_t link)
157 {
158   // add a tracklet on the given link
159
160   fTracklets[link]->Add(tracklet);
161   return kTRUE;
162 }
163
164
165 Bool_t AliTRDgtuTMU::RunTMU(TList *ListOfTracks, AliESDEvent *esd)
166 {
167   // performs the analysis as in a TMU module of the GTU, i. e.
168   // track matching
169   // calculation of track parameteres (pt, deflection, ???)
170
171   if (fStack < 0 || fSector < 0) {
172     AliError("No valid stack/sector set for this TMU! No tracking!");
173     return kFALSE;
174   }
175
176   // ----- Input units -----
177   AliDebug(1,"--------- Running Input units ----------");
178   for (Int_t layer = 0; layer < fGtuParam->GetNLayers(); layer++) {
179     if (!RunInputUnit(layer)) {
180       AliError(Form("Input unit in layer %i failed", layer));
181       return kFALSE;
182     }
183   }
184
185   // ----- Z-channel units -----
186   AliDebug(1,"--------- Running Z-channel units ----------");
187   for (Int_t layer = 0;  layer <  fGtuParam->GetNLayers(); layer++) {
188     fZChannelTracklets[layer] = new TList[fGtuParam->GetNZChannels()];
189     if (!RunZChannelUnit(layer)) {
190       AliError(Form("Z-Channel unit in layer %i failed", layer));
191       return kFALSE;
192     }
193   }
194
195   // ----- track finding -----
196   AliDebug(1,"--------- Running tracking units ----------");
197   for (Int_t zch = 0; zch < fGtuParam->GetNZChannels(); zch++) {
198     if (!RunTrackFinder(zch, ListOfTracks)) {
199       AliError(Form("Track Finder in z-channel %i failed", zch));
200       return kFALSE;
201     }
202   }
203
204   // ----- Track Merging -----
205   if (!RunTrackMerging(ListOfTracks)) {
206     AliError("Track merging failed");
207     return kFALSE;
208   }
209
210   // ----- track reconstruction -----
211   if (!RunTrackReconstruction(ListOfTracks)) {
212     AliError("Track reconstruction failed");
213     return kFALSE;
214   }
215
216   // ----- label calculation and ESD storage -----
217   TIter next(ListOfTracks);
218   while (AliTRDtrackGTU *trk = (AliTRDtrackGTU*) next()) {
219     trk->CookLabel();
220     if (esd) {
221       AliESDTrdTrack *trdtrack = trk->CreateTrdTrack();
222       esd->AddTrdTrack(trdtrack);
223       delete trdtrack;
224     }
225   }
226
227   return kTRUE;
228 }
229
230 Bool_t AliTRDgtuTMU::RunInputUnit(Int_t layer)
231 {
232   // precalculations for the tracking and reconstruction
233
234   Int_t iTrkl0 = 0; // A-side tracklet
235   Int_t iTrkl1 = 0; // B-side tracklet
236
237   while ((iTrkl0 < fTracklets[2*layer + 0]->GetEntriesFast()) || (iTrkl1 < fTracklets[2*layer + 1]->GetEntriesFast())) {
238     if (iTrkl1 >= fTracklets[2*layer + 1]->GetEntriesFast()) {
239       fTrackletsPostInput[layer]->AddLast(fTracklets[2*layer + 0]->At(iTrkl0));
240       iTrkl0++;
241     }
242     else {
243       if (iTrkl0 >= fTracklets[2*layer + 0]->GetEntriesFast()) {
244         fTrackletsPostInput[layer]->AddLast(fTracklets[2*layer + 1]->At(iTrkl1));
245         iTrkl1++;
246       }
247       else {
248         if (((AliTRDtrackletGTU*) fTracklets[2*layer + 1]->At(iTrkl1))->GetZbin() <
249             ((AliTRDtrackletGTU*) fTracklets[2*layer + 0]->At(iTrkl0))->GetZbin()) {
250           fTrackletsPostInput[layer]->AddLast(fTracklets[2*layer + 1]->At(iTrkl1));
251           iTrkl1++;
252
253         }
254         else {
255           fTrackletsPostInput[layer]->AddLast(fTracklets[2*layer + 0]->At(iTrkl0));
256           iTrkl0++;
257         }
258       }
259     }
260   }
261
262   TIter next(fTrackletsPostInput[layer]);
263
264   while ( AliTRDtrackletGTU *trk = (AliTRDtrackletGTU*) next() ) {
265     trk->SetIndex(fTrackletsPostInput[layer]->IndexOf(trk));
266
267     Int_t alpha = (trk->GetYbin() >> fGtuParam->GetBitExcessY()) * fGtuParam->GetCiAlpha(layer);
268     alpha = ( 2 * trk->GetdY() - (alpha >> fGtuParam->GetBitExcessAlpha()) + 1 ) >> 1;
269     // wrapping projected alpha as in hardware
270     if ((alpha < -64) || (alpha > 63))
271       AliError(Form("alpha out of range: %i", alpha));
272     alpha = alpha & 0x7f;
273     if (alpha & 0x40)
274       trk->SetAlpha(0xffffffc0 | alpha);
275     else
276       trk->SetAlpha(alpha);
277
278     Int_t yproj = trk->GetdY() * (fGtuParam->GetCiYProj(layer)); //??? sign?
279     yproj = ( ( ( (yproj >> fGtuParam->GetBitExcessYProj()) + trk->GetYbin() ) >> 2) + 1) >> 1;
280     trk->SetYProj(yproj);
281
282     trk->SetYPrime(trk->GetYbin() + fGtuParam->GetYt(fStack, layer, trk->GetZbin()));
283
284     AliDebug(10, Form("0x%08x: idx: %3i, z: %2i, y: %5i, dy: %3i, y': %5i, y_proj: %5i, alpha: %3i, pid: %3i",
285                       trk->GetTrackletWord(), trk->GetIndex(), trk->GetZbin(), trk->GetYbin(), trk->GetdY(), trk->GetYPrime(),
286                       trk->GetYProj(), trk->GetAlpha(), trk->GetPID() ));
287   }
288   return kTRUE;
289 }
290
291 Bool_t AliTRDgtuTMU::RunZChannelUnit(Int_t layer)
292 {
293   // run the z-channel unit
294
295   TIter next(fTrackletsPostInput[layer]);
296
297   while (AliTRDtrackletGTU *trk = (AliTRDtrackletGTU*) next()) {
298     for (Int_t zch = 0; zch < fGtuParam->GetNZChannels(); zch++) {
299       if (fGtuParam->IsInZChannel(fStack, layer, zch, trk->GetZbin()) ) {
300         trk->SetSubChannel(zch, fGtuParam->GetZSubchannel(fStack, layer, zch, trk->GetZbin()) );
301
302         TIter nexttrkl(&fZChannelTracklets[layer][zch], kIterBackward);
303         AliTRDtrackletGTU *t = 0x0;
304         while ((t = (AliTRDtrackletGTU*) nexttrkl.Next())) {
305           if (t->GetSubChannel(zch) < trk->GetSubChannel(zch) ||
306               ((t->GetSubChannel(zch) == trk->GetSubChannel(zch)) && (t->GetYProj() < trk->GetYProj())) ) {
307             break;
308           }
309         }
310         if (t)
311           fZChannelTracklets[layer][zch].AddAfter(t, trk);
312         else
313           fZChannelTracklets[layer][zch].AddFirst(trk);
314       }
315     }
316     AliDebug(10, Form("stack %d, layer %2d: 0x%08x Z(2..0)=%i/%i/%i",
317                       fStack, layer, trk->GetTrackletWord(),
318                       fGtuParam->IsInZChannel(fStack, layer, 2, trk->GetZbin()) ? trk->GetSubChannel(2) : -1,
319                       fGtuParam->IsInZChannel(fStack, layer, 1, trk->GetZbin()) ? trk->GetSubChannel(1) : -1,
320                       fGtuParam->IsInZChannel(fStack, layer, 0, trk->GetZbin()) ? trk->GetSubChannel(0) : -1
321                       ));
322   }
323   return kTRUE;
324 }
325
326 Bool_t AliTRDgtuTMU::RunTrackFinder(Int_t zch, TList* /* ListOfTracks */)
327 {
328   // run the track finding
329
330    Int_t        *notr = new Int_t[fGtuParam->GetNLayers()];
331    Int_t        *ptrA = new Int_t[fGtuParam->GetNLayers()];
332    Int_t        *ptrB = new Int_t[fGtuParam->GetNLayers()];
333
334    Bool_t       *bHitA     = new Bool_t[fGtuParam->GetNLayers()];
335    Bool_t       *bHitB     = new Bool_t[fGtuParam->GetNLayers()];
336    Bool_t       *bAligned  = new Bool_t[fGtuParam->GetNLayers()];
337    Bool_t       *bAlignedA = new Bool_t[fGtuParam->GetNLayers()];
338    Bool_t       *bAlignedB = new Bool_t[fGtuParam->GetNLayers()];
339    Bool_t       *bDone     = new Bool_t[fGtuParam->GetNLayers()];
340    Bool_t        ready;
341
342    Int_t        *inc      = new Int_t[fGtuParam->GetNLayers()];
343    Int_t        *incprime = new Int_t[fGtuParam->GetNLayers()];
344
345 // ----- signals within current layer -----
346    Int_t        yPlus;
347    Int_t        yMinus;
348    Int_t        ybPlus;
349    Int_t        ybMinus;
350    Int_t        alphaPlus;
351    Int_t        alphaMinus;
352    Int_t        nHits;
353    Int_t        nUnc;
354    Int_t        nWayBeyond;
355
356    AliTRDtrackletGTU    *trkRA  = 0x0;  // reference tracklet A
357    AliTRDtrackletGTU    *trkRB  = 0x0;  // reference tracklet B
358    AliTRDtrackletGTU    *trkA   = 0x0;  // tracklet A in current layer
359    AliTRDtrackletGTU    *trkB   = 0x0;  // tracklet B in current layer
360 /*
361    AliTRDtrackletGTU    *trkEnd = new AliTRDtrackletGTU();
362    for (Int_t i = 0; i < fGtuParam->GetNZChannels(); i++)
363        trkEnd->SetSubChannel(i, 7);
364    trkEnd->SetYProj(0);
365    trkEnd->SetAlpha(0);
366 */
367
368    for (Int_t refLayerIdx = 0; refLayerIdx < fGtuParam->GetNRefLayers(); refLayerIdx++) {
369      Int_t reflayer = fGtuParam->GetRefLayer(refLayerIdx);
370      AliDebug(5,Form("Tracking for z-channel: %i, reflayer: %i", zch, reflayer));
371
372      ready = kFALSE; // ready if all channels done
373
374 //      for (Int_t iLayer = 0; iLayer < fGtuParam->GetNLayers(); iLayer++) {
375 //        for (Int_t iTrkl = 0; iTrkl < fZChannelTracklets[iLayer][zch].GetEntries(); iTrkl++) {
376 //       printf("%4i/%4i(%i,%i) ",
377 //              ((AliTRDtrackletGTU*) fZChannelTracklets[iLayer][zch].At(iTrkl))->GetYProj(),
378 //              ((AliTRDtrackletGTU*) fZChannelTracklets[iLayer][zch].At(iTrkl))->GetAlpha(),
379 //              ((AliTRDtrackletGTU*) fZChannelTracklets[iLayer][zch].At(iTrkl))->GetIndex(),
380 //              ((AliTRDtrackletGTU*) fZChannelTracklets[iLayer][zch].At(iTrkl))->GetSubChannel(zch)
381 //              );
382 //        }
383 //        printf("\n");
384 //      }
385
386      for (Int_t layer = 0; layer < fGtuParam->GetNLayers(); layer++) {
387        notr[layer] = fZChannelTracklets[layer][zch].GetEntries();
388        ptrA[layer] = 0; // notr[layer] > 0 ? 0 : -1;
389        ptrB[layer] = 1; // notr[layer] > 1 ? 1 : -1;
390
391 // not necessary here
392 //       bDone[layer] = (ptrB[layer] >= notr[layer] - 1); // trkB is last one
393 //       bDone[layer] = (notr[layer] <= 0);
394 //       ready = ready && bDone[layer];
395
396        if (reflayer == 1)
397          AliDebug(5,Form("in layer: %i (zchannel = %i) there are: %i tracklets", layer, zch, notr[layer]));
398      }
399
400      if (ptrA[reflayer] < 0 && ptrB[reflayer] < 0)
401        continue;
402
403      while (!ready) {
404        // ----- reference tracklets -----
405        trkRA = 0x0;
406        trkRB = 0x0;
407        if (0 <= ptrA[reflayer] && ptrA[reflayer] < notr[reflayer])
408          trkRA = (AliTRDtrackletGTU*) fZChannelTracklets[reflayer][zch].At(ptrA[reflayer]);
409        else  {
410          AliDebug(10,Form("No valid tracklet in the reference at ptr: %i! Nothing done!", ptrA[reflayer]));
411          break;
412        }
413
414        if (0 <= ptrB[reflayer] && ptrB[reflayer] < notr[reflayer])
415          trkRB = (AliTRDtrackletGTU*) fZChannelTracklets[reflayer][zch].At(ptrB[reflayer]);
416
417        yPlus      = trkRA->GetYProj() + fGtuParam->GetDeltaY();
418        yMinus     = trkRA->GetYProj() - fGtuParam->GetDeltaY();
419        alphaPlus  = trkRA->GetAlpha() + fGtuParam->GetDeltaAlpha();
420        alphaMinus = trkRA->GetAlpha() - fGtuParam->GetDeltaAlpha();
421        if (trkRB) {
422          ybPlus           = trkRB->GetYProj() + fGtuParam->GetDeltaY();
423          ybMinus          = trkRB->GetYProj() - fGtuParam->GetDeltaY();
424        }
425        else { // irrelevant (should be, is it?)
426          ybPlus           = trkRA->GetYProj() + fGtuParam->GetDeltaY();
427          ybMinus          = trkRA->GetYProj() - fGtuParam->GetDeltaY();
428        }
429
430        nHits      = 0;
431        nUnc       = 0;
432        nWayBeyond = 0;
433
434        for (Int_t layer = 0; layer < fGtuParam->GetNLayers(); layer++) {
435          bHitA[layer] = bHitB[layer] = bAligned[layer] = kFALSE;
436          inc[layer] = incprime[layer] = 0;
437
438          if (layer == reflayer) {
439            bHitA[layer] = kTRUE;
440            bAligned[layer] = kTRUE;
441            nHits++;
442            continue;
443          }
444
445          trkA = 0x0;
446          trkB = 0x0;
447          if (0 <= ptrA[layer] && ptrA[layer] < notr[layer])
448            trkA = (AliTRDtrackletGTU*) fZChannelTracklets[layer][zch].At(ptrA[layer]);
449          if (0 <= ptrB[layer] && ptrB[layer] < notr[layer])
450            trkB = (AliTRDtrackletGTU*) fZChannelTracklets[layer][zch].At(ptrB[layer]);
451
452          bAlignedA[layer] = kFALSE;
453          bAlignedB[layer] = kFALSE;
454
455          if (trkA) {
456            bHitA[layer] = ( !(trkA->GetSubChannel(zch) < trkRA->GetSubChannel(zch) || (trkA->GetSubChannel(zch) == trkRA->GetSubChannel(zch) && trkA->GetYProj() < yMinus) ) &&
457                             !(trkA->GetSubChannel(zch) > trkRA->GetSubChannel(zch) || (trkA->GetSubChannel(zch) == trkRA->GetSubChannel(zch) && trkA->GetYProj() > yPlus)  ) &&
458                             !(trkA->GetAlpha() < alphaMinus) &&
459                             !(trkA->GetAlpha() > alphaPlus) );
460            bAlignedA[layer] = !(trkA->GetSubChannel(zch) < trkRA->GetSubChannel(zch) || (trkA->GetSubChannel(zch) == trkRA->GetSubChannel(zch) && trkA->GetYProj() < yMinus) );
461          }
462          else {
463            bHitA[layer] = 0;
464            bAlignedA[layer] = kTRUE;
465          }
466
467          if (trkB) {
468            bHitB[layer] = ( !(trkB->GetSubChannel(zch) < trkRA->GetSubChannel(zch) || (trkB->GetSubChannel(zch) == trkRA->GetSubChannel(zch) && trkB->GetYProj() < yMinus) ) &&
469                             !(trkB->GetSubChannel(zch) > trkRA->GetSubChannel(zch) || (trkB->GetSubChannel(zch) == trkRA->GetSubChannel(zch) && trkB->GetYProj() > yPlus) ) &&
470                             !(trkB->GetAlpha() < alphaMinus) &&
471                             !(trkB->GetAlpha() > alphaPlus) );
472            bAlignedB[layer] = (trkB->GetSubChannel(zch) > trkRA->GetSubChannel(zch) || (trkB->GetSubChannel(zch) == trkRA->GetSubChannel(zch) && trkB->GetYProj() > yPlus) );
473          }
474          else {
475            bHitB[layer] = 0;
476            bAlignedB[layer] = kTRUE;
477          }
478
479          bAligned[layer] = bAlignedA[layer] || bAlignedB[layer]; //???
480 //       bAligned[layer] = bAlignedA[layer]; //???
481
482          if (bAligned[layer] && (bHitA[layer] || bHitB[layer]) )
483            nHits++;
484          else if (!bAligned[layer] )
485            nUnc++;
486          if (trkRB) {
487            if (trkA) {
488              if ((trkA->GetSubChannel(zch) > trkRB->GetSubChannel(zch)) || (trkA->GetSubChannel(zch) == trkRB->GetSubChannel(zch) && trkA->GetYProj() > ybPlus) )
489                nWayBeyond++;
490            }
491            else
492              nWayBeyond++;
493          }
494
495          //  pre-calculation for the layer shifting (alignment w. r. t. trkRB)
496          if (trkA) {
497              if(trkRB) {
498                  if ((trkA->GetSubChannel(zch) < trkRB->GetSubChannel(zch)) || (trkA->GetSubChannel(zch) == trkRB->GetSubChannel(zch) && trkA->GetYProj() < ybMinus )) // could trkA be aligned for trkRB
499                      incprime[layer] = 1;
500                  else
501                      incprime[layer] = 0;
502              }
503              else
504                  incprime[layer] = 1;
505
506              if (trkB) {
507                  if (trkRB) {
508                      if ((trkB->GetSubChannel(zch) < trkRB->GetSubChannel(zch)) || (trkB->GetSubChannel(zch) == trkRB->GetSubChannel(zch) && trkB->GetYProj() < ybMinus )) // could trkB be aligned for trkRB
509                          incprime[layer] = 2;
510                  }
511                  else
512                      incprime[layer] = 2;
513              }
514          }
515        } // end of loop over layers
516
517        AliDebug(5,Form("logic calculation finished, Nhits: %i %s",
518                        nHits, (nHits >= 4) ? "-> track found" : ""));
519
520        if (nHits >= 4) {
521          // ----- track registration -----
522          AliTRDtrackGTU *track = new AliTRDtrackGTU();
523          track->SetSector(fSector);
524          track->SetStack(fStack);
525          for (Int_t layer = 0; layer < fGtuParam->GetNLayers(); layer++ ) {
526            if (bHitA[layer] || layer == reflayer)
527              track->AddTracklet((AliTRDtrackletGTU* ) fZChannelTracklets[layer][zch].At(ptrA[layer]), layer );
528            else if (bHitB[layer])
529              track->AddTracklet((AliTRDtrackletGTU* ) fZChannelTracklets[layer][zch].At(ptrB[layer]), layer );
530          }
531
532          Bool_t registerTrack = kTRUE;
533          for (Int_t layerIdx = refLayerIdx-1; layerIdx >= 0; layerIdx--) {
534            if (track->IsTrackletInLayer(fGtuParam->GetRefLayer(layerIdx))) {
535              if ((track->GetTracklet(fGtuParam->GetRefLayer(layerIdx)))->GetSubChannel(zch) > 0) {
536                AliDebug(1,"Not registered");
537                  registerTrack = kFALSE;
538              }
539            }
540          }
541          if (registerTrack) {
542              track->SetZChannel(zch);
543              track->SetRefLayerIdx(refLayerIdx);
544              fTracks[zch][refLayerIdx].Add(track);
545          }
546        }
547
548        if ( (nUnc != 0) && (nUnc + nHits >= 4) ) // could this position of the reference layer give some track //??? special check in case of hit?
549           inc[reflayer] = 0;
550        else if (nWayBeyond > 2) // no track possible for both reference tracklets
551          inc[reflayer] = 2;
552        else
553          inc[reflayer] = 1;
554
555        if (inc[reflayer] != 0) // reflayer is shifted
556          for (Int_t layer = 0; layer < fGtuParam->GetNLayers(); layer++) {
557            if (layer == reflayer)
558              continue;
559            inc[layer] = incprime[layer];
560          }
561        else { // reflayer is not shifted
562          for (Int_t layer = 0; layer < fGtuParam->GetNLayers(); layer++) {
563            if (layer == reflayer || notr[layer] == 0)
564              continue;
565            inc[layer] = 0;
566            trkA = 0x0;
567            trkB = 0x0;
568            if (0 <= ptrA[layer] && ptrA[layer] < notr[layer])
569              trkA = (AliTRDtrackletGTU*) fZChannelTracklets[layer][zch].At(ptrA[layer]);
570
571            if (0 <= ptrB[layer] && ptrB[layer] < notr[layer])
572              trkB = (AliTRDtrackletGTU*) fZChannelTracklets[layer][zch].At(ptrB[layer]);
573
574            if (trkA) {
575              if ( !(trkA->GetSubChannel(zch) < trkRA->GetSubChannel(zch) || (trkA->GetSubChannel(zch) == trkRA->GetSubChannel(zch) && trkA->GetYProj() < yMinus) ) &&
576                   !(trkA->GetSubChannel(zch) > trkRA->GetSubChannel(zch) || (trkA->GetSubChannel(zch) == trkRA->GetSubChannel(zch) && trkA->GetYProj() > yPlus ) ) )  // trkA could hit trkRA
577              // if ( !(trkA->GetSubChannel(zch) < trkRA->GetSubChannel(zch) || (trkA->GetSubChannel(zch) == trkRA->GetSubChannel(zch) && trkA->GetYProj() < yMinus) ) )
578                    inc[layer] = 0;
579                else if (trkB) {
580                    if ( trkB->GetSubChannel(zch) < trkRA->GetSubChannel(zch) || (trkB->GetSubChannel(zch) == trkRA->GetSubChannel(zch) && trkB->GetYProj() < yMinus) )
581                        inc[layer] = 2;
582                    else if ( !(trkB->GetSubChannel(zch) > trkRA->GetSubChannel(zch) || (trkB->GetSubChannel(zch) == trkRA->GetSubChannel(zch) && trkB->GetYProj() > yPlus) ) )
583                        inc[layer] = 1;
584                    else
585                        inc[layer] = incprime[layer];
586                }
587                else
588                    inc[layer] = incprime[layer];
589            }
590          }
591        }
592
593        ready = kTRUE;
594        for (Int_t layer = fGtuParam->GetNLayers()-1; layer >= 0; layer--) {
595
596          bDone[layer] = ptrB[layer] < 0 || ptrB[layer] >= notr[layer];
597          ready = ready && bDone[layer];
598
599          if (inc[layer] != 0 && ptrA[layer] >= notr[layer])
600            AliError(Form("Invalid increment: %i at ptrA: %i, notr: %i", inc[layer], ptrA[layer], notr[layer]));
601
602 //       AliInfo(Form("Shifting layer: %i, notr: %i, ptrA: %i, ptrB: %i, inc: %i", layer, notr[layer], ptrA[layer], ptrB[layer], inc[layer]));
603          AliDebug(10,Form("    Layer: %i   %2i(%2i, %2i, %4i, %3i)%s%s   %2i(%2i, %2i, %4i, %3i)%s%s   +%i   %s  (no: %i)",
604                           layer,
605                           ptrA[layer],
606                           (0 <= ptrA[layer] && ptrA[layer] < notr[layer]) ? ((AliTRDtrackletGTU*) fZChannelTracklets[layer][zch].At(ptrA[layer]))->GetIndex() : -1,
607                           (0 <= ptrA[layer] && ptrA[layer] < notr[layer]) ? ((AliTRDtrackletGTU*) fZChannelTracklets[layer][zch].At(ptrA[layer]))->GetSubChannel(zch) : -1,
608                           (0 <= ptrA[layer] && ptrA[layer] < notr[layer]) ? ((AliTRDtrackletGTU*) fZChannelTracklets[layer][zch].At(ptrA[layer]))->GetYProj() : -1,
609                           (0 <= ptrA[layer] && ptrA[layer] < notr[layer]) ? ((AliTRDtrackletGTU*) fZChannelTracklets[layer][zch].At(ptrA[layer]))->GetAlpha() : -1,
610                           bHitA[layer] ? "*" : " ", bAlignedA[layer] ? "+" : " ",
611                           ptrB[layer],
612                           (0 <= ptrB[layer] && ptrB[layer] < notr[layer]) ? ((AliTRDtrackletGTU*) fZChannelTracklets[layer][zch].At(ptrB[layer]))->GetIndex() : -1,
613                           (0 <= ptrB[layer] && ptrB[layer] < notr[layer]) ? ((AliTRDtrackletGTU*) fZChannelTracklets[layer][zch].At(ptrB[layer]))->GetSubChannel(zch) : -1,
614                           (0 <= ptrB[layer] && ptrB[layer] < notr[layer]) ? ((AliTRDtrackletGTU*) fZChannelTracklets[layer][zch].At(ptrB[layer]))->GetYProj() : -1,
615                           (0 <= ptrB[layer] && ptrB[layer] < notr[layer]) ? ((AliTRDtrackletGTU*) fZChannelTracklets[layer][zch].At(ptrB[layer]))->GetAlpha() : -1,
616                           bHitB[layer] ? "*" : " ", bAlignedB[layer] ? "+" : " ",
617                           inc[layer], bDone[layer] ? "done" : "    ", notr[layer]));
618          ptrA[layer] += inc[layer];
619          ptrB[layer] += inc[layer];
620        }
621
622      } // end of while
623    } // end of loop over reflayer
624
625    delete [] bHitA;
626    delete [] bHitB;
627    delete [] bAligned;
628    delete [] bDone;
629    delete [] inc;
630    delete [] incprime;
631    delete [] bAlignedA;
632    delete [] bAlignedB;
633    delete [] notr;
634    delete [] ptrA;
635    delete [] ptrB;
636
637    return kTRUE;
638 }
639
640 Bool_t AliTRDgtuTMU::RunTrackMerging(TList* ListOfTracks)
641 {
642     TList **tracksRefMerged = new TList*[fGtuParam->GetNZChannels()];
643     TList **tracksRefUnique = new TList*[fGtuParam->GetNZChannels()];
644
645     for (Int_t zch = 0; zch < fGtuParam->GetNZChannels(); zch++) {
646         tracksRefMerged[zch] = new TList;
647         tracksRefUnique[zch] = new TList;
648     }
649
650     TList *tracksZMergedStage0 = new TList;
651     TList *tracksZUniqueStage0 = new TList;
652
653     TList **tracksZSplitted = new TList*[2];
654     for (Int_t i = 0; i < 2; i++)
655         tracksZSplitted[i] = new TList;
656
657     TList *tracksZMergedStage1 = new TList;
658
659     AliTRDtrackGTU **trkInRefLayer = new AliTRDtrackGTU*[fGtuParam->GetNRefLayers()];
660
661     // Bool_t done = kFALSE;
662     Int_t minIdx = 0;
663     AliTRDtrackGTU *trkStage0 = 0x0;
664
665     for (Int_t zch = 0; zch < fGtuParam->GetNZChannels(); zch++) {
666         // ----- Merging and Unification in Reflayers (seed_merger) -----
667         do {
668           // done = kTRUE;
669             trkStage0 = 0x0;
670             for (Int_t refLayerIdx = 0; refLayerIdx < fGtuParam->GetNRefLayers(); refLayerIdx++) {
671                 trkInRefLayer[refLayerIdx] = (AliTRDtrackGTU*) fTracks[zch][refLayerIdx].First();
672                 if (trkInRefLayer[refLayerIdx] == 0) {
673                     continue;
674                 }
675                 else if (trkStage0 == 0x0 ) {
676                     trkStage0 = trkInRefLayer[refLayerIdx];
677                     minIdx = refLayerIdx;
678                     // done = kFALSE;
679                 }
680                 else if (trkInRefLayer[refLayerIdx]->GetZSubChannel() < trkStage0->GetZSubChannel() ||
681                          (trkInRefLayer[refLayerIdx]->GetZSubChannel() == trkStage0->GetZSubChannel() && trkInRefLayer[refLayerIdx]->GetYapprox() < trkStage0->GetYapprox()) ) {
682                     minIdx = refLayerIdx;
683                     trkStage0 = trkInRefLayer[refLayerIdx];
684                     // done = kFALSE;
685                 }
686             }
687             if (!trkStage0)
688               break;
689             tracksRefMerged[zch]->Add(trkStage0);
690             fTracks[zch][minIdx].RemoveFirst();
691         } while (trkStage0 != 0);
692
693         Uniquifier(tracksRefMerged[zch], tracksRefUnique[zch]);
694
695         AliDebug(2, Form("zchannel %i", zch));
696         TIter trackRefMerged(tracksRefMerged[zch]);
697         while (AliTRDtrackGTU *trk = (AliTRDtrackGTU*) trackRefMerged())
698           AliDebug(2, Form("track ref layer %i : %i %i %i %i %i %i, y=%i, z_idx=%i",
699                            AliTRDgtuParam::GetRefLayer(trk->GetRefLayerIdx()),
700                            trk->GetTrackletIndex(5),
701                            trk->GetTrackletIndex(4),
702                            trk->GetTrackletIndex(3),
703                            trk->GetTrackletIndex(2),
704                            trk->GetTrackletIndex(1),
705                            trk->GetTrackletIndex(0),
706                            trk->GetYapprox() >> 3,
707                            trk->GetZSubChannel()));
708         AliDebug(2, "uniquified:");
709         TIter trackRefUnique(tracksRefUnique[zch]);
710         while (AliTRDtrackGTU *trk = (AliTRDtrackGTU*) trackRefUnique())
711           AliDebug(2, Form("track ref layer %i : %i %i %i %i %i %i, y=%i, z_idx=%i",
712                            AliTRDgtuParam::GetRefLayer(trk->GetRefLayerIdx()),
713                            trk->GetTrackletIndex(5),
714                            trk->GetTrackletIndex(4),
715                            trk->GetTrackletIndex(3),
716                            trk->GetTrackletIndex(2),
717                            trk->GetTrackletIndex(1),
718                            trk->GetTrackletIndex(0),
719                            trk->GetYapprox() >> 3,
720                            trk->GetZSubChannel()));
721     }
722
723 // ----- Merging in zchannels - 1st stage -----
724
725     do {
726       // done = kTRUE;
727         trkStage0 = 0x0;
728         for (Int_t zch = fGtuParam->GetNZChannels() - 1; zch > -1; zch--) {
729             AliTRDtrackGTU *trk = (AliTRDtrackGTU*) tracksRefUnique[zch]->First();
730             if (trk == 0) {
731                 continue;
732             }
733             else if (trkStage0 == 0x0 ) {
734                 trkStage0 = trk;
735                 minIdx = zch;
736                 // done = kFALSE;
737             }
738             else if ( ((trk->GetZChannel() + 3 * trk->GetZSubChannel()) / 2 - 1) <  ((trkStage0->GetZChannel() + 3 * trkStage0->GetZSubChannel()) / 2 -1 ) ||
739                       (((trk->GetZChannel() + 3 * trk->GetZSubChannel()) / 2 - 1) == ((trkStage0->GetZChannel() + 3 * trkStage0->GetZSubChannel()) / 2 -1 ) && (trk->GetYapprox() < trkStage0->GetYapprox()) ) ) {
740                 minIdx = zch;
741                 trkStage0 = trk;
742                 // done = kFALSE;
743             }
744         }
745
746         if (!trkStage0)
747           break;
748         tracksZMergedStage0->Add(trkStage0);
749         tracksRefUnique[minIdx]->RemoveFirst();
750     } while (trkStage0 != 0);
751
752     Uniquifier(tracksZMergedStage0, tracksZUniqueStage0);
753
754     AliDebug(2, "stage 0:");
755     TIter trackZMergedStage0(tracksZMergedStage0);
756     while (AliTRDtrackGTU *trk = (AliTRDtrackGTU*) trackZMergedStage0())
757       AliDebug(2, Form("track ref layer %i : %i %i %i %i %i %i, y=%i, zch=%i idx=%i",
758                        AliTRDgtuParam::GetRefLayer(trk->GetRefLayerIdx()),
759                        trk->GetTrackletIndex(5),
760                        trk->GetTrackletIndex(4),
761                        trk->GetTrackletIndex(3),
762                        trk->GetTrackletIndex(2),
763                        trk->GetTrackletIndex(1),
764                        trk->GetTrackletIndex(0),
765                        trk->GetYapprox() >> 3,
766                        trk->GetZChannel(),
767                        trk->GetZSubChannel()));
768
769     AliDebug(2, "uniquified:");
770     TIter trackZUniqueStage0(tracksZUniqueStage0);
771     while (AliTRDtrackGTU *trk = (AliTRDtrackGTU*) trackZUniqueStage0())
772       AliDebug(2, Form("track ref layer %i : %i %i %i %i %i %i, y=%i, zch=%i idx=%i",
773                        AliTRDgtuParam::GetRefLayer(trk->GetRefLayerIdx()),
774                        trk->GetTrackletIndex(5),
775                        trk->GetTrackletIndex(4),
776                        trk->GetTrackletIndex(3),
777                        trk->GetTrackletIndex(2),
778                        trk->GetTrackletIndex(1),
779                        trk->GetTrackletIndex(0),
780                        trk->GetYapprox() >> 3,
781                        trk->GetZChannel(),
782                        trk->GetZSubChannel()));
783
784 // ----- Splitting in z -----
785
786     TIter next(tracksZUniqueStage0);
787     while (AliTRDtrackGTU *trk = (AliTRDtrackGTU*) next()) {
788         tracksZSplitted[(trk->GetZChannel() + 3 * (trk->GetZSubChannel() - 1)) % 2]->Add(trk);
789     }
790
791     for (Int_t i = 0; i < 2; i++) {
792       AliDebug(2, Form("split %i:", i));
793       TIter trackZSplit(tracksZSplitted[i]);
794       while (AliTRDtrackGTU *trk = (AliTRDtrackGTU*) trackZSplit())
795         AliDebug(2, Form("track ref layer %i : %i %i %i %i %i %i, y=%i, zch=%i idx=%i",
796                          AliTRDgtuParam::GetRefLayer(trk->GetRefLayerIdx()),
797                          trk->GetTrackletIndex(5),
798                          trk->GetTrackletIndex(4),
799                          trk->GetTrackletIndex(3),
800                          trk->GetTrackletIndex(2),
801                          trk->GetTrackletIndex(1),
802                          trk->GetTrackletIndex(0),
803                          trk->GetYapprox() >> 3,
804                          trk->GetZChannel(),
805                          trk->GetZSubChannel()));
806     }
807
808 // ----- Merging in zchanels - 2nd stage -----
809
810     do {
811       // done = kTRUE;
812         trkStage0 = 0x0;
813         for (Int_t i = 1; i >= 0; i--) {
814             AliTRDtrackGTU *trk = (AliTRDtrackGTU*) tracksZSplitted[i]->First();
815             if (trk == 0) {
816                 continue;
817             }
818             else if (trkStage0 == 0x0 ) {
819                 trkStage0 = trk;
820                 minIdx = i;
821                 // done = kFALSE;
822             }
823             else if ( (((trk->GetZChannel() + 3 * (trk->GetZSubChannel() - 1)) / 2) <  ((trkStage0->GetZChannel() + 3 * (trkStage0->GetZSubChannel() - 1)) / 2)) ||
824                       ((((trk->GetZChannel() + 3 * (trk->GetZSubChannel() - 1)) / 2) == ((trkStage0->GetZChannel() + 3 * (trkStage0->GetZSubChannel() - 1)) / 2)) && (trk->GetYapprox() < trkStage0->GetYapprox()) ) ) {
825                 minIdx = i;
826                 trkStage0 = trk;
827                 // done = kFALSE;
828             }
829         }
830
831         if (!trkStage0)
832           break;
833         tracksZMergedStage1->Add(trkStage0);
834         tracksZSplitted[minIdx]->RemoveFirst();
835     } while (trkStage0 != 0);
836
837     Uniquifier(tracksZMergedStage1, ListOfTracks);
838
839     AliDebug(2, "merged:");
840     TIter trackZMergedStage1(tracksZMergedStage1);
841     while (AliTRDtrackGTU *trk = (AliTRDtrackGTU*) trackZMergedStage1())
842       AliDebug(2, Form("track ref layer %i : %i %i %i %i %i %i, y=%i, zch=%i idx=%i",
843                        AliTRDgtuParam::GetRefLayer(trk->GetRefLayerIdx()),
844                        trk->GetTrackletIndex(5),
845                        trk->GetTrackletIndex(4),
846                        trk->GetTrackletIndex(3),
847                        trk->GetTrackletIndex(2),
848                        trk->GetTrackletIndex(1),
849                        trk->GetTrackletIndex(0),
850                        trk->GetYapprox() >> 3,
851                        trk->GetZChannel(),
852                        trk->GetZSubChannel()));
853
854     AliDebug(2, "uniquified:");
855     TIter track(ListOfTracks);
856     while (AliTRDtrackGTU *trk = (AliTRDtrackGTU*) track())
857       AliDebug(2, Form("track ref layer %i : %i %i %i %i %i %i, y=%i, zch=%i idx=%i",
858                        AliTRDgtuParam::GetRefLayer(trk->GetRefLayerIdx()),
859                        trk->GetTrackletIndex(5),
860                        trk->GetTrackletIndex(4),
861                        trk->GetTrackletIndex(3),
862                        trk->GetTrackletIndex(2),
863                        trk->GetTrackletIndex(1),
864                        trk->GetTrackletIndex(0),
865                        trk->GetYapprox() >> 3,
866                        trk->GetZChannel(),
867                        trk->GetZSubChannel()));
868
869     // cleaning up
870     for (Int_t zch = 0; zch < fGtuParam->GetNZChannels(); zch++) {
871       delete tracksRefMerged[zch];
872       delete tracksRefUnique[zch];
873     }
874     delete [] tracksRefMerged;
875     delete [] tracksRefUnique;
876
877
878     delete tracksZMergedStage0;
879     delete tracksZUniqueStage0;
880
881     for (Int_t i = 0; i < 2; i++)
882       delete tracksZSplitted[i];
883     delete [] tracksZSplitted;
884
885     delete tracksZMergedStage1;
886
887     delete [] trkInRefLayer;
888
889     return kTRUE;
890 }
891
892 Bool_t AliTRDgtuTMU::RunTrackReconstruction(TList* ListOfTracks)
893 {
894   // run the track reconstruction for all tracks in the list
895
896   TIter next(ListOfTracks);
897   while (AliTRDtrackGTU *track = (AliTRDtrackGTU*) next()) {
898     CalculateTrackParams(track);
899     CalculatePID(track);
900     AliDebug(1, Form("track with pid = %i", track->GetPID()));
901   }
902   return kTRUE;
903 }
904
905 Bool_t AliTRDgtuTMU::CalculatePID(AliTRDtrackGTU *track)
906 {
907   // calculate PID for the given track
908   if (!track) {
909     AliError("No track to calculate!");
910     return kFALSE;
911   }
912
913   if (AliTRDgtuParam::GetUseGTUconst()) {
914     // averaging as in GTU
915     ULong64_t coeff;
916
917     // selection of coefficient for averaging
918     switch(track->GetTrackletMask()){
919     case 0x3f:
920       // 6 tracklets
921       coeff=0x5558; // approx. 1/6
922       break;
923
924     case 0x3e:
925     case 0x3d:
926     case 0x3b:
927     case 0x37:
928     case 0x2f:
929     case 0x1f:
930       // 5 tracklets
931       coeff=0x6666; // approx. 1/5
932       break;
933
934     default:
935       // 4 tracklets
936       coeff=0x8000; // = 0.25
937     }
938     coeff &= 0x1ffff; // 17-bit constant
939
940     ULong64_t sum = 0;
941     for (Int_t iLayer = 0; iLayer < fGtuParam->GetNLayers(); iLayer++) {
942       if ((track->GetTrackletMask() >> iLayer) & 1) {
943         sum += track->GetTracklet(iLayer)->GetPID();
944       }
945     }
946
947     sum = sum & 0x7ff;
948     ULong64_t prod   = (sum * coeff) & 0xfffffffffull;
949     ULong64_t prodFinal = ((prod >> 17) + ((prod >> 16) & 1)) & 0xff;
950
951     track->SetPID(prodFinal & 0xff);
952
953     return kTRUE;
954   }
955   else {
956
957     // simple averaging
958     Int_t nTracklets = 0;
959     Int_t pidSum = 0;
960     for (Int_t layer = 0; layer < fGtuParam->GetNLayers(); layer++) {
961       if (!track->IsTrackletInLayer(layer)) {
962         continue;
963       }
964       AliTRDtrackletGTU *trk = track->GetTracklet(layer);
965       pidSum += trk->GetPID();
966       nTracklets++;
967     }
968
969     if (nTracklets>0)
970       track->SetPID(pidSum/nTracklets);
971     else
972       AliError("Track without contributing tracklets, no PID assigned");
973
974     return kTRUE;
975   }
976 }
977
978 Bool_t AliTRDgtuTMU::CalculateTrackParams(AliTRDtrackGTU *track)
979 {
980   // calculate the track parameters
981
982   if (!track) {
983     AliError("No track to calculate!");
984     return kFALSE;
985   }
986
987   Int_t a = 0;
988   Float_t b = 0;
989   Float_t c = 0;
990   Float_t x1;
991   Float_t x2;
992
993   AliDebug(5,Form("There are %i tracklets in this track.", track->GetNTracklets()));
994
995   for (Int_t layer = 0; layer < fGtuParam->GetNLayers(); layer++) {
996     if (!track->IsTrackletInLayer(layer)) {
997       continue;
998     }
999     AliTRDtrackletGTU *trk = track->GetTracklet(layer);
1000     if (!trk) {
1001       AliError(Form("Could not get tracklet in layer %i\n", layer));
1002       continue;
1003     }
1004     AliDebug(10,Form("  layer %i trk yprime: %6i, aki: %6i", layer, trk->GetYPrime(), (Int_t) (2048 * fGtuParam->GetAki(track->GetTrackletMask(), layer))));
1005     a += (((Int_t) (2048 * fGtuParam->GetAki(track->GetTrackletMask(), layer))) * trk->GetYPrime() + 1) >> 8;
1006     b += fGtuParam->GetBki(track->GetTrackletMask(), layer) * trk->GetYPrime() * fGtuParam->GetBinWidthY();
1007     c += fGtuParam->GetCki(track->GetTrackletMask(), layer) * trk->GetYPrime() * fGtuParam->GetBinWidthY();
1008   }
1009   if (a < 0)
1010       a += 3;
1011   a = a >> 2;
1012
1013   track->SetFitParams(a, b, c);
1014
1015   fGtuParam->GetIntersectionPoints(track->GetTrackletMask(), x1, x2);
1016
1017   AliDebug(5,Form("  Track parameters: a = %i, b = %f, c = %f, x1 = %f, x2 = %f, pt = %f (trkl mask: %i)", a, b, c, x1, x2, track->GetPt(), track->GetTrackletMask()));
1018
1019   return kTRUE;
1020 }
1021
1022
1023 Bool_t AliTRDgtuTMU::Uniquifier(const TList *inlist, TList *outlist)
1024 {
1025   // remove multiple occurences of the same track
1026
1027     TIter next(inlist);
1028     AliTRDtrackGTU *trkStage0 = 0x0;
1029     AliTRDtrackGTU *trkStage1 = 0x0;
1030
1031     do {
1032         trkStage0 = (AliTRDtrackGTU*) next();
1033
1034         Bool_t tracksEqual = kFALSE;
1035         if (trkStage0 != 0 && trkStage1 != 0) {
1036             for (Int_t layer = 0; layer < fGtuParam->GetNLayers(); layer++) {
1037                 if (trkStage0->GetTrackletIndex(layer) != -1 && trkStage0->GetTrackletIndex(layer) == trkStage1->GetTrackletIndex(layer)) {
1038                     tracksEqual = kTRUE;
1039                     break;
1040                 }
1041             }
1042         }
1043
1044         if (tracksEqual) {
1045             if (trkStage0->GetNTracklets() >= trkStage1->GetNTracklets())
1046                 trkStage1 = trkStage0;
1047         }
1048         else {
1049             if (trkStage1 != 0x0)
1050                 outlist->Add(trkStage1);
1051             trkStage1 = trkStage0;
1052         }
1053
1054     } while (trkStage1 != 0x0);
1055     return kTRUE;
1056 }