]> git.uio.no Git - u/mrichter/AliRoot.git/blame_incremental - TRD/AliTRDgtuTMU.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame_incremental
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (non-incremental) | history | HEAD
Removing centrality check from default:
[u/mrichter/AliRoot.git] / TRD / AliTRDgtuTMU.cxx
... / ...
CommitLineData
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),
47 fTrackletsPostInput(0x0),
48 fZChannelTracklets(0x0),
49 fTrackArray(new TClonesArray("AliTRDtrackGTU", 50)),
50 fTracks(0x0),
51 fGtuParam(0x0),
52 fStack(-1),
53 fSector(-1)
54{
55 // constructor which initializes the position information of the TMU
56
57 fGtuParam = AliTRDgtuParam::Instance();
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()];
67 fZChannelTracklets = new TList*[fGtuParam->GetNLayers()];
68 for (Int_t layer = 0; layer < fGtuParam->GetNLayers(); layer++) {
69 fZChannelTracklets[layer] = new TList[fGtuParam->GetNZChannels()];
70 fTrackletsPostInput[layer] = new TObjArray();
71 }
72
73 fTracks = new TList*[fGtuParam->GetNZChannels()];
74 for (Int_t zch = 0; zch < fGtuParam->GetNZChannels(); zch++) {
75 fTracks[zch] = new TList[fGtuParam->GetNRefLayers()];
76 }
77
78 if (stack > -1)
79 SetStack(stack);
80 if (sector > -1)
81 SetSector(sector);
82}
83
84AliTRDgtuTMU::~AliTRDgtuTMU()
85{
86 // destructor
87
88 for (Int_t zch = 0; zch < fGtuParam->GetNZChannels(); zch++) {
89 delete [] fTracks[zch];
90 }
91 delete [] fTracks;
92 delete fTrackArray;
93
94 for (Int_t layer = 0; layer < fGtuParam->GetNLayers(); layer++) {
95 delete [] fZChannelTracklets[layer];
96 delete fTrackletsPostInput[layer];
97 }
98 delete [] fZChannelTracklets;
99 delete [] fTrackletsPostInput;
100
101 for (Int_t iLink = 0; iLink < fGtuParam->GetNLinks(); iLink++) {
102 delete fTracklets[iLink];
103 }
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() ) {
112 fSector = sector;
113 return kTRUE;
114 }
115
116 AliError(Form("Invalid sector given: %i", sector));
117 return kFALSE;
118}
119
120Bool_t AliTRDgtuTMU::SetStack(Int_t stack)
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
133Bool_t AliTRDgtuTMU::Reset()
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
143 fTrackArray->Delete();
144
145 // delete all added tracklets
146 for (Int_t iLink = 0; iLink < fGtuParam->GetNLinks(); iLink++) {
147 fTracklets[iLink]->Clear();
148 }
149 for (Int_t layer = 0; layer < fGtuParam->GetNLinks()/2; layer++) {
150 fTrackletsPostInput[layer]->Clear();
151 for (Int_t zch = 0; zch < fGtuParam->GetNZChannels(); zch++)
152 fZChannelTracklets[layer][zch].Clear();
153 }
154
155 fSector = -1;
156 fStack = -1;
157
158 return kTRUE;
159}
160
161Bool_t AliTRDgtuTMU::AddTracklet(AliTRDtrackletGTU *tracklet, Int_t link)
162{
163 // add a tracklet on the given link
164
165 fTracklets[link]->Add(tracklet);
166 return kTRUE;
167}
168
169
170Bool_t AliTRDgtuTMU::RunTMU(TList *ListOfTracks, AliESDEvent *esd, Int_t outLabel)
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) {
177 AliError("No valid stack/sector set for this TMU! No tracking!");
178 return kFALSE;
179 }
180
181 // ----- Input units -----
182 AliDebug(1,"--------- Running Input units ----------");
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 -----
191 AliDebug(1,"--------- Running Z-channel units ----------");
192 for (Int_t layer = 0; layer < fGtuParam->GetNLayers(); layer++) {
193 if (!RunZChannelUnit(layer)) {
194 AliError(Form("Z-Channel unit in layer %i failed", layer));
195 return kFALSE;
196 }
197 }
198
199 // ----- track finding -----
200 AliDebug(1,"--------- Running tracking units ----------");
201 for (Int_t zch = 0; zch < fGtuParam->GetNZChannels(); zch++) {
202 if (!RunTrackFinder(zch, ListOfTracks)) {
203 AliError(Form("Track Finder in z-channel %i failed", zch));
204 return kFALSE;
205 }
206 }
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
220 // ----- label calculation and ESD storage -----
221 TIter next(ListOfTracks);
222 while (AliTRDtrackGTU *trk = (AliTRDtrackGTU*) next()) {
223 if (outLabel == -1)
224 trk->CookLabel();
225 else
226 trk->SetLabel(outLabel);
227 if (esd) {
228 AliESDTrdTrack *trdtrack = trk->CreateTrdTrack();
229 esd->AddTrdTrack(trdtrack);
230 delete trdtrack;
231 }
232 }
233
234 return kTRUE;
235}
236
237Bool_t AliTRDgtuTMU::RunInputUnit(Int_t layer)
238{
239 // precalculations for the tracking and reconstruction
240
241 Int_t iTrkl0 = 0; // A-side tracklet
242 Int_t iTrkl1 = 0; // B-side tracklet
243
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()))) {
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 }
271 ++nTracklets;
272 }
273
274 TIter next(fTrackletsPostInput[layer]);
275
276 while ( AliTRDtrackletGTU *trk = (AliTRDtrackletGTU*) next() ) {
277 trk->SetIndex(fTrackletsPostInput[layer]->IndexOf(trk));
278
279 Int_t alpha = (trk->GetYbin() >> fGtuParam->GetBitExcessY()) * fGtuParam->GetCiAlpha(layer);
280 alpha = ( 2 * trk->GetdY() - (alpha >> fGtuParam->GetBitExcessAlpha()) + 1 ) >> 1;
281 // wrapping projected alpha as in hardware
282 if ((alpha < -64) || (alpha > 63))
283 AliDebug(1, Form("alpha out of range: %i", alpha));
284 alpha = alpha & 0x7f;
285 if (alpha & 0x40)
286 trk->SetAlpha(0xffffffc0 | alpha);
287 else
288 trk->SetAlpha(alpha);
289
290 Int_t yproj = trk->GetdY() * (fGtuParam->GetCiYProj(layer));
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
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",
297 trk->GetTrackletWord(), trk->GetIndex(), trk->GetZbin(), trk->GetYbin(), trk->GetdY(), trk->GetYPrime(),
298 trk->GetYProj(), trk->GetAlpha(), trk->GetPID(),
299 fGtuParam->GetCiYProj(layer), fGtuParam->GetYt(fStack, layer, trk->GetZbin()) ));
300 }
301 return kTRUE;
302}
303
304Bool_t AliTRDgtuTMU::RunZChannelUnit(Int_t layer)
305{
306 // run the z-channel unit
307
308 TIter next(fTrackletsPostInput[layer]);
309
310 while (AliTRDtrackletGTU *trk = (AliTRDtrackletGTU*) next()) {
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()) );
314
315 TIter nexttrkl(&fZChannelTracklets[layer][zch], kIterBackward);
316 AliTRDtrackletGTU *t = 0x0;
317 while ((t = (AliTRDtrackletGTU*) nexttrkl.Next())) {
318 if (t->GetSubChannel(zch) < trk->GetSubChannel(zch) ||
319 ((t->GetSubChannel(zch) == trk->GetSubChannel(zch)) && (t->GetYProj() < trk->GetYProj())) ) {
320 break;
321 }
322 }
323 if (t)
324 fZChannelTracklets[layer][zch].AddAfter(t, trk);
325 else
326 fZChannelTracklets[layer][zch].AddFirst(trk);
327 }
328 }
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 ));
335 }
336 return kTRUE;
337}
338
339Bool_t AliTRDgtuTMU::RunTrackFinder(Int_t zch, TList* /* ListOfTracks */)
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
347 Bool_t *bHitA = new Bool_t[fGtuParam->GetNLayers()];
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 -----
359 Int_t yPlus;
360 Int_t yMinus;
361 Int_t ybPlus;
362 Int_t ybMinus;
363 Int_t alphaPlus;
364 Int_t alphaMinus;
365 Int_t nHits;
366 Int_t nUnc;
367 Int_t nWayBeyond;
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();
375 for (Int_t i = 0; i < fGtuParam->GetNZChannels(); i++)
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);
383 AliDebug(5,Form("Tracking for z-channel: %i, reflayer: %i", zch, reflayer));
384
385 ready = kFALSE; // ready if all channels done
386
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
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)
410 AliDebug(5,Form("in layer: %i (zchannel = %i) there are: %i tracklets", layer, zch, notr[layer]));
411 }
412
413 if (ptrA[reflayer] < 0 && ptrB[reflayer] < 0)
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 {
423 AliDebug(10,Form("No valid tracklet in the reference at ptr: %i! Nothing done!", ptrA[reflayer]));
424 break;
425 }
426
427 if (0 <= ptrB[reflayer] && ptrB[reflayer] < notr[reflayer])
428 trkRB = (AliTRDtrackletGTU*) fZChannelTracklets[reflayer][zch].At(ptrB[reflayer]);
429
430 yPlus = trkRA->GetYProj() + fGtuParam->GetDeltaY();
431 yMinus = trkRA->GetYProj() - fGtuParam->GetDeltaY();
432 alphaPlus = trkRA->GetAlpha() + fGtuParam->GetDeltaAlpha();
433 alphaMinus = trkRA->GetAlpha() - fGtuParam->GetDeltaAlpha();
434 if (trkRB) {
435 ybPlus = trkRB->GetYProj() + fGtuParam->GetDeltaY();
436 ybMinus = trkRB->GetYProj() - fGtuParam->GetDeltaY();
437 }
438 else { // irrelevant (should be, is it?)
439 ybPlus = trkRA->GetYProj() + fGtuParam->GetDeltaY();
440 ybMinus = trkRA->GetYProj() - fGtuParam->GetDeltaY();
441 }
442
443 nHits = 0;
444 nUnc = 0;
445 nWayBeyond = 0;
446
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;
454 nHits++;
455 continue;
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
468 if (trkA) {
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) );
473 bAlignedA[layer] = !(trkA->GetSubChannel(zch) < trkRA->GetSubChannel(zch) || (trkA->GetSubChannel(zch) == trkRA->GetSubChannel(zch) && trkA->GetYProj() < yMinus) );
474 }
475 else {
476 bHitA[layer] = 0;
477 bAlignedA[layer] = kTRUE;
478 }
479
480 if (trkB) {
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) ) &&
483 !(trkB->GetAlpha() < alphaMinus) &&
484 !(trkB->GetAlpha() > alphaPlus) );
485 bAlignedB[layer] = (trkB->GetSubChannel(zch) > trkRA->GetSubChannel(zch) || (trkB->GetSubChannel(zch) == trkRA->GetSubChannel(zch) && trkB->GetYProj() > yPlus) );
486 }
487 else {
488 bHitB[layer] = 0;
489 bAlignedB[layer] = kTRUE;
490 }
491
492 bAligned[layer] = bAlignedA[layer] || bAlignedB[layer]; //???
493// bAligned[layer] = bAlignedA[layer]; //???
494
495 if (bAligned[layer] && (bHitA[layer] || bHitB[layer]) )
496 nHits++;
497 else if (!bAligned[layer] )
498 nUnc++;
499 if (trkRB) {
500 if (trkA) {
501 if ((trkA->GetSubChannel(zch) > trkRB->GetSubChannel(zch)) || (trkA->GetSubChannel(zch) == trkRB->GetSubChannel(zch) && trkA->GetYProj() > ybPlus) )
502 nWayBeyond++;
503 }
504 else
505 nWayBeyond++;
506 }
507
508 // pre-calculation for the layer shifting (alignment w. r. t. trkRB)
509 if (trkA) {
510 if(trkRB) {
511 if ((trkA->GetSubChannel(zch) < trkRB->GetSubChannel(zch)) || (trkA->GetSubChannel(zch) == trkRB->GetSubChannel(zch) && trkA->GetYProj() < ybMinus )) // could trkA be aligned for trkRB
512 incprime[layer] = 1;
513 else
514 incprime[layer] = 0;
515 }
516 else
517 incprime[layer] = 1;
518
519 if (trkB) {
520 if (trkRB) {
521 if ((trkB->GetSubChannel(zch) < trkRB->GetSubChannel(zch)) || (trkB->GetSubChannel(zch) == trkRB->GetSubChannel(zch) && trkB->GetYProj() < ybMinus )) // could trkB be aligned for trkRB
522 incprime[layer] = 2;
523 }
524 else
525 incprime[layer] = 2;
526 }
527 }
528 } // end of loop over layers
529
530 AliDebug(5,Form("logic calculation finished, Nhits: %i %s",
531 nHits, (nHits >= 4) ? "-> track found" : ""));
532
533 if (nHits >= 4) {
534 // ----- track registration -----
535 AliTRDtrackGTU *track = new ((*fTrackArray)[fTrackArray->GetEntriesFast()]) AliTRDtrackGTU();
536 track->SetSector(fSector);
537 track->SetStack(fStack);
538 for (Int_t layer = 0; layer < fGtuParam->GetNLayers(); layer++ ) {
539 if (bHitA[layer] || layer == reflayer)
540 track->AddTracklet((AliTRDtrackletGTU* ) fZChannelTracklets[layer][zch].At(ptrA[layer]), layer );
541 else if (bHitB[layer])
542 track->AddTracklet((AliTRDtrackletGTU* ) fZChannelTracklets[layer][zch].At(ptrB[layer]), layer );
543 }
544
545 Bool_t registerTrack = kTRUE;
546 for (Int_t layerIdx = refLayerIdx-1; layerIdx >= 0; layerIdx--) {
547 if (track->IsTrackletInLayer(fGtuParam->GetRefLayer(layerIdx))) {
548 if ((track->GetTracklet(fGtuParam->GetRefLayer(layerIdx)))->GetSubChannel(zch) > 0) {
549 AliDebug(1,"Not registered");
550 registerTrack = kFALSE;
551 }
552 }
553 }
554 if (registerTrack) {
555 track->SetZChannel(zch);
556 track->SetRefLayerIdx(refLayerIdx);
557 fTracks[zch][refLayerIdx].Add(track);
558 }
559 }
560
561 if ( (nUnc != 0) && (nUnc + nHits >= 4) ) // could this position of the reference layer give some track //??? special check in case of hit?
562 inc[reflayer] = 0;
563 else if (nWayBeyond > 2) // no track possible for both reference tracklets
564 inc[reflayer] = 2;
565 else
566 inc[reflayer] = 1;
567
568 if (inc[reflayer] != 0) // reflayer is shifted
569 for (Int_t layer = 0; layer < fGtuParam->GetNLayers(); layer++) {
570 if (layer == reflayer)
571 continue;
572 inc[layer] = incprime[layer];
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) {
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) ) )
591 inc[layer] = 0;
592 else if (trkB) {
593 if ( trkB->GetSubChannel(zch) < trkRA->GetSubChannel(zch) || (trkB->GetSubChannel(zch) == trkRA->GetSubChannel(zch) && trkB->GetYProj() < yMinus) )
594 inc[layer] = 2;
595 else if ( !(trkB->GetSubChannel(zch) > trkRA->GetSubChannel(zch) || (trkB->GetSubChannel(zch) == trkRA->GetSubChannel(zch) && trkB->GetYProj() > yPlus) ) )
596 inc[layer] = 1;
597 else
598 inc[layer] = incprime[layer];
599 }
600 else
601 inc[layer] = incprime[layer];
602 }
603 }
604 }
605
606 ready = kTRUE;
607 for (Int_t layer = fGtuParam->GetNLayers()-1; layer >= 0; layer--) {
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]));
616 AliDebug(10,Form(" Layer: %i %2i(%2i, %2i, %4i, %3i)%s%s %2i(%2i, %2i, %4i, %3i)%s%s +%i %s (no: %i)",
617 layer,
618 ptrA[layer],
619 (0 <= ptrA[layer] && ptrA[layer] < notr[layer]) ? ((AliTRDtrackletGTU*) fZChannelTracklets[layer][zch].At(ptrA[layer]))->GetIndex() : -1,
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,
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,
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,
629 bHitB[layer] ? "*" : " ", bAlignedB[layer] ? "+" : " ",
630 inc[layer], bDone[layer] ? "done" : " ", notr[layer]));
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
653Bool_t AliTRDgtuTMU::RunTrackMerging(TList* ListOfTracks)
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
672 AliTRDtrackGTU **trkInRefLayer = new AliTRDtrackGTU*[fGtuParam->GetNRefLayers()];
673
674 // Bool_t done = kFALSE;
675 Int_t minIdx = 0;
676 AliTRDtrackGTU *trkStage0 = 0x0;
677
678 for (Int_t zch = 0; zch < fGtuParam->GetNZChannels(); zch++) {
679 // ----- Merging and Unification in Reflayers (seed_merger) -----
680 do {
681 // done = kTRUE;
682 trkStage0 = 0x0;
683 for (Int_t refLayerIdx = 0; refLayerIdx < fGtuParam->GetNRefLayers(); refLayerIdx++) {
684 trkInRefLayer[refLayerIdx] = (AliTRDtrackGTU*) fTracks[zch][refLayerIdx].First();
685 if (trkInRefLayer[refLayerIdx] == 0) {
686 continue;
687 }
688 else if (trkStage0 == 0x0 ) {
689 trkStage0 = trkInRefLayer[refLayerIdx];
690 minIdx = refLayerIdx;
691 // done = kFALSE;
692 }
693 else if ( (trkInRefLayer[refLayerIdx]->GetZSubChannel() < trkStage0->GetZSubChannel()) ||
694 ((trkInRefLayer[refLayerIdx]->GetZSubChannel() == trkStage0->GetZSubChannel()) &&
695 ((trkInRefLayer[refLayerIdx]->GetYapprox() >> 3) < (trkStage0->GetYapprox() >> 3)) ) ) {
696 minIdx = refLayerIdx;
697 trkStage0 = trkInRefLayer[refLayerIdx];
698 // done = kFALSE;
699 }
700 }
701 if (!trkStage0)
702 break;
703 tracksRefMerged[zch]->Add(trkStage0);
704 fTracks[zch][minIdx].RemoveFirst();
705 } while (trkStage0 != 0);
706
707 Uniquifier(tracksRefMerged[zch], tracksRefUnique[zch]);
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",
713 AliTRDgtuParam::GetRefLayer(trk->GetRefLayerIdx()),
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",
726 AliTRDgtuParam::GetRefLayer(trk->GetRefLayerIdx()),
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()));
735 }
736
737// ----- Merging in zchannels - 1st stage -----
738
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;
760 }
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;
860 }
861 }
862 }
863
864 if (!trkStage0)
865 break;
866 tracksZMergedStage0->Add(trkStage0);
867 tracksRefUnique[minIdx]->RemoveFirst();
868 } while (trkStage0 != 0);
869 }
870
871 Uniquifier(tracksZMergedStage0, tracksZUniqueStage0);
872
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
903// ----- Splitting in z -----
904
905 TIter next(tracksZUniqueStage0);
906 while (AliTRDtrackGTU *trk = (AliTRDtrackGTU*) next()) {
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);
926 }
927
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
945// ----- Merging in zchanels - 2nd stage -----
946
947 do {
948 // done = kTRUE;
949 trkStage0 = 0x0;
950 for (Int_t i = 1; i >= 0; i--) {
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;
958 // done = kFALSE;
959 }
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)) ) ) {
963 minIdx = i;
964 trkStage0 = trk;
965 // done = kFALSE;
966 }
967 }
968
969 if (!trkStage0)
970 break;
971 tracksZMergedStage1->Add(trkStage0);
972 tracksZSplitted[minIdx]->RemoveFirst();
973 } while (trkStage0 != 0);
974
975 Uniquifier(tracksZMergedStage1, ListOfTracks);
976
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
1007 // cleaning up
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];
1021 delete [] tracksZSplitted;
1022
1023 delete tracksZMergedStage1;
1024
1025 delete [] trkInRefLayer;
1026
1027 return kTRUE;
1028}
1029
1030Bool_t AliTRDgtuTMU::RunTrackReconstruction(TList* ListOfTracks)
1031{
1032 // run the track reconstruction for all tracks in the list
1033
1034 TIter next(ListOfTracks);
1035 while (AliTRDtrackGTU *track = (AliTRDtrackGTU*) next()) {
1036 CalculateTrackParams(track);
1037 CalculatePID(track);
1038 AliDebug(1, Form("track with pid = %i", track->GetPID()));
1039 }
1040 return kTRUE;
1041}
1042
1043Bool_t AliTRDgtuTMU::CalculatePID(AliTRDtrackGTU *track)
1044{
1045 // calculate PID for the given track
1046 if (!track) {
1047 AliError("No track to calculate!");
1048 return kFALSE;
1049 }
1050
1051 if (AliTRDgtuParam::GetUseGTUconst()) {
1052 // averaging as in GTU
1053 AliDebug(1, "using GTU constants for PID calculation");
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
1076 }
1077 coeff &= 0x1ffff; // 17-bit constant
1078
1079 ULong64_t sum = 0;
1080 Int_t i = 0;
1081 for (Int_t iLayer = 0; iLayer < fGtuParam->GetNLayers(); iLayer++) {
1082 if ((track->GetTrackletMask() >> iLayer) & 1) {
1083 sum += track->GetTracklet(iLayer)->GetPID();
1084 ++i;
1085 }
1086 }
1087
1088 Float_t av = 1./i * sum;
1089 sum = sum & 0x7ff;
1090 ULong64_t prod = (sum * coeff) & 0xfffffffffull;
1091 ULong64_t prodFinal = ((prod >> 17) + ((prod >> 16) & 1)) & 0xff;
1092
1093 if (TMath::Abs((prodFinal & 0xff) - av) > 0.5)
1094 AliError(Form("problem with PID averaging (hw <-> ar): %3lli <-> %4.1f", prodFinal & 0xff, av));
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 }
1112
1113 if (nTracklets>0)
1114 track->SetPID(pidSum/nTracklets);
1115 else
1116 AliError("Track without contributing tracklets, no PID assigned");
1117
1118 return kTRUE;
1119 }
1120}
1121
1122Bool_t AliTRDgtuTMU::CalculateTrackParams(AliTRDtrackGTU *track)
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
1137 AliDebug(5,Form("There are %i tracklets in this track.", track->GetNTracklets()));
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 }
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;
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
1158 track->SetFitParams(a << 2, b, c);
1159
1160 fGtuParam->GetIntersectionPoints(track->GetTrackletMask(), x1, x2);
1161
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()));
1164
1165 return kTRUE;
1166}
1167
1168
1169Bool_t AliTRDgtuTMU::Uniquifier(const TList *inlist, TList *outlist)
1170{
1171 // remove multiple occurences of the same track
1172
1173 TIter next(inlist);
1174 AliTRDtrackGTU *trkStage0 = 0x0;
1175 AliTRDtrackGTU *trkStage1 = 0x0;
1176
1177 do {
1178 trkStage0 = (AliTRDtrackGTU*) next();
1179
1180 Bool_t tracksEqual = kFALSE;
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)) {
1184 tracksEqual = kTRUE;
1185 break;
1186 }
1187 }
1188 }
1189
1190 if (tracksEqual) {
1191 if (trkStage0->GetNTracklets() >= trkStage1->GetNTracklets())
1192 trkStage1 = trkStage0;
1193 }
1194 else {
1195 if (trkStage1 != 0x0)
1196 outlist->Add(trkStage1);
1197 trkStage1 = trkStage0;
1198 }
1199
1200 } while (trkStage1 != 0x0);
1201 return kTRUE;
1202}