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