Weird inline removed
[u/mrichter/AliRoot.git] / TRD / AliTRDtracker.cxx
CommitLineData
46d29e70 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/*
17$Log$
bdbd0f7a 18Revision 1.10 2001/05/07 08:08:05 cblume
19Update of TRD code
20
a2b90f83 21Revision 1.9 2001/02/14 18:22:26 cblume
22Change in the geometry of the padplane
23
71d9fa7b 24Revision 1.8 2000/12/20 13:00:44 cblume
25Modifications for the HP-compiler
26
d1b06c24 27Revision 1.7 2000/12/08 16:07:02 cblume
28Update of the tracking by Sergei
29
bbf92647 30Revision 1.6 2000/11/30 17:38:08 cblume
31Changes to get in line with new STEER and EVGEN
32
1819f4bb 33Revision 1.5 2000/11/14 14:40:27 cblume
34Correction for the Sun compiler (kTRUE and kFALSE)
35
57527628 36Revision 1.4 2000/11/10 14:57:52 cblume
37Changes in the geometry constants for the DEC compiler
38
dd56b762 39Revision 1.3 2000/10/15 23:40:01 cblume
40Remove AliTRDconst
41
0e9c2ad5 42Revision 1.2 2000/10/06 16:49:46 cblume
43Made Getters const
44
46d29e70 45Revision 1.1.2.2 2000/10/04 16:34:58 cblume
46Replace include files by forward declarations
47
48Revision 1.1.2.1 2000/09/22 14:47:52 cblume
49Add the tracking code
50
51*/
bbf92647 52
1819f4bb 53#include <iostream.h>
46d29e70 54
55#include <TFile.h>
56#include <TROOT.h>
57#include <TBranch.h>
58#include <TTree.h>
59
60#include "AliRun.h"
46d29e70 61#include "AliTRD.h"
46d29e70 62#include "AliTRDgeometry.h"
63#include "AliTRDrecPoint.h"
64#include "AliTRDcluster.h"
65#include "AliTRDtrack.h"
66#include "AliTRDtrackingSector.h"
67#include "AliTRDtimeBin.h"
68
69#include "AliTRDtracker.h"
70
71ClassImp(AliTRDtracker)
72
bbf92647 73 const Int_t AliTRDtracker::fSeedGap = 35;
74 const Int_t AliTRDtracker::fSeedStep = 5;
75
76
77 const Float_t AliTRDtracker::fMinClustersInTrack = 0.5;
78 const Float_t AliTRDtracker::fMinClustersInSeed = 0.5;
79 const Float_t AliTRDtracker::fSeedDepth = 0.5;
80 const Float_t AliTRDtracker::fSkipDepth = 0.2;
81 const Float_t AliTRDtracker::fMaxSeedDeltaZ = 30.;
82 const Float_t AliTRDtracker::fMaxSeedC = 0.01;
83 const Float_t AliTRDtracker::fMaxSeedTan = 1.2;
84 const Float_t AliTRDtracker::fMaxSeedVertexZ = 200.;
85 const Float_t AliTRDtracker::fLabelFraction = 0.5;
86 const Float_t AliTRDtracker::fWideRoad = 30.;
87
88 const Double_t AliTRDtracker::fMaxChi2 = 12.;
89 const Double_t AliTRDtracker::fSeedErrorSY = 0.1;
90 const Double_t AliTRDtracker::fSeedErrorSY3 = 2.5;
91 const Double_t AliTRDtracker::fSeedErrorSZ = 0.1;
92
46d29e70 93//____________________________________________________________________
94AliTRDtracker::AliTRDtracker()
95{
96 //
97 // Default constructor
98 //
99
46d29e70 100 fEvent = 0;
46d29e70 101 fGeom = NULL;
bbf92647 102
46d29e70 103 fNclusters = 0;
104 fClusters = NULL;
105 fNseeds = 0;
106 fSeeds = NULL;
107 fNtracks = 0;
108 fTracks = NULL;
109
110}
111
112//____________________________________________________________________
113AliTRDtracker::AliTRDtracker(const Text_t* name, const Text_t* title)
114 :TNamed(name, title)
115{
46d29e70 116 fEvent = 0;
46d29e70 117 fGeom = NULL;
bbf92647 118
46d29e70 119 fNclusters = 0;
120 fClusters = new TObjArray(2000);
121 fNseeds = 0;
122 fSeeds = new TObjArray(20000);
123 fNtracks = 0;
124 fTracks = new TObjArray(10000);
125
126}
127
46d29e70 128//___________________________________________________________________
129AliTRDtracker::~AliTRDtracker()
130{
46d29e70 131 delete fClusters;
132 delete fTracks;
133 delete fSeeds;
134 delete fGeom;
135}
136
bbf92647 137//___________________________________________________________________
138void AliTRDtracker::Clusters2Tracks()
139{
140 Int_t inner, outer;
141 Int_t fTotalNofTimeBins = fGeom->GetTimeMax() * AliTRDgeometry::Nplan();
142 Int_t nSteps = (Int_t) (fTotalNofTimeBins * fSeedDepth)/fSeedStep;
143
144 for(Int_t i=0; i<nSteps; i++) {
145 printf("step %d out of %d \n", i+1, nSteps);
146 outer=fTotalNofTimeBins-1-i*fSeedStep; inner=outer-fSeedGap;
147 MakeSeeds(inner,outer);
148 FindTracks();
149 }
150}
46d29e70 151
152//_____________________________________________________________________
bbf92647 153Double_t AliTRDtracker::ExpectedSigmaY2(Double_t r, Double_t tgl, Double_t pt)
46d29e70 154{
155 // Parametrised "expected" error of the cluster reconstruction in Y
156
bbf92647 157 Double_t s = 0.2;
46d29e70 158 return s;
159}
160
161//_____________________________________________________________________
bbf92647 162Double_t AliTRDtracker::ExpectedSigmaZ2(Double_t r, Double_t tgl)
46d29e70 163{
164 // Parametrised "expected" error of the cluster reconstruction in Z
165
71d9fa7b 166 // Take an example pad size for the time being, check back
167 // again with Sergei
168 Double_t s, pad = fGeom->GetRowPadSize(0,2,0);
46d29e70 169 s = pad * pad /12.;
170 return s;
171}
172
173//_____________________________________________________________________
bdbd0f7a 174Double_t f1trd(Double_t x1,Double_t y1,
46d29e70 175 Double_t x2,Double_t y2,
176 Double_t x3,Double_t y3)
177{
178 // Initial approximation of the track curvature
179 // Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch
180
181 Double_t d=(x2-x1)*(y3-y2)-(x3-x2)*(y2-y1);
182 Double_t a=0.5*((y3-y2)*(y2*y2-y1*y1+x2*x2-x1*x1)-
183 (y2-y1)*(y3*y3-y2*y2+x3*x3-x2*x2));
184 Double_t b=0.5*((x2-x1)*(y3*y3-y2*y2+x3*x3-x2*x2)-
185 (x3-x2)*(y2*y2-y1*y1+x2*x2-x1*x1));
186
187 Double_t xr=TMath::Abs(d/(d*x1-a)), yr=d/(d*y1-b);
188
189 return -xr*yr/sqrt(xr*xr+yr*yr);
190}
191
192//_____________________________________________________________________
bdbd0f7a 193Double_t f2trd(Double_t x1,Double_t y1,
46d29e70 194 Double_t x2,Double_t y2,
195 Double_t x3,Double_t y3)
196{
197 // Initial approximation of the track curvature times center of curvature
198 // Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch
199
200 Double_t d=(x2-x1)*(y3-y2)-(x3-x2)*(y2-y1);
201 Double_t a=0.5*((y3-y2)*(y2*y2-y1*y1+x2*x2-x1*x1)-
202 (y2-y1)*(y3*y3-y2*y2+x3*x3-x2*x2));
203 Double_t b=0.5*((x2-x1)*(y3*y3-y2*y2+x3*x3-x2*x2)-
204 (x3-x2)*(y2*y2-y1*y1+x2*x2-x1*x1));
205
206 Double_t xr=TMath::Abs(d/(d*x1-a)), yr=d/(d*y1-b);
207
208 return -a/(d*y1-b)*xr/sqrt(xr*xr+yr*yr);
209}
210
211//_____________________________________________________________________
bdbd0f7a 212Double_t f3trd(Double_t x1,Double_t y1,
46d29e70 213 Double_t x2,Double_t y2,
214 Double_t z1,Double_t z2)
215{
216 // Initial approximation of the tangent of the track dip angle
217 // Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch
218
219 return (z1 - z2)/sqrt((x1-x2)*(x1-x2)+(y1-y2)*(y1-y2));
220}
221
222
223//___________________________________________________________________
224
bbf92647 225Int_t AliTRDtracker::FindProlongation(AliTRDtrack& t, AliTRDtrackingSector *sec,
d1b06c24 226 Int_t s, Int_t rf)
46d29e70 227{
228 // Starting from current position on track=t this function tries
229 // to extrapolate the track up to timeBin=rf and to confirm prolongation
230 // if a close cluster is found. *sec is a pointer to allocated
231 // array of sectors, in which the initial sector has index=s.
232
bbf92647 233 const Int_t TIME_BINS_TO_SKIP=Int_t(fSkipDepth*sec->GetNtimeBins());
46d29e70 234 Int_t try_again=TIME_BINS_TO_SKIP;
235
0e9c2ad5 236 Double_t alpha=AliTRDgeometry::GetAlpha();
46d29e70 237
238 Int_t ns=Int_t(2*TMath::Pi()/alpha+0.5);
239
240 for (Int_t nr=sec->GetTimeBinNumber(t.GetX())-1; nr>=rf; nr--) {
241 Double_t x=sec->GetX(nr), ymax=sec->GetMaxY(nr);
242 if (!t.PropagateTo(x)) {
243 cerr<<"Can't propagate to x = "<<x<<endl;
244 return 0;
245 }
246
247 AliTRDcluster *cl=0;
248 UInt_t index=0;
249
bbf92647 250 Double_t max_chi2=fMaxChi2;
251
46d29e70 252 AliTRDtimeBin& time_bin=sec[s][nr];
bbf92647 253 Double_t sy2=ExpectedSigmaY2(t.GetX(),t.GetTgl(),t.GetPt());
254 Double_t sz2=ExpectedSigmaZ2(t.GetX(),t.GetTgl());
46d29e70 255 Double_t road=5.*sqrt(t.GetSigmaY2() + sy2), y=t.GetY(), z=t.GetZ();
256
bbf92647 257 if (road>fWideRoad) {
46d29e70 258 if (t.GetNclusters() > 4) {
259 cerr<<t.GetNclusters()<<" FindProlongation: Road is too wide !\n";
260 }
261 return 0;
262 }
263
264 if (time_bin) {
265 for (Int_t i=time_bin.Find(y-road); i<time_bin; i++) {
266 AliTRDcluster* c=(AliTRDcluster*)(time_bin[i]);
267 if (c->GetY() > y+road) break;
268 if (c->IsUsed() > 0) continue;
269
270 if((c->GetZ()-z)*(c->GetZ()-z) > 25. + sz2) continue;
271
272 Double_t chi2=t.GetPredictedChi2(c);
273
274 if (chi2 > max_chi2) continue;
275 max_chi2=chi2;
276 cl=c;
277 index=time_bin.GetIndex(i);
278 }
279 }
280 if (cl) {
281
282 // Float_t l=sec->GetPitch();
283 // t.SetSampledEdx(cl->fQ/l,Int_t(t));
284
285 t.Update(cl,max_chi2,index);
286
287 try_again=TIME_BINS_TO_SKIP;
288 } else {
289 if (try_again==0) break;
290 if (y > ymax) {
291 cerr<<"y > ymax: "<<y<<" > "<<ymax<<endl;
292 s = (s+1) % ns;
293 if (!t.Rotate(alpha)) {
294 cerr<<"Failed to rotate, alpha = "<<alpha<<endl;
295 return 0;
296 }
297 } else if (y <-ymax) {
298 cerr<<"y < -ymax: "<<y<<" < "<<-ymax<<endl;
299 s = (s-1+ns) % ns;
300 if (!t.Rotate(-alpha)) {
301 cerr<<"Failed to rotate, alpha = "<<alpha<<endl;
302 return 0;
303 }
304 }
305 try_again--;
306 }
307 }
308
309 return 1;
310}
311
312
313
314//_____________________________________________________________________________
bbf92647 315void AliTRDtracker::GetEvent(const Char_t *hitfile, const Char_t *clusterfile)
46d29e70 316{
317 // Opens a ROOT-file with TRD-clusters and reads in the cluster-tree
46d29e70 318
bbf92647 319 ReadClusters(fClusters, clusterfile);
320
321 // get geometry from the file with hits
322
323 TFile *fInputFile = (TFile*) gROOT->GetListOfFiles()->FindObject(hitfile);
46d29e70 324 if (!fInputFile) {
325 printf("AliTRDtracker::Open -- ");
bbf92647 326 printf("Open the ALIROOT-file %s.\n",hitfile);
327 fInputFile = new TFile(hitfile);
46d29e70 328 }
329 else {
330 printf("AliTRDtracker::Open -- ");
bbf92647 331 printf("%s is already open.\n",hitfile);
46d29e70 332 }
333
334 // Get AliRun object from file or create it if not on file
46d29e70 335
bbf92647 336 gAlice = (AliRun*) fInputFile->Get("gAlice");
337 if (gAlice) {
338 printf("AliTRDtracker::GetEvent -- ");
339 printf("AliRun object found on file.\n");
340 }
341 else {
342 printf("AliTRDtracker::GetEvent -- ");
343 printf("Could not find AliRun object.\n");
344 }
46d29e70 345
346 // Import the Trees for the event nEvent in the file
347 Int_t nparticles = gAlice->GetEvent(fEvent);
348 cerr<<"nparticles = "<<nparticles<<endl;
349
350 if (nparticles <= 0) {
351 printf("AliTRDtracker::GetEvent -- ");
352 printf("No entries in the trees for event %d.\n",fEvent);
353 }
354
355 AliTRD *TRD = (AliTRD*) gAlice->GetDetector("TRD");
356 fGeom = TRD->GetGeometry();
357
46d29e70 358}
359
360
361//_____________________________________________________________________________
362void AliTRDtracker::SetUpSectors(AliTRDtrackingSector *sec)
363{
364 // Fills clusters into TRD tracking_sectors
365 // Note that the numbering scheme for the TRD tracking_sectors
366 // differs from that of TRD sectors
367
bbf92647 368 for (Int_t i=0; i<AliTRDgeometry::Nsect(); i++) sec[i].SetUp();
46d29e70 369
370 // Sort clusters into AliTRDtimeBin's within AliTRDtrackSector's
371
372 cerr<<"MakeSeeds: sorting clusters"<<endl;
373
374 Int_t ncl=fClusters->GetEntriesFast();
375 UInt_t index;
376 while (ncl--) {
377 printf("\r %d left",ncl);
378 AliTRDcluster *c=(AliTRDcluster*)fClusters->UncheckedAt(ncl);
379 Int_t detector=c->GetDetector(), local_time_bin=c->GetLocalTimeBin();
380 Int_t sector=fGeom->GetSector(detector);
381
382 Int_t tracking_sector=sector;
dd56b762 383 if(sector > 0) tracking_sector = AliTRDgeometry::kNsect - sector;
46d29e70 384
385 Int_t tb=sec[sector].GetTimeBin(detector,local_time_bin);
386 index=ncl;
387 sec[tracking_sector][tb].InsertCluster(c,index);
388 }
389 printf("\r\n");
390}
391
392
393//_____________________________________________________________________________
394void AliTRDtracker::MakeSeeds(Int_t inner, Int_t outer)
395{
396 // Creates track seeds using clusters in timeBins=i1,i2
397
bbf92647 398 Int_t i2 = inner, i1 = outer;
46d29e70 399
400 if (!fClusters) return;
401
dd56b762 402 AliTRDtrackingSector fTrSec[AliTRDgeometry::kNsect];
46d29e70 403 SetUpSectors(fTrSec);
404
405 // find seeds
406
407 Double_t x[5], c[15];
dd56b762 408 Int_t max_sec=AliTRDgeometry::kNsect;
46d29e70 409
0e9c2ad5 410 Double_t alpha=AliTRDgeometry::GetAlpha();
411 Double_t shift=AliTRDgeometry::GetAlpha()/2.;
46d29e70 412 Double_t cs=cos(alpha), sn=sin(alpha);
413
414 Double_t x1 =fTrSec[0].GetX(i1);
415 Double_t xx2=fTrSec[0].GetX(i2);
416
417 for (Int_t ns=0; ns<max_sec; ns++) {
418
419 Int_t nl=fTrSec[(ns-1+max_sec)%max_sec][i2];
420 Int_t nm=fTrSec[ns][i2];
421 Int_t nu=fTrSec[(ns+1)%max_sec][i2];
422
423 AliTRDtimeBin& r1=fTrSec[ns][i1];
424
425 for (Int_t is=0; is < r1; is++) {
426 Double_t y1=r1[is]->GetY(), z1=r1[is]->GetZ();
427
428 for (Int_t js=0; js < nl+nm+nu; js++) {
429 const AliTRDcluster *cl;
430 Double_t x2, y2, z2;
431 Double_t x3=0.,y3=0.;
432
433 if (js<nl) {
434 AliTRDtimeBin& r2=fTrSec[(ns-1+max_sec)%max_sec][i2];
435 cl=r2[js];
436 y2=cl->GetY(); z2=cl->GetZ();
437
438 x2= xx2*cs+y2*sn;
439 y2=-xx2*sn+y2*cs;
440
441 } else
442 if (js<nl+nm) {
443 AliTRDtimeBin& r2=fTrSec[ns][i2];
444 cl=r2[js-nl];
445 x2=xx2; y2=cl->GetY(); z2=cl->GetZ();
446 } else {
447 AliTRDtimeBin& r2=fTrSec[(ns+1)%max_sec][i2];
448 cl=r2[js-nl-nm];
449 y2=cl->GetY(); z2=cl->GetZ();
450
451 x2=xx2*cs-y2*sn;
452 y2=xx2*sn+y2*cs;
453
454 }
455
456 Double_t zz=z1 - z1/x1*(x1-x2);
457
bbf92647 458 if (TMath::Abs(zz-z2)>fMaxSeedDeltaZ) continue;
46d29e70 459
460 Double_t d=(x2-x1)*(0.-y2)-(0.-x2)*(y2-y1);
461 if (d==0.) {cerr<<"TRD MakeSeeds: Straight seed !\n"; continue;}
462
463 x[0]=y1;
464 x[1]=z1;
465 x[2]=f1trd(x1,y1,x2,y2,x3,y3);
466
bbf92647 467 if (TMath::Abs(x[2]) >= fMaxSeedC) continue;
46d29e70 468
469 x[3]=f2trd(x1,y1,x2,y2,x3,y3);
470
471 if (TMath::Abs(x[2]*x1-x[3]) >= 0.99999) continue;
472
473 x[4]=f3trd(x1,y1,x2,y2,z1,z2);
474
bbf92647 475 if (TMath::Abs(x[4]) > fMaxSeedTan) continue;
46d29e70 476
477 Double_t a=asin(x[3]);
478 Double_t zv=z1 - x[4]/x[2]*(a+asin(x[2]*x1-x[3]));
bbf92647 479 if (TMath::Abs(zv)>fMaxSeedVertexZ) continue;
46d29e70 480
bbf92647 481 Double_t sy1=r1[is]->GetSigmaY2(), sz1=r1[is]->GetSigmaZ2();
482 Double_t sy2=cl->GetSigmaY2(), sz2=cl->GetSigmaZ2();
483 Double_t sy3=fSeedErrorSY3, sy=fSeedErrorSY, sz=fSeedErrorSZ;
46d29e70 484
485 Double_t f20=(f1trd(x1,y1+sy,x2,y2,x3,y3)-x[2])/sy;
486 Double_t f22=(f1trd(x1,y1,x2,y2+sy,x3,y3)-x[2])/sy;
487 Double_t f24=(f1trd(x1,y1,x2,y2,x3,y3+sy)-x[2])/sy;
488 Double_t f30=(f2trd(x1,y1+sy,x2,y2,x3,y3)-x[3])/sy;
489 Double_t f32=(f2trd(x1,y1,x2,y2+sy,x3,y3)-x[3])/sy;
490 Double_t f34=(f2trd(x1,y1,x2,y2,x3,y3+sy)-x[3])/sy;
491 Double_t f40=(f3trd(x1,y1+sy,x2,y2,z1,z2)-x[4])/sy;
492 Double_t f41=(f3trd(x1,y1,x2,y2,z1+sz,z2)-x[4])/sz;
493 Double_t f42=(f3trd(x1,y1,x2,y2+sy,z1,z2)-x[4])/sy;
494 Double_t f43=(f3trd(x1,y1,x2,y2,z1,z2+sz)-x[4])/sz;
495
496 c[0]=sy1;
497 c[1]=0.; c[2]=sz1;
498 c[3]=f20*sy1; c[4]=0.; c[5]=f20*sy1*f20+f22*sy2*f22+f24*sy3*f24;
499 c[6]=f30*sy1; c[7]=0.; c[8]=f30*sy1*f20+f32*sy2*f22+f34*sy3*f24;
500 c[9]=f30*sy1*f30+f32*sy2*f32+f34*sy3*f34;
501 c[10]=f40*sy1; c[11]=f41*sz1; c[12]=f40*sy1*f20+f42*sy2*f22;
502 c[13]=f40*sy1*f30+f42*sy2*f32;
503 c[14]=f40*sy1*f40+f41*sz1*f41+f42*sy2*f42+f43*sz2*f43;
504
505 UInt_t index=r1.GetIndex(is);
bbf92647 506 AliTRDtrack *track=new AliTRDtrack(index, x, c, x1, ns*alpha+shift);
46d29e70 507
508 Int_t rc=FindProlongation(*track,fTrSec,ns,i2);
509
bbf92647 510 if ((rc < 0) ||
511 (track->GetNclusters() < (i1-i2)*fMinClustersInSeed)) delete track;
46d29e70 512 else {
513 fSeeds->AddLast(track); fNseeds++;
514 cerr<<"found seed "<<fNseeds<<endl;
515 }
516 }
517 }
518 }
519
46d29e70 520 fSeeds->Sort();
521}
522
523//___________________________________________________________________
524void AliTRDtracker::FindTracks()
525{
526 if (!fClusters) return;
527
dd56b762 528 AliTRDtrackingSector fTrSec[AliTRDgeometry::kNsect];
46d29e70 529 SetUpSectors(fTrSec);
530
531 // find tracks
532
533 Int_t num_of_time_bins = fTrSec[0].GetNtimeBins();
534 Int_t nseed=fSeeds->GetEntriesFast();
535
536 Int_t nSeedClusters;
537 for (Int_t i=0; i<nseed; i++) {
538 cerr<<"FindTracks: seed "<<i+1<<" out of "<<nseed<<endl;
539
bbf92647 540 AliTRDtrack& t=*((AliTRDtrack*)fSeeds->UncheckedAt(i));
46d29e70 541
542 nSeedClusters = t.GetNclusters();
bbf92647 543 Double_t alpha=t.GetAlpha();
46d29e70 544
545 if (alpha > 2.*TMath::Pi()) alpha -= 2.*TMath::Pi();
546 if (alpha < 0. ) alpha += 2.*TMath::Pi();
dd56b762 547 Int_t ns=Int_t(alpha/AliTRDgeometry::GetAlpha())%AliTRDgeometry::kNsect;
46d29e70 548
549 if (FindProlongation(t,fTrSec,ns)) {
550 cerr<<"No of clusters in the track = "<<t.GetNclusters()<<endl;
bbf92647 551 if (t.GetNclusters() >= Int_t(fMinClustersInTrack*num_of_time_bins)) {
46d29e70 552 Int_t label = GetTrackLabel(t);
553 t.SetLabel(label);
46d29e70 554 UseClusters(t);
bbf92647 555
556 AliTRDtrack *pt = new AliTRDtrack(t);
557 fTracks->AddLast(pt); fNtracks++;
558
46d29e70 559 cerr<<"found track "<<fNtracks<<endl;
560 }
561 }
bbf92647 562 delete fSeeds->RemoveAt(i);
46d29e70 563 }
564}
565
566//__________________________________________________________________
bbf92647 567void AliTRDtracker::UseClusters(AliTRDtrack t) {
46d29e70 568 Int_t ncl=t.GetNclusters();
569 for (Int_t i=0; i<ncl; i++) {
570 Int_t index = t.GetClusterIndex(i);
571 AliTRDcluster *c=(AliTRDcluster*)fClusters->UncheckedAt(index);
572 c->Use();
573 }
574}
575
576//__________________________________________________________________
bbf92647 577Int_t AliTRDtracker::GetTrackLabel(AliTRDtrack t) {
46d29e70 578
579 Int_t label=123456789, index, i, j;
580 Int_t ncl=t.GetNclusters();
bbf92647 581 const Int_t range = AliTRDgeometry::kNplan * fGeom->GetTimeMax();
46d29e70 582 Bool_t label_added;
583
d1b06c24 584 // Int_t s[range][2];
585 Int_t **s = new Int_t* [range];
586 for (i=0; i<range; i++) {
587 s[i] = new Int_t[2];
588 }
46d29e70 589 for (i=0; i<range; i++) {
590 s[i][0]=-1;
591 s[i][1]=0;
592 }
593
594 Int_t t0,t1,t2;
595 for (i=0; i<ncl; i++) {
596 index=t.GetClusterIndex(i);
597 AliTRDcluster *c=(AliTRDcluster*)fClusters->UncheckedAt(index);
598 t0=c->GetTrackIndex(0);
599 t1=c->GetTrackIndex(1);
600 t2=c->GetTrackIndex(2);
601 }
602
603 for (i=0; i<ncl; i++) {
604 index=t.GetClusterIndex(i);
605 AliTRDcluster *c=(AliTRDcluster*)fClusters->UncheckedAt(index);
606 for (Int_t k=0; k<3; k++) {
607 label=c->GetTrackIndex(k);
57527628 608 label_added=kFALSE; j=0;
46d29e70 609 if (label >= 0) {
610 while ( (!label_added) && ( j < range ) ) {
611 if (s[j][0]==label || s[j][1]==0) {
612 s[j][0]=label;
613 s[j][1]=s[j][1]+1;
57527628 614 label_added=kTRUE;
46d29e70 615 }
616 j++;
617 }
618 }
619 }
620 }
621
46d29e70 622 Int_t max=0;
623 label = -123456789;
624
625 for (i=0; i<range; i++) {
626 if (s[i][1]>max) {
627 max=s[i][1]; label=s[i][0];
628 }
629 }
d1b06c24 630 delete []s;
bbf92647 631 if(max > ncl*fLabelFraction) return label;
46d29e70 632 else return -1;
633}
634
635//___________________________________________________________________
636
bbf92647 637Int_t AliTRDtracker::WriteTracks(const Char_t *filename) {
46d29e70 638
639 TDirectory *savedir=gDirectory;
640
bbf92647 641 TFile *out=TFile::Open(filename,"RECREATE");
46d29e70 642
643 TTree tracktree("TreeT","Tree with TRD tracks");
644
645 AliTRDtrack *iotrack=0;
646 tracktree.Branch("tracks","AliTRDtrack",&iotrack,32000,0);
647
648 Int_t ntracks=fTracks->GetEntriesFast();
649
650 for (Int_t i=0; i<ntracks; i++) {
651 AliTRDtrack *pt=(AliTRDtrack*)fTracks->UncheckedAt(i);
652 iotrack=pt;
653 tracktree.Fill();
654 cerr<<"WriteTracks: put track "<<i<<" in the tree"<<endl;
655 }
656
657 tracktree.Write();
658 out->Close();
659
660 savedir->cd();
661
662 cerr<<"WriteTracks: done"<<endl;
663 return 0;
664}
665
46d29e70 666//_____________________________________________________________________________
bbf92647 667void AliTRDtracker::ReadClusters(TObjArray *array, const Char_t *filename,
d1b06c24 668Int_t option)
46d29e70 669{
670 //
671 // Reads AliTRDclusters (option >= 0) or AliTRDrecPoints (option < 0)
672 // from the file. The names of the cluster tree and branches
673 // should match the ones used in AliTRDclusterizer::WriteClusters()
674 //
675
676 TDirectory *savedir=gDirectory;
677
678 TFile *file = TFile::Open(filename);
679 if (!file->IsOpen()) {printf("Can't open file %s !\n",filename); return;}
680
bbf92647 681 TTree *tree = (TTree*)file->Get("ClusterTree");
46d29e70 682 Int_t nentr = (Int_t) tree->GetEntries();
bbf92647 683 printf("found %d entries in %s.\n",nentr,tree->GetName());
684
685 TBranch *branch;
686 TObjArray *ioArray = new TObjArray(400);
687
688 if( option < 0 ) {
a2b90f83 689 //branch = tree->GetBranch("RecPoints");
690 // changed CBL
691 branch = tree->GetBranch("TRDrecPoints");
bbf92647 692
693 for (Int_t i=0; i<nentr; i++) {
694 branch->SetAddress(&ioArray);
695 tree->GetEvent(i);
696 Int_t npoints = ioArray->GetEntriesFast();
697 printf("Read %d rec. points from entry %d \n", npoints, i);
698
699 for(Int_t j=0; j<npoints; j++) {
700 AliTRDrecPoint *p=(AliTRDrecPoint*)ioArray->UncheckedAt(j);
701 array->AddLast(p);
702 ioArray->RemoveAt(j);
703 }
704 }
705 }
706 else {
a2b90f83 707 branch = tree->GetBranch("TRDcluster");
bbf92647 708
709 for (Int_t i=0; i<nentr; i++) {
710 branch->SetAddress(&ioArray);
711 tree->GetEvent(i);
712 Int_t npoints = ioArray->GetEntriesFast();
713 printf("Read %d clusters from entry %d \n", npoints, i);
714
715 for(Int_t j=0; j<npoints; j++) {
716 AliTRDcluster *c=(AliTRDcluster*)ioArray->UncheckedAt(j);
717 array->AddLast(c);
718 ioArray->RemoveAt(j);
719 }
46d29e70 720 }
721 }
722
723 file->Close();
724 savedir->cd();
bbf92647 725
46d29e70 726}
727