Making online tracklets usable in offline reconstruction
[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),
47 fZChannelTracklets(0x0),
48 fTracks(0x0),
49 fGtuParam(0x0),
50 fStack(-1),
51 fSector(-1)
52{
53 fGtuParam = AliTRDgtuParam::Instance();
54 fTracklets = new TObjArray*[fGtuParam->GetNLayers()];
55 fZChannelTracklets = new TList*[fGtuParam->GetNLayers()];
56 for (Int_t layer = 0; layer < fGtuParam->GetNLayers(); layer++) {
57 fTracklets[layer] = new TObjArray();
58 fZChannelTracklets[layer] = new TList[fGtuParam->GetNZChannels()];
59 }
60 fTracks = new TList*[fGtuParam->GetNZChannels()];
61 for (Int_t zch = 0; zch < fGtuParam->GetNZChannels(); zch++) {
62 fTracks[zch] = new TList[fGtuParam->GetNRefLayers()];
63 }
64 if (stack > -1)
65 SetStack(stack);
66 if (sector > -1)
67 SetSector(sector);
68}
69
70AliTRDgtuTMU::~AliTRDgtuTMU()
71{
72 for (Int_t zch = 0; zch < fGtuParam->GetNZChannels(); zch++) {
73 delete [] fTracks[zch];
74 }
75 delete [] fTracks;
76 for (Int_t layer = 0; layer < fGtuParam->GetNLayers(); layer++) {
77 fTracklets[layer]->Delete();
78 delete [] fZChannelTracklets[layer];
79 delete fTracklets[layer];
80 }
81 delete [] fZChannelTracklets;
82 delete [] fTracklets;
83}
84
85Bool_t AliTRDgtuTMU::SetSector(Int_t sector)
86{
87 // set the sector
88
89 if (sector > -1 && sector < fGtuParam->GetGeo()->Nsector() ) {
90 fSector = sector;
91 return kTRUE;
92 }
93
94 AliError(Form("Invalid sector given: %i", sector));
95 return kFALSE;
96}
97
98Bool_t AliTRDgtuTMU::SetStack(Int_t stack)
99{
100 // Set the stack (necessary for tracking)
101
102 if (stack > -1 && stack < fGtuParam->GetGeo()->Nstack() ) {
103 fStack = stack;
104 return kTRUE;
105 }
106
107 AliError(Form("Invalid stack given: %i", stack));
108 return kFALSE;
109}
110
111Bool_t AliTRDgtuTMU::Reset()
112{
113 // delete all tracks
114
115 for (Int_t zch = 0; zch < fGtuParam->GetNZChannels(); zch++) {
116 for (Int_t reflayeridx = 0; reflayeridx < fGtuParam->GetNRefLayers(); reflayeridx++) {
117 fTracks[zch][reflayeridx].Clear();
118 }
119 }
120
121 // delete all added tracklets
122 for (Int_t layer = 0; layer < fGtuParam->GetNLinks()/2; layer++) {
123 for (Int_t zch = 0; zch < fGtuParam->GetNZChannels(); zch++)
124 fZChannelTracklets[layer][zch].Clear();
125 fTracklets[layer]->Delete();
126 }
127
128 fSector = -1;
129 fStack = -1;
130
131 return kTRUE;
132}
133
134Bool_t AliTRDgtuTMU::AddTracklet(AliTRDtrackletBase *tracklet, Int_t link)
135{
136 // add a tracklet on the given link
137
138 AliTRDtrackletGTU *trkl = new AliTRDtrackletGTU(tracklet);
139 fTracklets[(Int_t) link/2]->Add(trkl);
140 return kTRUE;
141}
142
143
144Bool_t AliTRDgtuTMU::RunTMU(TList *ListOfTracks, AliESDEvent *esd)
145{
146 // performs the analysis as in a TMU module of the GTU, i. e.
147 // track matching
148 // calculation of track parameteres (pt, deflection, ???)
149
150 if (fStack < 0 || fSector < 0) {
151 AliError("No valid stack/sector set for this TMU! Tracking aborted!");
152 return kFALSE;
153 }
154
155 // ----- Input units -----
c8b1590d 156 AliDebug(1,"--------- Running Input units ----------");
52c19022 157 for (Int_t layer = 0; layer < fGtuParam->GetNLayers(); layer++) {
158 if (!RunInputUnit(layer)) {
159 AliError(Form("Input unit in layer %i failed", layer));
160 return kFALSE;
161 }
162 }
163
164 // ----- Z-channel units -----
c8b1590d 165 AliDebug(1,"--------- Running Z-channel units ----------");
52c19022 166 for (Int_t layer = 0; layer < fGtuParam->GetNLayers(); layer++) {
167 fZChannelTracklets[layer] = new TList[fGtuParam->GetNZChannels()];
168 if (!RunZChannelUnit(layer)) {
169 AliError(Form("Z-Channel unit in layer %i failed", layer));
170 return kFALSE;
171 }
172 }
173
174 // ----- track finding -----
c8b1590d 175 AliDebug(1,"--------- Running tracking units ----------");
52c19022 176 for (Int_t zch = 0; zch < fGtuParam->GetNZChannels(); zch++) {
c8b1590d 177 AliDebug(2,Form("Track finder for Zchannel: %i", zch));
52c19022 178 if (!RunTrackFinder(zch, ListOfTracks)) {
179 AliError(Form("Track Finder in z-channel %i failed", zch));
180 return kFALSE;
181 }
182 }
183
184 // ----- Track Merging -----
185 if (!RunTrackMerging(ListOfTracks)) {
186 AliError("Track merging failed");
187 return kFALSE;
188 }
189
190 // ----- track reconstruction -----
191 if (!RunTrackReconstruction(ListOfTracks)) {
192 AliError("Track reconstruction failed");
193 return kFALSE;
194 }
195
196 if (esd) {
197 TIter next(ListOfTracks);
198 while (AliTRDtrackGTU *trk = (AliTRDtrackGTU*) next()) {
199 AliESDTrdTrack *trdtrack = trk->CreateTrdTrack();
200 esd->AddTrdTrack(trdtrack);
201 delete trdtrack;
202 }
203 }
204
205 return kTRUE;
206}
207
208Bool_t AliTRDgtuTMU::RunInputUnit(Int_t layer)
209{
210 // precalculations for the tracking and reconstruction
211
212 fTracklets[layer]->Sort();
213 TIter next(fTracklets[layer]);
214
215 while ( AliTRDtrackletGTU *trk = (AliTRDtrackletGTU*) next() ) {
216 trk->SetIndex(fTracklets[layer]->IndexOf(trk));
217
218 Int_t alpha = (trk->GetYbin() >> fGtuParam->GetBitExcessY()) * fGtuParam->GetCiAlpha(layer);
219 alpha = ( 2 * trk->GetdY() - (alpha >> fGtuParam->GetBitExcessAlpha()) + 1 ) >> 1;
220 trk->SetAlpha(alpha);
221
80f93426 222 Int_t yproj = -1 * trk->GetdY() * (fGtuParam->GetCiYProj(layer)); //??? sign?
52c19022 223 yproj = ( ( ( (yproj >> fGtuParam->GetBitExcessYProj()) + trk->GetYbin() ) >> 2) + 1) >> 1;
224 trk->SetYProj(yproj);
225
226 trk->SetYPrime(trk->GetYbin() + fGtuParam->GetYt(fStack, layer, trk->GetZbin()));
227
228// printf("InputUnit : GetIndex(): %3i, GetZbin(): %2i, GetY() : %5i, GetdY() : %3i, GetYPrime() : %5i, GetYProj() : %5i, GetAlpha() : %3i \n",
229// trk->GetIndex(), trk->GetZbin(), trk->GetYbin(), trk->GetdY(), trk->GetYPrime(), trk->GetYProj(), trk->GetAlpha() );
230 }
231 return kTRUE;
232}
233
234Bool_t AliTRDgtuTMU::RunZChannelUnit(Int_t layer)
235{
236 // run the z-channel unit
237
238 TIter next(fTracklets[layer]);
239
240 while (AliTRDtrackletGTU *trk = (AliTRDtrackletGTU*) next()) {
c8b1590d 241// printf("*TMU* Tracklet in stack %d, layer %2d: 0x%08x ", fStack, layer, trk->GetTrackletWord());
52c19022 242 for (Int_t zch = 0; zch < fGtuParam->GetNZChannels(); zch++) {
243 if (fGtuParam->IsInZChannel(fStack, layer, zch, trk->GetZbin()) ) {
244 trk->SetSubChannel(zch, fGtuParam->GetZSubchannel(fStack, layer, zch, trk->GetZbin()) );
c8b1590d 245// printf("Z%i(%i) ", zch, trk->GetSubChannel(zch));
52c19022 246
247 TIter nexttrkl(&fZChannelTracklets[layer][zch], kIterBackward);
248 AliTRDtrackletGTU *t = 0x0;
4cc89512 249 while ((t = (AliTRDtrackletGTU*) nexttrkl.Next())) {
52c19022 250 if (t->GetSubChannel(zch) < trk->GetSubChannel(zch) ||
251 (t->GetSubChannel(zch) == trk->GetSubChannel(zch) && t->GetYProj() < trk->GetYProj()) )
252 break;
253 }
254 fZChannelTracklets[layer][zch].AddAfter(t, trk);
255 }
c8b1590d 256// else
257// printf(" ");
52c19022 258 }
c8b1590d 259// printf("\n");
52c19022 260 }
261 return kTRUE;
262}
263
4cc89512 264Bool_t AliTRDgtuTMU::RunTrackFinder(Int_t zch, TList* /* ListOfTracks */)
52c19022 265{
266 // run the track finding
267
268 Int_t *notr = new Int_t[fGtuParam->GetNLayers()];
269 Int_t *ptrA = new Int_t[fGtuParam->GetNLayers()];
270 Int_t *ptrB = new Int_t[fGtuParam->GetNLayers()];
271
272 Bool_t *bHitA = new Bool_t[fGtuParam->GetNLayers()];
273 Bool_t *bHitB = new Bool_t[fGtuParam->GetNLayers()];
274 Bool_t *bAligned = new Bool_t[fGtuParam->GetNLayers()];
275 Bool_t *bAlignedA = new Bool_t[fGtuParam->GetNLayers()];
276 Bool_t *bAlignedB = new Bool_t[fGtuParam->GetNLayers()];
277 Bool_t *bDone = new Bool_t[fGtuParam->GetNLayers()];
278 Bool_t ready;
279
280 Int_t *inc = new Int_t[fGtuParam->GetNLayers()];
281 Int_t *incprime = new Int_t[fGtuParam->GetNLayers()];
282
283// ----- signals within current layer -----
284 Int_t Yplus;
285 Int_t Yminus;
286 Int_t YBplus;
287 Int_t YBminus;
288 Int_t Alphaplus;
289 Int_t Alphaminus;
290 Int_t NHits;
291 Int_t NUnc;
292 Int_t NWayBeyond;
293
294 AliTRDtrackletGTU *trkRA = 0x0; // reference tracklet A
295 AliTRDtrackletGTU *trkRB = 0x0; // reference tracklet B
296 AliTRDtrackletGTU *trkA = 0x0; // tracklet A in current layer
297 AliTRDtrackletGTU *trkB = 0x0; // tracklet B in current layer
298/*
299 AliTRDtrackletGTU *trkEnd = new AliTRDtrackletGTU();
300 for (Int_t i = 0; i < fGtuParam->GetNZChannels(); i++)
301 trkEnd->SetSubChannel(i, 7);
302 trkEnd->SetYProj(0);
303 trkEnd->SetAlpha(0);
304*/
305
306 for (Int_t refLayerIdx = 0; refLayerIdx < fGtuParam->GetNRefLayers(); refLayerIdx++) {
307 Int_t reflayer = fGtuParam->GetRefLayer(refLayerIdx);
c8b1590d 308 AliDebug(5,Form("~~~~~ Reflayer: %i", reflayer));
52c19022 309
310 ready = kFALSE; // ready if all channels done
311
312 for (Int_t layer = 0; layer < fGtuParam->GetNLayers(); layer++) {
313 notr[layer] = fZChannelTracklets[layer][zch].GetEntries();
314 ptrA[layer] = 0; // notr[layer] > 0 ? 0 : -1;
315 ptrB[layer] = 1; // notr[layer] > 1 ? 1 : -1;
316
317// not necessary here
318// bDone[layer] = (ptrB[layer] >= notr[layer] - 1); // trkB is last one
319// bDone[layer] = (notr[layer] <= 0);
320// ready = ready && bDone[layer];
321
322 if (reflayer == 1)
c8b1590d 323 AliDebug(5,Form("in layer: %i (zchannel = %i) there are: %i tracklets", layer, zch, notr[layer]));
52c19022 324 }
325
326 if (ptrA[reflayer] < 0 && ptrB[reflayer] < 0)
327 continue;
328
329 while (!ready) {
330 // ----- reference tracklets -----
331 trkRA = 0x0;
332 trkRB = 0x0;
333 if (0 <= ptrA[reflayer] && ptrA[reflayer] < notr[reflayer])
334 trkRA = (AliTRDtrackletGTU*) fZChannelTracklets[reflayer][zch].At(ptrA[reflayer]);
335 else {
c8b1590d 336 AliDebug(10,Form("No valid tracklet in the reference at ptr: %i! Aborting!", ptrA[reflayer]));
52c19022 337 break;
338 }
339
340 if (0 <= ptrB[reflayer] && ptrB[reflayer] < notr[reflayer])
341 trkRB = (AliTRDtrackletGTU*) fZChannelTracklets[reflayer][zch].At(ptrB[reflayer]);
342
c8b1590d 343 AliDebug(10,Form("ptrRA: %i, ptrRB: %i", ptrA[reflayer], ptrB[reflayer]));
52c19022 344 Yplus = trkRA->GetYProj() + fGtuParam->GetDeltaY();
345 Yminus = trkRA->GetYProj() - fGtuParam->GetDeltaY();
346 Alphaplus = trkRA->GetAlpha() + fGtuParam->GetDeltaAlpha();
347 Alphaminus = trkRA->GetAlpha() - fGtuParam->GetDeltaAlpha();
348 if (trkRB) {
349 YBplus = trkRB->GetYProj() + fGtuParam->GetDeltaY();
350 YBminus = trkRB->GetYProj() - fGtuParam->GetDeltaY();
351 }
352 else { // irrelevant (should be, is it?)
353 YBplus = trkRA->GetYProj() + fGtuParam->GetDeltaY();
354 YBminus = trkRA->GetYProj() - fGtuParam->GetDeltaY();
355 }
356
357 NHits = 0;
358 NUnc = 0;
359 NWayBeyond = 0;
360
361 for (Int_t layer = 0; layer < fGtuParam->GetNLayers(); layer++) {
362 bHitA[layer] = bHitB[layer] = bAligned[layer] = kFALSE;
363 inc[layer] = incprime[layer] = 0;
364
365 if (layer == reflayer) {
366 bHitA[layer] = kTRUE;
367 bAligned[layer] = kTRUE;
368 NHits++;
369 continue;
370 }
371
372 trkA = 0x0;
373 trkB = 0x0;
374 if (0 <= ptrA[layer] && ptrA[layer] < notr[layer])
375 trkA = (AliTRDtrackletGTU*) fZChannelTracklets[layer][zch].At(ptrA[layer]);
376 if (0 <= ptrB[layer] && ptrB[layer] < notr[layer])
377 trkB = (AliTRDtrackletGTU*) fZChannelTracklets[layer][zch].At(ptrB[layer]);
378
379 bAlignedA[layer] = kFALSE;
380 bAlignedB[layer] = kFALSE;
381
382 if (trkA) {
383 bHitA[layer] = ( !(trkA->GetSubChannel(zch) < trkRA->GetSubChannel(zch) || (trkA->GetSubChannel(zch) == trkRA->GetSubChannel(zch) && trkA->GetYProj() < Yminus) ) &&
384 !(trkA->GetSubChannel(zch) > trkRA->GetSubChannel(zch) || (trkA->GetSubChannel(zch) == trkRA->GetSubChannel(zch) && trkA->GetYProj() > Yplus) ) &&
385 !(trkA->GetAlpha() < Alphaminus) &&
386 !(trkA->GetAlpha() > Alphaplus) );
387 bAlignedA[layer] = !(trkA->GetSubChannel(zch) < trkRA->GetSubChannel(zch) || (trkA->GetSubChannel(zch) == trkRA->GetSubChannel(zch) && trkA->GetYProj() < Yminus) );
388 }
389 else {
390 bHitA[layer] = 0;
391 bAlignedA[layer] = kTRUE;
392 }
393
394 if (trkB) {
395 bHitB[layer] = ( !(trkB->GetSubChannel(zch) < trkRA->GetSubChannel(zch) || (trkB->GetSubChannel(zch) == trkRA->GetSubChannel(zch) && trkB->GetYProj() < Yminus) ) &&
396 !(trkB->GetSubChannel(zch) > trkRA->GetSubChannel(zch) || (trkB->GetSubChannel(zch) == trkRA->GetSubChannel(zch) && trkB->GetYProj() > Yplus) ) &&
397 !(Alphaminus > trkB->GetAlpha()) &&
398 !(Alphaplus > trkB->GetAlpha()) );
399 bAlignedB[layer] = (trkB->GetSubChannel(zch) > trkRA->GetSubChannel(zch) || (trkB->GetSubChannel(zch) == trkRA->GetSubChannel(zch) && trkB->GetYProj() > Yplus) );
400 }
401 else {
402 bHitB[layer] = 0;
403 bAlignedB[layer] = kTRUE;
404 }
405
406 bAligned[layer] = bAlignedA[layer] || bAlignedB[layer]; //???
407// bAligned[layer] = bAlignedA[layer]; //???
408
409 if (bAligned[layer] && (bHitA[layer] || bHitB[layer]) )
410 NHits++;
411 else if (!bAligned[layer] )
412 NUnc++;
413 if (trkRB) {
414 if (trkA) {
415 if ((trkA->GetSubChannel(zch) > trkRB->GetSubChannel(zch)) || (trkA->GetSubChannel(zch) == trkRB->GetSubChannel(zch) && trkA->GetYProj() > YBplus) )
416 NWayBeyond++;
417 }
418 else
419 NWayBeyond++;
420 }
421
422 // pre-calculation for the layer shifting (alignment w. r. t. trkRB)
423 if (trkA) {
424 if(trkRB) {
425 if ((trkA->GetSubChannel(zch) < trkRB->GetSubChannel(zch)) || (trkA->GetSubChannel(zch) == trkRB->GetSubChannel(zch) && trkA->GetYProj() < YBminus )) // could trkA be aligned for trkRB
426 incprime[layer] = 1;
427 else
428 incprime[layer] = 0;
429 }
430 else
431 incprime[layer] = 1;
432
433 if (trkB) {
434 if (trkRB) {
435 if ((trkB->GetSubChannel(zch) < trkRB->GetSubChannel(zch)) || (trkB->GetSubChannel(zch) == trkRB->GetSubChannel(zch) && trkB->GetYProj() < YBminus )) // could trkB be aligned for trkRB
436 incprime[layer] = 2;
437 }
438 else
439 incprime[layer] = 2;
440 }
441 }
442 } // end of loop over layers
443
c8b1590d 444 AliDebug(5,Form("logic calculation finished, Nhits: %i", NHits));
52c19022 445
446 if (NHits >= 4) {
447 // ----- track registration -----
c8b1590d 448 AliDebug(1,"***** TMU: Track found *****");
52c19022 449 AliTRDtrackGTU *track = new AliTRDtrackGTU();
450 track->SetSector(fSector);
451 track->SetStack(fStack);
452 for (Int_t layer = 0; layer < fGtuParam->GetNLayers(); layer++ ) {
453 if (bHitA[layer] || layer == reflayer)
454 track->AddTracklet((AliTRDtrackletGTU* ) fZChannelTracklets[layer][zch].At(ptrA[layer]), layer );
455 else if (bHitB[layer])
456 track->AddTracklet((AliTRDtrackletGTU* ) fZChannelTracklets[layer][zch].At(ptrB[layer]), layer );
457 }
458
459 Bool_t registerTrack = kTRUE;
460 for (Int_t layerIdx = refLayerIdx; layerIdx > 0; layerIdx--) {
80f93426 461 if (track->IsTrackletInLayer(fGtuParam->GetRefLayer(layerIdx))) {
462 if ((track->GetTracklet(fGtuParam->GetRefLayer(layerIdx)))->GetSubChannel(zch) > 0) {
c8b1590d 463 AliDebug(1,"Not registered");
52c19022 464 registerTrack = kFALSE;
80f93426 465 }
466 }
52c19022 467 }
468 if (registerTrack) {
469 track->SetZChannel(zch);
470 track->SetRefLayerIdx(refLayerIdx);
471 fTracks[zch][refLayerIdx].Add(track);
472 }
473 }
474
475 if ( (NUnc != 0) && (NUnc + NHits >= 4) ) // could this position of the reference layer give some track //??? special check in case of hit?
476 inc[reflayer] = 0;
477 else if (NWayBeyond > 2) // no track possible for both reference tracklets
478 inc[reflayer] = 2;
479 else
480 inc[reflayer] = 1;
481
482 if (inc[reflayer] != 0) // reflayer is shifted
483 for (Int_t layer = 0; layer < fGtuParam->GetNLayers(); layer++) {
484 if (layer == reflayer)
485 continue;
486 inc[layer] = incprime[layer];
487 }
488 else { // reflayer is not shifted
489 for (Int_t layer = 0; layer < fGtuParam->GetNLayers(); layer++) {
490 if (layer == reflayer || notr[layer] == 0)
491 continue;
492 inc[layer] = 0;
493 trkA = 0x0;
494 trkB = 0x0;
495 if (0 <= ptrA[layer] && ptrA[layer] < notr[layer])
496 trkA = (AliTRDtrackletGTU*) fZChannelTracklets[layer][zch].At(ptrA[layer]);
497
498 if (0 <= ptrB[layer] && ptrB[layer] < notr[layer])
499 trkB = (AliTRDtrackletGTU*) fZChannelTracklets[layer][zch].At(ptrB[layer]);
500
501 if (trkA) {
502 if ( !(trkA->GetSubChannel(zch) < trkRA->GetSubChannel(zch) || (trkA->GetSubChannel(zch) == trkRA->GetSubChannel(zch) && trkA->GetYProj() < Yminus) ) &&
503 !(trkA->GetSubChannel(zch) > trkRA->GetSubChannel(zch) || (trkA->GetSubChannel(zch) == trkRA->GetSubChannel(zch) && trkA->GetYProj() > Yplus ) ) ) // trkA could hit trkRA
504 inc[layer] = 0;
505 else if (trkB) {
506 if ( trkB->GetSubChannel(zch) < trkRA->GetSubChannel(zch) || (trkB->GetSubChannel(zch) == trkRA->GetSubChannel(zch) && trkB->GetYProj() < Yminus) )
507 inc[layer] = 2;
508 else if ( !(trkB->GetSubChannel(zch) > trkRA->GetSubChannel(zch) || (trkB->GetSubChannel(zch) == trkRA->GetSubChannel(zch) && trkB->GetYProj() > Yplus) ) )
509 inc[layer] = 1;
510 else
511 inc[layer] = incprime[layer];
512 }
513 else
514 inc[layer] = incprime[layer];
515 }
516 }
517 }
518
519 ready = kTRUE;
520 for (Int_t layer = 0; layer < fGtuParam->GetNLayers(); layer++) {
521
522 bDone[layer] = ptrB[layer] < 0 || ptrB[layer] >= notr[layer];
523 ready = ready && bDone[layer];
524
525 if (inc[layer] != 0 && ptrA[layer] >= notr[layer])
526 AliError(Form("Invalid increment: %i at ptrA: %i, notr: %i", inc[layer], ptrA[layer], notr[layer]));
527
528// AliInfo(Form("Shifting layer: %i, notr: %i, ptrA: %i, ptrB: %i, inc: %i", layer, notr[layer], ptrA[layer], ptrB[layer], inc[layer]));
c8b1590d 529 AliDebug(10,Form(" -- Layer: %i %i %i +%i %s (no: %i)", layer, ptrA[layer], ptrB[layer], inc[layer], bDone[layer] ? "done" : " ", notr[layer]));
52c19022 530 ptrA[layer] += inc[layer];
531 ptrB[layer] += inc[layer];
532 }
533
534 } // end of while
535 } // end of loop over reflayer
536
537 delete [] bHitA;
538 delete [] bHitB;
539 delete [] bAligned;
540 delete [] bDone;
541 delete [] inc;
542 delete [] incprime;
543 delete [] bAlignedA;
544 delete [] bAlignedB;
545 delete [] notr;
546 delete [] ptrA;
547 delete [] ptrB;
548
549 return kTRUE;
550}
551
552Bool_t AliTRDgtuTMU::RunTrackMerging(TList* ListOfTracks)
553{
554 TList **tracksRefMerged = new TList*[fGtuParam->GetNZChannels()];
555 TList **tracksRefUnique = new TList*[fGtuParam->GetNZChannels()];
556
557 for (Int_t zch = 0; zch < fGtuParam->GetNZChannels(); zch++) {
558 tracksRefMerged[zch] = new TList;
559 tracksRefUnique[zch] = new TList;
560 }
561
562 TList *tracksZMergedStage0 = new TList;
563 TList *tracksZUniqueStage0 = new TList;
564
565 TList **tracksZSplitted = new TList*[2];
566 for (Int_t i = 0; i < 2; i++)
567 tracksZSplitted[i] = new TList;
568
569 TList *tracksZMergedStage1 = new TList;
570
4cc89512 571 AliTRDtrackGTU **trkInRefLayer = new AliTRDtrackGTU*[fGtuParam->GetNRefLayers()];
52c19022 572
573 Bool_t done = kFALSE;
574 Int_t minIdx = 0;
575 AliTRDtrackGTU *trkStage0 = 0x0;
576
577 for (Int_t zch = 0; zch < fGtuParam->GetNZChannels(); zch++) {
578 // ----- Merging and Unification in Reflayers -----
579 do {
580 done = kTRUE;
581 trkStage0 = 0x0;
582 for (Int_t refLayerIdx = 0; refLayerIdx < fGtuParam->GetNRefLayers(); refLayerIdx++) {
4cc89512 583 trkInRefLayer[refLayerIdx] = (AliTRDtrackGTU*) fTracks[zch][refLayerIdx].First();
584 if (trkInRefLayer[refLayerIdx] == 0) {
52c19022 585 continue;
586 }
587 else if (trkStage0 == 0x0 ) {
4cc89512 588 trkStage0 = trkInRefLayer[refLayerIdx];
52c19022 589 minIdx = refLayerIdx;
590 done = kFALSE;
591 }
4cc89512 592 else if (trkInRefLayer[refLayerIdx]->GetZSubChannel() < trkStage0->GetZSubChannel() ||
593 (trkInRefLayer[refLayerIdx]->GetZSubChannel() == trkStage0->GetZSubChannel() && trkInRefLayer[refLayerIdx]->GetYapprox() < trkStage0->GetYapprox()) ) {
52c19022 594 minIdx = refLayerIdx;
4cc89512 595 trkStage0 = trkInRefLayer[refLayerIdx];
52c19022 596 done = kFALSE;
597 }
598 }
c8b1590d 599 if (!trkStage0)
600 break;
52c19022 601 tracksRefMerged[zch]->Add(trkStage0);
602 fTracks[zch][minIdx].RemoveFirst();
603 } while (trkStage0 != 0);
604
605 Uniquifier(tracksRefMerged[zch], tracksRefUnique[zch]);
606 }
607
608// ----- Merging in zchannels - 1st stage -----
609
610 do {
611 done = kTRUE;
612 trkStage0 = 0x0;
613 for (Int_t zch = 0; zch < fGtuParam->GetNZChannels(); zch++) {
614 AliTRDtrackGTU *trk = (AliTRDtrackGTU*) tracksRefUnique[zch]->First();
615 if (trk == 0) {
616 continue;
617 }
618 else if (trkStage0 == 0x0 ) {
619 trkStage0 = trk;
620 minIdx = zch;
621 done = kFALSE;
622 }
623 else if ( ((trk->GetZChannel() + 3 * trk->GetZSubChannel()) / 2 - 1) < ((trkStage0->GetZChannel() + 3 * trkStage0->GetZSubChannel()) / 2 -1 ) ) {
624 minIdx = zch;
625 trkStage0 = trk;
626 done = kFALSE;
627 }
628 }
629
c8b1590d 630 if (!trkStage0)
631 break;
52c19022 632 tracksZMergedStage0->Add(trkStage0);
633 tracksRefUnique[minIdx]->RemoveFirst();
634 } while (trkStage0 != 0);
635
636 Uniquifier(tracksZMergedStage0, tracksZUniqueStage0);
637
638// ----- Splitting in z -----
639
640 TIter next(tracksZUniqueStage0);
641 while (AliTRDtrackGTU *trk = (AliTRDtrackGTU*) next()) {
642 tracksZSplitted[(trk->GetZChannel() + 3 * (trk->GetZSubChannel() - 1)) % 2]->Add(trk);
643 }
644
645// ----- Merging in zchanels - 2nd stage -----
646
647 do {
648 done = kTRUE;
649 trkStage0 = 0x0;
650 for (Int_t i = 0; i < 2; i++) {
651 AliTRDtrackGTU *trk = (AliTRDtrackGTU*) tracksZSplitted[i]->First();
652 if (trk == 0) {
653 continue;
654 }
655 else if (trkStage0 == 0x0 ) {
656 trkStage0 = trk;
657 minIdx = i;
658 done = kFALSE;
659 }
660 else if ( ((trk->GetZChannel() + 3 * (trk->GetZSubChannel() - 1)) / 2) < ((trkStage0->GetZChannel() + 3 * (trkStage0->GetZSubChannel() - 1)) / 2) ) {
661 minIdx = i;
662 trkStage0 = trk;
663 done = kFALSE;
664 }
665 }
666
c8b1590d 667 if (!trkStage0)
668 break;
52c19022 669 tracksZMergedStage1->Add(trkStage0);
670 tracksZSplitted[minIdx]->RemoveFirst();
671 } while (trkStage0 != 0);
672
673 Uniquifier(tracksZMergedStage1, ListOfTracks);
674
675 return kTRUE;
676}
677
678Bool_t AliTRDgtuTMU::RunTrackReconstruction(TList* ListOfTracks)
679{
680 TIter next(ListOfTracks);
681 while (AliTRDtrackGTU *track = (AliTRDtrackGTU*) next()) {
682 CalculateTrackParams(track);
683 }
684 return kTRUE;
685}
686
687Bool_t AliTRDgtuTMU::CalculateTrackParams(AliTRDtrackGTU *track)
688{
689 // calculate the track parameters
690
691 if (!track) {
692 AliError("No track to calculate!");
693 return kFALSE;
694 }
695
696 Int_t a = 0;
697 Float_t b = 0;
698 Float_t c = 0;
699 Float_t x1;
700 Float_t x2;
701
c8b1590d 702 AliDebug(5,Form("There are %i tracklets in this track.", track->GetNTracklets()));
52c19022 703
704 for (Int_t layer = 0; layer < fGtuParam->GetNLayers(); layer++) {
705 if (!track->IsTrackletInLayer(layer)) {
706 continue;
707 }
708 AliTRDtrackletGTU *trk = track->GetTracklet(layer);
709 if (!trk) {
710 AliError(Form("Could not get tracklet in layer %i\n", layer));
711 continue;
712 }
c8b1590d 713 AliDebug(10,Form("trk yprime: %i", trk->GetYPrime()));
52c19022 714 a += (((Int_t) (2048 * fGtuParam->GetAki(track->GetTrackletMask(), layer))) * trk->GetYPrime() + 1) >> 8;
715 b += fGtuParam->GetBki(track->GetTrackletMask(), layer) * trk->GetYPrime() * fGtuParam->GetBinWidthY();
716 c += fGtuParam->GetCki(track->GetTrackletMask(), layer) * trk->GetYPrime() * fGtuParam->GetBinWidthY();
717 }
718 if (a < 0)
719 a += 3;
720 a = a >> 2;
721
722 fGtuParam->GetIntersectionPoints(track->GetTrackletMask(), x1, x2);
c8b1590d 723 AliDebug(10,Form("Intersection points: %f, %f", x1, x2));
724 AliDebug(10,Form("Sum: a = %5i, b = %9.2f, c = %9.2f\n", a, b, c));
52c19022 725 track->SetFitParams(a, b, c);
726
727 Float_t r = fGtuParam->GetRadius(a, b, x1, x2);
728 Int_t pt = (Int_t) (2 * r);
729 if (pt >= 0)
730 pt += 32;
731 else
732 pt -= 29;
733 pt /= 2;
734 track->SetPtInt(pt);
c8b1590d 735 AliDebug(5,Form("Track parameters: a = %i, b = %f, c = %f, x1 = %f, x2 = %f, r = %f, pt = %f (trkl mask: %i)", a, b, c, x1, x2, r, track->GetPt(), track->GetTrackletMask()));
52c19022 736 return kTRUE;
737}
738
739Bool_t AliTRDgtuTMU::WriteTrackletsToTree(TTree *trklTree)
740{
741 if (!trklTree) {
742 AliError("No tree given");
743 return kFALSE;
744 }
745 AliTRDtrackletGTU *trkl = 0x0;
746 TBranch *branch = trklTree->GetBranch("gtutracklets");
747 if (!branch) {
748 branch = trklTree->Branch("gtutracklets", "AliTRDtrackletGTU", &trkl, 32000, 99);
749 }
750
c8b1590d 751 AliDebug(5,Form("---------- Writing tracklets to tree (not yet) ----------"));
52c19022 752 for (Int_t layer = 0; layer < fGtuParam->GetNLayers(); layer++) {
753 TIter next(fTracklets[layer]);
edf78ab4 754 while ((trkl = (AliTRDtrackletGTU*) next())) {
c8b1590d 755 AliDebug(10,Form("InputUnit : GetIndex(): %3i, GetZbin(): %2i, GetY() : %5i, GetdY() : %3i, GetYPrime() : %5i, GetYProj() : %5i, GetAlpha() : %3i, Zidx(2..0): %i %i %i", trkl->GetIndex(), trkl->GetZbin(), trkl->GetYbin(), trkl->GetdY(), trkl->GetYPrime(), trkl->GetYProj(), trkl->GetAlpha(), trkl->GetSubChannel(2), trkl->GetSubChannel(1), trkl->GetSubChannel(0) ));
52c19022 756 branch->SetAddress(&trkl);
757 trklTree->Fill();
758 }
759 }
760 return kTRUE;
761}
762
763Bool_t AliTRDgtuTMU::Uniquifier(TList *inlist, TList *outlist)
764{
765 TIter next(inlist);
766 AliTRDtrackGTU *trkStage0 = 0x0;
767 AliTRDtrackGTU *trkStage1 = 0x0;
768
769 do {
770 trkStage0 = (AliTRDtrackGTU*) next();
771
772 Bool_t tracks_equal = kFALSE;
773 if (trkStage0 != 0 && trkStage1 != 0) {
774 for (Int_t layer = 0; layer < fGtuParam->GetNLayers(); layer++) {
775 if (trkStage0->GetTrackletIndex(layer) != -1 && trkStage0->GetTrackletIndex(layer) == trkStage1->GetTrackletIndex(layer)) {
776 tracks_equal = kTRUE;
777 break;
778 }
779 }
780 }
781
782 if (tracks_equal) {
783 if (trkStage0->GetNTracklets() > trkStage1->GetNTracklets())
784 trkStage1 = trkStage0;
785 }
786 else {
787 if (trkStage1 != 0x0)
788 outlist->Add(trkStage1);
789 trkStage1 = trkStage0;
790 }
791
792 } while (trkStage1 != 0x0);
793 return kTRUE;
794}