1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
18 Revision 1.38 2003/10/17 12:01:16 kowal2
19 Removed compiler warning.
21 Revision 1.37 2003/07/22 15:56:14 hristov
22 Implementing ESD functionality in the NewIO (Yu.Belikov)
24 Revision 1.35.2.3 2003/07/15 09:58:03 hristov
25 Corrected back-propagation (Yu.Belikov)
27 Revision 1.35.2.2 2003/07/14 09:19:33 hristov
28 TOF included in the combined PID (Yu.Belikov)
30 Revision 1.35.2.1 2003/07/11 10:53:01 hristov
31 Inward refit for TPC and TRD in the ESD schema (T.Kuhr)
33 Revision 1.35 2003/05/23 10:08:51 hristov
34 SetLabel replaced by SetNumber (Yu.Belikov)
36 Revision 1.34 2003/05/22 13:57:48 hristov
37 First implementation of ESD classes (Yu.Belikov)
39 Revision 1.32 2003/04/10 10:36:54 hristov
40 Code for unified TPC/TRD tracking (S.Radomski)
42 Revision 1.31 2003/03/19 17:14:11 hristov
43 Load/UnloadClusters added to the base class and the derived classes changed correspondingly. Possibility to give 2 input files for ITS and TPC tracks in PropagateBack. TRD tracker uses fEventN from the base class (T.Kuhr)
45 Revision 1.30 2003/02/28 16:13:32 hristov
48 Revision 1.29 2003/02/28 15:18:16 hristov
49 Corrections suggested by J.Chudoba
51 Revision 1.28 2003/02/27 16:15:52 hristov
52 Code for inward refitting (S.Radomski)
54 Revision 1.27 2003/02/25 16:47:58 hristov
55 allow back propagation for more than 1 event (J.Chudoba)
57 Revision 1.26 2003/02/19 08:49:46 hristov
58 Track time measurement (S.Radomski)
60 Revision 1.25 2003/01/28 16:43:35 hristov
61 Additional protection: to be discussed with the Root team (M.Ivanov)
63 Revision 1.24 2002/11/19 16:13:24 hristov
64 stdlib.h included to declare exit() on HP
66 Revision 1.23 2002/11/19 11:50:08 hristov
67 Removing CONTAINERS (Yu.Belikov)
69 Revision 1.19 2002/07/19 07:31:40 kowal2
70 Improvement in tracking by J. Belikov
72 Revision 1.18 2002/05/13 07:33:52 kowal2
73 Added protection in Int_t AliTPCtracker::AliTPCRow::Find(Double_t y) const
74 in the case of defined region of interests.
76 Revision 1.17 2002/03/18 17:59:13 kowal2
77 Chnges in the pad geometry - 3 pad lengths introduced.
79 Revision 1.16 2001/11/08 16:39:03 hristov
80 Additional protection (M.Masera)
82 Revision 1.15 2001/11/08 16:36:33 hristov
83 Updated V2 stream of tracking (Yu.Belikov). The new long waited features are: 1) Possibility to pass the primary vertex position to the trackers (both for the TPC and the ITS) 2) Possibility to specify the number of tracking passes together with applying (or not applying) the vertex constraint (ITS only) 3) Possibility to make some use of partial PID provided by the TPC when doing tracking in the ITS (ITS only) 4) V0 reconstruction with a helix minimisation of the DCA. (new macros: AliV0FindVertices.C and AliV0Comparison.C) 4a) ( Consequence of the 4) ) All the efficiencies and resolutions are from now on calculated including *secondary*particles* too. (Don't be surprised by the drop in efficiency etc)
85 Revision 1.14 2001/10/21 19:04:55 hristov
86 Several patches were done to adapt the barel reconstruction to the multi-event case. Some memory leaks were corrected. (Yu.Belikov)
88 Revision 1.13 2001/07/23 12:01:30 hristov
91 Revision 1.12 2001/07/20 14:32:44 kowal2
92 Processing of many events possible now
94 Revision 1.11 2001/05/23 08:50:10 hristov
97 Revision 1.10 2001/05/16 14:57:25 alibrary
98 New files for folders and Stack
100 Revision 1.9 2001/05/11 07:16:56 hristov
101 Fix needed on Sun and Alpha
103 Revision 1.8 2001/05/08 15:00:15 hristov
104 Corrections for tracking in arbitrary magnenetic field. Changes towards a concept of global Alice track. Back propagation of reconstructed tracks (Yu.Belikov)
106 Revision 1.5 2000/12/20 07:51:59 kowal2
107 Changes suggested by Alessandra and Paolo to avoid overlapped
108 data fields in encapsulated classes.
110 Revision 1.4 2000/11/02 07:27:16 kowal2
113 Revision 1.2 2000/06/30 12:07:50 kowal2
114 Updated from the TPC-PreRelease branch
116 Revision 1.1.2.1 2000/06/25 08:53:55 kowal2
117 Splitted from AliTPCtracking
121 //-------------------------------------------------------
122 // Implementation of the TPC tracker
124 // Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch
125 //-------------------------------------------------------
126 #include <TObjArray.h>
134 #include "AliTPCtracker.h"
135 #include "AliTPCcluster.h"
136 #include "AliTPCParam.h"
137 #include "AliClusters.h"
139 extern Double_t SigmaY2(Double_t, Double_t, Double_t);
140 extern Double_t SigmaZ2(Double_t, Double_t);
142 ClassImp(AliTPCtracker)
144 //_____________________________________________________________________________
145 AliTPCtracker::AliTPCtracker(const AliTPCParam *par):
146 AliTracker(), fkNIS(par->GetNInnerSector()/2), fkNOS(par->GetNOuterSector()/2)
148 //---------------------------------------------------------------------
149 // The main TPC tracker constructor
150 //---------------------------------------------------------------------
151 fInnerSec=new AliTPCSector[fkNIS];
152 fOuterSec=new AliTPCSector[fkNOS];
155 for (i=0; i<fkNIS; i++) fInnerSec[i].Setup(par,0);
156 for (i=0; i<fkNOS; i++) fOuterSec[i].Setup(par,1);
158 fParam = (AliTPCParam*) par;
170 //_____________________________________________________________________________
171 AliTPCtracker::~AliTPCtracker() {
172 //------------------------------------------------------------------
173 // TPC tracker destructor
174 //------------------------------------------------------------------
184 fBarrelFile->Close();
189 //_____________________________________________________________________________
191 void AliTPCtracker::SetBarrelTree(const char *mode) {
193 // Creates a tree for BarrelTracks
194 // mode = "back" or "refit"
199 if (!IsStoringBarrel()) return;
201 TDirectory *sav = gDirectory;
202 if (!fBarrelFile) fBarrelFile = new TFile("AliBarrelTracks.root", "UPDATE");
205 sprintf(buff, "BarrelTPC_%d_%s", GetEventNumber(), mode);
208 fBarrelTree = new TTree(buff, "Barrel TPC tracks");
210 if (!fBarrelArray) fBarrelArray = new TClonesArray("AliBarrelTrack", 4);
211 for(Int_t i=0; i<4; i++) new((*fBarrelArray)[i]) AliBarrelTrack();
213 fBarrelTree->Branch("tracks", &fBarrelArray);
217 //_____________________________________________________________________________
219 void AliTPCtracker::StoreBarrelTrack(AliTPCtrack *ps, Int_t refPlane, Int_t isIn) {
221 // Stores Track at a given reference plane
224 // isIn: 1 - backward, 2 - refit
227 if (!IsStoringBarrel()) return;
228 if (refPlane < 0 || refPlane > 4) return;
229 if (isIn > 2) return;
231 static Int_t nClusters;
233 static Double_t chi2;
236 Int_t newClusters, newWrong;
239 if ( (refPlane == 1 && isIn == kTrackBack) ||
240 (refPlane == 4 && isIn == kTrackRefit) ) {
242 fBarrelArray->Clear();
243 nClusters = nWrong = 0;
250 if (refPlane == 1) refX = fParam->GetInnerRadiusLow();
251 if (refPlane == 2) refX = fParam->GetInnerRadiusUp();
252 if (refPlane == 3) refX = fParam->GetOuterRadiusLow();
253 if (refPlane == 4) refX = fParam->GetOuterRadiusUp();
255 ps->PropagateTo(refX);
257 fBarrelTrack = (AliBarrelTrack*)(*fBarrelArray)[index++];
258 ps->GetBarrelTrack(fBarrelTrack);
260 newClusters = ps->GetNumberOfClusters() - nClusters;
261 newWrong = ps->GetNWrong() - nWrong;
262 newChi2 = ps->GetChi2() - chi2;
264 nClusters = ps->GetNumberOfClusters();
265 nWrong = ps->GetNWrong();
266 chi2 = ps->GetChi2();
268 fBarrelTrack->SetNClusters(newClusters, newChi2);
269 fBarrelTrack->SetNWrongClusters(newWrong);
270 fBarrelTrack->SetRefPlane(refPlane, isIn);
274 //_____________________________________________________________________________
275 Double_t f1(Double_t x1,Double_t y1,
276 Double_t x2,Double_t y2,
277 Double_t x3,Double_t y3)
279 //-----------------------------------------------------------------
280 // Initial approximation of the track curvature
281 //-----------------------------------------------------------------
282 Double_t d=(x2-x1)*(y3-y2)-(x3-x2)*(y2-y1);
283 Double_t a=0.5*((y3-y2)*(y2*y2-y1*y1+x2*x2-x1*x1)-
284 (y2-y1)*(y3*y3-y2*y2+x3*x3-x2*x2));
285 Double_t b=0.5*((x2-x1)*(y3*y3-y2*y2+x3*x3-x2*x2)-
286 (x3-x2)*(y2*y2-y1*y1+x2*x2-x1*x1));
288 Double_t xr=TMath::Abs(d/(d*x1-a)), yr=d/(d*y1-b);
290 return -xr*yr/sqrt(xr*xr+yr*yr);
294 //_____________________________________________________________________________
295 Double_t f2(Double_t x1,Double_t y1,
296 Double_t x2,Double_t y2,
297 Double_t x3,Double_t y3)
299 //-----------------------------------------------------------------
300 // Initial approximation of the track curvature times center of curvature
301 //-----------------------------------------------------------------
302 Double_t d=(x2-x1)*(y3-y2)-(x3-x2)*(y2-y1);
303 Double_t a=0.5*((y3-y2)*(y2*y2-y1*y1+x2*x2-x1*x1)-
304 (y2-y1)*(y3*y3-y2*y2+x3*x3-x2*x2));
305 Double_t b=0.5*((x2-x1)*(y3*y3-y2*y2+x3*x3-x2*x2)-
306 (x3-x2)*(y2*y2-y1*y1+x2*x2-x1*x1));
308 Double_t xr=TMath::Abs(d/(d*x1-a)), yr=d/(d*y1-b);
310 return -a/(d*y1-b)*xr/sqrt(xr*xr+yr*yr);
313 //_____________________________________________________________________________
314 Double_t f3(Double_t x1,Double_t y1,
315 Double_t x2,Double_t y2,
316 Double_t z1,Double_t z2)
318 //-----------------------------------------------------------------
319 // Initial approximation of the tangent of the track dip angle
320 //-----------------------------------------------------------------
321 return (z1 - z2)/sqrt((x1-x2)*(x1-x2)+(y1-y2)*(y1-y2));
324 //_____________________________________________________________________________
325 Int_t AliTPCtracker::LoadClusters(TTree *cTree) {
326 //-----------------------------------------------------------------
327 // This function loads TPC clusters.
328 //-----------------------------------------------------------------
329 TBranch *branch=cTree->GetBranch("Segment");
331 Error("LoadClusters","Can't get the branch !");
335 AliClusters carray, *addr=&carray;
336 carray.SetClass("AliTPCcluster");
338 branch->SetAddress(&addr);
340 Int_t nentr=(Int_t)cTree->GetEntries();
342 for (Int_t i=0; i<nentr; i++) {
345 Int_t ncl=carray.GetArray()->GetEntriesFast();
347 Int_t nir=fInnerSec->GetNRows(), nor=fOuterSec->GetNRows();
348 Int_t id=carray.GetID();
349 if ((id<0) || (id>2*(fkNIS*nir + fkNOS*nor))) {
350 Error("LoadClusters","Wrong index !");
353 Int_t outindex = 2*fkNIS*nir;
356 Int_t row = id - sec*nir;
358 AliTPCRow &padrow=fInnerSec[sec][row];
360 AliTPCcluster *c=(AliTPCcluster*)carray[ncl];
361 padrow.InsertCluster(c,sec,row);
366 Int_t row = id - sec*nor;
368 AliTPCRow &padrow=fOuterSec[sec][row];
370 AliTPCcluster *c=(AliTPCcluster*)carray[ncl];
371 padrow.InsertCluster(c,sec+fkNIS,row);
374 carray.GetArray()->Clear();
380 //_____________________________________________________________________________
381 void AliTPCtracker::UnloadClusters() {
382 //-----------------------------------------------------------------
383 // This function unloads TPC clusters.
384 //-----------------------------------------------------------------
386 for (i=0; i<fkNIS; i++) {
387 Int_t nr=fInnerSec->GetNRows();
388 for (Int_t n=0; n<nr; n++) fInnerSec[i][n].ResetClusters();
390 for (i=0; i<fkNOS; i++) {
391 Int_t nr=fOuterSec->GetNRows();
392 for (Int_t n=0; n<nr; n++) fOuterSec[i][n].ResetClusters();
396 //_____________________________________________________________________________
397 Int_t AliTPCtracker::FollowProlongation(AliTPCseed& t, Int_t rf) {
398 //-----------------------------------------------------------------
399 // This function tries to find a track prolongation.
400 //-----------------------------------------------------------------
401 Double_t xt=t.GetX();
402 const Int_t kSKIP=(t.GetNumberOfClusters()<10) ? kRowsToSkip :
403 Int_t(0.5*fSectors->GetNRows());
404 Int_t tryAgain=kSKIP;
406 Double_t alpha=t.GetAlpha() - fSectors->GetAlphaShift();
407 if (alpha > 2.*TMath::Pi()) alpha -= 2.*TMath::Pi();
408 if (alpha < 0. ) alpha += 2.*TMath::Pi();
409 Int_t s=Int_t(alpha/fSectors->GetAlpha())%fN;
411 Int_t nrows=fSectors->GetRowNumber(xt)-1;
412 for (Int_t nr=nrows; nr>=rf; nr--) {
413 Double_t x=fSectors->GetX(nr), ymax=fSectors->GetMaxY(nr);
414 if (!t.PropagateTo(x)) return 0;
418 Double_t maxchi2=kMaxCHI2;
419 const AliTPCRow &krow=fSectors[s][nr];
420 Double_t pt=t.GetConvConst()/(100/0.299792458/0.2)/t.Get1Pt();
421 Double_t sy2=SigmaY2(t.GetX(),t.GetTgl(),pt);
422 Double_t sz2=SigmaZ2(t.GetX(),t.GetTgl());
423 Double_t road=4.*sqrt(t.GetSigmaY2() + sy2), y=t.GetY(), z=t.GetZ();
426 if (t.GetNumberOfClusters()>4)
427 Warning("FindProlongation","Too broad road !");
432 for (Int_t i=krow.Find(y-road); i<krow; i++) {
433 AliTPCcluster *c=(AliTPCcluster*)(krow[i]);
434 if (c->GetY() > y+road) break;
435 if (c->IsUsed()) continue;
436 if ((c->GetZ()-z)*(c->GetZ()-z) > 16.*(t.GetSigmaZ2()+sz2)) continue;
437 Double_t chi2=t.GetPredictedChi2(c);
438 if (chi2 > maxchi2) continue;
441 index=krow.GetIndex(i);
445 Float_t l=fSectors->GetPadPitchWidth();
446 Float_t corr=1.; if (nr>63) corr=0.67; // new (third) pad response !
447 t.SetSampledEdx(cl->GetQ()/l*corr,t.GetNumberOfClusters());
448 if (!t.Update(cl,maxchi2,index)) {
449 if (!tryAgain--) return 0;
450 } else tryAgain=kSKIP;
452 if (tryAgain==0) break;
455 if (!t.Rotate(fSectors->GetAlpha())) return 0;
456 } else if (y <-ymax) {
458 if (!t.Rotate(-fSectors->GetAlpha())) return 0;
466 //_____________________________________________________________________________
468 Int_t AliTPCtracker::FollowRefitInward(AliTPCseed *seed, AliTPCtrack *track) {
470 // This function propagates seed inward TPC using old clusters
473 // Sylwester Radomski, GSI
479 Int_t nRows = fSectors->GetNRows();
480 for (Int_t iRow = nRows; iRow >= 0; iRow--) {
482 Double_t x = fSectors->GetX(iRow);
483 if (!seed->PropagateTo(x)) return 0;
485 // try to find an assigned cluster in this row
487 AliTPCcluster* cluster = NULL;
490 for (Int_t iCluster = 0; iCluster < track->GetNumberOfClusters(); iCluster++) {
491 idx = track->GetClusterIndex(iCluster);
492 sec = (idx&0xff000000)>>24;
493 Int_t row = (idx&0x00ff0000)>>16;
494 if (((fSectors == fInnerSec) && (sec >= fkNIS)) ||
495 ((fSectors == fOuterSec) && (sec < fkNIS))) continue;
497 cluster = (AliTPCcluster*) GetCluster(idx);
502 // update the track seed with the found cluster
505 Double_t dAlpha = fParam->GetAngle(sec) - seed->GetAlpha();
506 if (TMath::Abs(dAlpha) > 0.0001) {
507 if (!seed->Rotate(dAlpha)) return 0;
508 if (!seed->PropagateTo(x)) return 0;
511 seed->Update(cluster, seed->GetPredictedChi2(cluster), idx);
518 //_____________________________________________________________________________
519 Int_t AliTPCtracker::FollowBackProlongation
520 (AliTPCseed& seed, const AliTPCtrack &track) {
521 //-----------------------------------------------------------------
522 // This function propagates tracks back through the TPC
523 //-----------------------------------------------------------------
524 Double_t alpha=seed.GetAlpha() - fSectors->GetAlphaShift();
525 if (alpha > 2.*TMath::Pi()) alpha -= 2.*TMath::Pi();
526 if (alpha < 0. ) alpha += 2.*TMath::Pi();
527 Int_t s=Int_t(alpha/fSectors->GetAlpha())%fN;
529 Int_t idx=-1, sec=-1, row=-1;
530 Int_t nc=seed.GetNumber();
533 idx=track.GetClusterIndex(nc);
534 sec=(idx&0xff000000)>>24; row=(idx&0x00ff0000)>>16;
536 if (fSectors==fInnerSec) { if (sec >= fkNIS) row=-1; }
537 else { if (sec < fkNIS) row=-1; }
539 Int_t nr=fSectors->GetNRows();
540 for (Int_t i=0; i<nr; i++) {
541 Double_t x=fSectors->GetX(i), ymax=fSectors->GetMaxY(i);
542 Double_t y=seed.GetYat(x);
546 if (!seed.Rotate(fSectors->GetAlpha())) return 0;
547 } else if (y <-ymax) {
549 if (!seed.Rotate(-fSectors->GetAlpha())) return 0;
552 if (!seed.PropagateTo(x)) return 0;
556 Double_t maxchi2=kMaxCHI2;
557 Double_t pt=seed.GetConvConst()/(100/0.299792458/0.2)/seed.Get1Pt();
558 Double_t sy2=SigmaY2(seed.GetX(),seed.GetTgl(),pt);
559 Double_t sz2=SigmaZ2(seed.GetX(),seed.GetTgl());
560 Double_t road=4.*sqrt(seed.GetSigmaY2() + sy2), z=seed.GetZ();
562 Warning("FollowBackProlongation","Too broad road !");
566 Int_t accepted=seed.GetNumberOfClusters();
568 //try to accept already found cluster
569 AliTPCcluster *c=(AliTPCcluster*)GetCluster(idx);
571 if ((chi2=seed.GetPredictedChi2(c))<maxchi2 || accepted<27) {
572 index=idx; cl=c; maxchi2=chi2;
573 } //else cerr<<"AliTPCtracker::FollowBackProlongation: oulier !\n";
576 idx=track.GetClusterIndex(nc);
577 sec=(idx&0xff000000)>>24; row=(idx&0x00ff0000)>>16;
579 if (fSectors==fInnerSec) { if (sec >= fkNIS) row=-1; }
580 else { if (sec < fkNIS) row=-1; }
584 //try to fill the gap
585 const AliTPCRow &krow=fSectors[s][i];
588 for (Int_t i=krow.Find(y-road); i<krow; i++) {
589 AliTPCcluster *c=(AliTPCcluster*)(krow[i]);
590 if (c->GetY() > y+road) break;
591 if (c->IsUsed()) continue;
592 if ((c->GetZ()-z)*(c->GetZ()-z)>16.*(seed.GetSigmaZ2()+sz2)) continue;
593 Double_t chi2=seed.GetPredictedChi2(c);
594 if (chi2 > maxchi2) continue;
597 index=krow.GetIndex(i);
603 Float_t l=fSectors->GetPadPitchWidth();
604 Float_t corr=1.; if (i>63) corr=0.67; // new (third) pad response !
605 seed.SetSampledEdx(cl->GetQ()/l*corr,seed.GetNumberOfClusters());
606 seed.Update(cl,maxchi2,index);
616 //_____________________________________________________________________________
617 void AliTPCtracker::MakeSeeds(Int_t i1, Int_t i2) {
618 //-----------------------------------------------------------------
619 // This function creates track seeds.
620 //-----------------------------------------------------------------
621 Double_t x[5], c[15];
623 Double_t alpha=fSectors->GetAlpha(), shift=fSectors->GetAlphaShift();
624 Double_t cs=cos(alpha), sn=sin(alpha);
626 Double_t x1 =fSectors->GetX(i1);
627 Double_t xx2=fSectors->GetX(i2);
629 for (Int_t ns=0; ns<fN; ns++) {
630 Int_t nl=fSectors[(ns-1+fN)%fN][i2];
631 Int_t nm=fSectors[ns][i2];
632 Int_t nu=fSectors[(ns+1)%fN][i2];
633 const AliTPCRow& kr1=fSectors[ns][i1];
634 for (Int_t is=0; is < kr1; is++) {
635 Double_t y1=kr1[is]->GetY(), z1=kr1[is]->GetZ();
636 for (Int_t js=0; js < nl+nm+nu; js++) {
637 const AliTPCcluster *kcl;
639 Double_t x3=GetX(), y3=GetY(), z3=GetZ();
642 const AliTPCRow& kr2=fSectors[(ns-1+fN)%fN][i2];
644 y2=kcl->GetY(); z2=kcl->GetZ();
649 const AliTPCRow& kr2=fSectors[ns][i2];
651 x2=xx2; y2=kcl->GetY(); z2=kcl->GetZ();
653 const AliTPCRow& kr2=fSectors[(ns+1)%fN][i2];
655 y2=kcl->GetY(); z2=kcl->GetZ();
660 Double_t zz=z1 - (z1-z3)/(x1-x3)*(x1-x2);
661 if (TMath::Abs(zz-z2)>5.) continue;
663 Double_t d=(x2-x1)*(0.-y2)-(0.-x2)*(y2-y1);
665 Warning("MakeSeeds","Straight seed !");
670 x[4]=f1(x1,y1,x2,y2,x3,y3);
671 if (TMath::Abs(x[4]) >= 0.0066) continue;
672 x[2]=f2(x1,y1,x2,y2,x3,y3);
673 //if (TMath::Abs(x[4]*x1-x[2]) >= 0.99999) continue;
674 x[3]=f3(x1,y1,x2,y2,z1,z2);
675 if (TMath::Abs(x[3]) > 1.2) continue;
676 Double_t a=asin(x[2]);
677 Double_t zv=z1 - x[3]/x[4]*(a+asin(x[4]*x1-x[2]));
678 if (TMath::Abs(zv-z3)>10.) continue;
680 Double_t sy1=kr1[is]->GetSigmaY2(), sz1=kr1[is]->GetSigmaZ2();
681 Double_t sy2=kcl->GetSigmaY2(), sz2=kcl->GetSigmaZ2();
682 //Double_t sy3=400*3./12., sy=0.1, sz=0.1;
683 Double_t sy3=25000*x[4]*x[4]+0.1, sy=0.1, sz=0.1;
685 Double_t f40=(f1(x1,y1+sy,x2,y2,x3,y3)-x[4])/sy;
686 Double_t f42=(f1(x1,y1,x2,y2+sy,x3,y3)-x[4])/sy;
687 Double_t f43=(f1(x1,y1,x2,y2,x3,y3+sy)-x[4])/sy;
688 Double_t f20=(f2(x1,y1+sy,x2,y2,x3,y3)-x[2])/sy;
689 Double_t f22=(f2(x1,y1,x2,y2+sy,x3,y3)-x[2])/sy;
690 Double_t f23=(f2(x1,y1,x2,y2,x3,y3+sy)-x[2])/sy;
691 Double_t f30=(f3(x1,y1+sy,x2,y2,z1,z2)-x[3])/sy;
692 Double_t f31=(f3(x1,y1,x2,y2,z1+sz,z2)-x[3])/sz;
693 Double_t f32=(f3(x1,y1,x2,y2+sy,z1,z2)-x[3])/sy;
694 Double_t f34=(f3(x1,y1,x2,y2,z1,z2+sz)-x[3])/sz;
698 c[3]=f20*sy1; c[4]=0.; c[5]=f20*sy1*f20+f22*sy2*f22+f23*sy3*f23;
699 c[6]=f30*sy1; c[7]=f31*sz1; c[8]=f30*sy1*f20+f32*sy2*f22;
700 c[9]=f30*sy1*f30+f31*sz1*f31+f32*sy2*f32+f34*sz2*f34;
701 c[10]=f40*sy1; c[11]=0.; c[12]=f40*sy1*f20+f42*sy2*f22+f43*sy3*f23;
702 c[13]=f30*sy1*f40+f32*sy2*f42;
703 c[14]=f40*sy1*f40+f42*sy2*f42+f43*sy3*f43;
705 UInt_t index=kr1.GetIndex(is);
706 AliTPCseed *track=new AliTPCseed(index, x, c, x1, ns*alpha+shift);
707 Float_t l=fSectors->GetPadPitchWidth();
708 track->SetSampledEdx(kr1[is]->GetQ()/l,0);
710 Int_t rc=FollowProlongation(*track, i2);
711 if (rc==0 || track->GetNumberOfClusters()<(i1-i2)/2) delete track;
712 else fSeeds->AddLast(track);
718 //_____________________________________________________________________________
719 Int_t AliTPCtracker::ReadSeeds(const TFile *inp) {
720 //-----------------------------------------------------------------
721 // This function reades track seeds.
722 //-----------------------------------------------------------------
723 TDirectory *savedir=gDirectory;
725 TFile *in=(TFile*)inp;
727 Error("ReadSeeds","Input file has not been open !");
732 TTree *seedTree=(TTree*)in->Get("Seeds");
734 Error("ReadSeeds","Can't get a tree with track seeds !");
737 AliTPCtrack *seed=new AliTPCtrack;
738 seedTree->SetBranchAddress("tracks",&seed);
740 if (fSeeds==0) fSeeds=new TObjArray(15000);
742 Int_t n=(Int_t)seedTree->GetEntries();
743 for (Int_t i=0; i<n; i++) {
744 seedTree->GetEvent(i);
745 fSeeds->AddLast(new AliTPCseed(*seed,seed->GetAlpha()));
750 delete seedTree; //Thanks to Mariana Bondila
757 //_____________________________________________________________________________
758 Int_t AliTPCtracker::Clusters2Tracks(AliESD *event) {
759 //-----------------------------------------------------------------
760 // This is a track finder.
761 // The clusters must be already loaded !
762 //-----------------------------------------------------------------
765 Int_t nup=fOuterSec->GetNRows(), nlow=fInnerSec->GetNRows();
766 Int_t nrows=nlow+nup;
768 Int_t gap=Int_t(0.125*nrows), shift=Int_t(0.5*gap);
769 fSectors=fOuterSec; fN=fkNOS;
770 fSeeds=new TObjArray(15000);
771 MakeSeeds(nup-1, nup-1-gap);
772 MakeSeeds(nup-1-shift, nup-1-shift-gap);
776 Int_t nseed=fSeeds->GetEntriesFast();
777 for (Int_t i=0; i<nseed; i++) {
778 //tracking in the outer sectors
779 fSectors=fOuterSec; fN=fkNOS;
781 AliTPCseed *pt=(AliTPCseed*)fSeeds->UncheckedAt(i), &t=*pt;
782 if (!FollowProlongation(t)) {
783 delete fSeeds->RemoveAt(i);
787 //tracking in the inner sectors
788 fSectors=fInnerSec; fN=fkNIS;
790 Double_t alpha=t.GetAlpha() - fInnerSec->GetAlphaShift();
791 if (alpha > 2.*TMath::Pi()) alpha -= 2.*TMath::Pi();
792 if (alpha < 0. ) alpha += 2.*TMath::Pi();
793 Int_t ns=Int_t(alpha/fInnerSec->GetAlpha())%fkNIS;
795 alpha=ns*fInnerSec->GetAlpha()+fInnerSec->GetAlphaShift()-t.GetAlpha();
797 if (t.Rotate(alpha)) {
798 if (FollowProlongation(t)) {
799 if (t.GetNumberOfClusters() >= Int_t(0.4*nrows)) {
801 CookLabel(pt,0.1); //For comparison only
802 pt->PropagateTo(fParam->GetInnerRadiusLow());
804 iotrack.UpdateTrackParams(pt,AliESDtrack::kTPCin);
806 event->AddTrack(&iotrack);
812 delete fSeeds->RemoveAt(i);
815 Info("Clusters2Tracks","Number of found tracks : %d",
816 event->GetNumberOfTracks());
818 fSeeds->Clear(); delete fSeeds; fSeeds=0;
823 //_____________________________________________________________________________
824 Int_t AliTPCtracker::Clusters2Tracks(TTree *cTree, TTree *tTree) {
825 //-----------------------------------------------------------------
826 // This is a track finder.
827 //-----------------------------------------------------------------
828 if (LoadClusters(cTree)) {
829 Error("Clusters2Tracks","Problem with loading the clusters...");
833 AliTPCtrack *iotrack=0;
834 TBranch *branch=tTree->GetBranch("tracks");
835 if (!branch) tTree->Branch("tracks","AliTPCtrack",&iotrack,32000,3);
836 else branch->SetAddress(&iotrack);
839 Int_t nup=fOuterSec->GetNRows(), nlow=fInnerSec->GetNRows();
840 Int_t nrows=nlow+nup;
842 Int_t gap=Int_t(0.125*nrows), shift=Int_t(0.5*gap);
843 fSectors=fOuterSec; fN=fkNOS;
844 fSeeds=new TObjArray(15000);
845 MakeSeeds(nup-1, nup-1-gap);
846 MakeSeeds(nup-1-shift, nup-1-shift-gap);
851 Int_t nseed=fSeeds->GetEntriesFast();
853 for (Int_t i=0; i<nseed; i++) {
854 //tracking in the outer sectors
855 fSectors=fOuterSec; fN=fkNOS;
857 AliTPCseed *pt=(AliTPCseed*)fSeeds->UncheckedAt(i), &t=*pt;
858 if (!FollowProlongation(t)) {
859 delete fSeeds->RemoveAt(i);
863 //tracking in the inner sectors
864 fSectors=fInnerSec; fN=fkNIS;
866 Double_t alpha=t.GetAlpha() - fInnerSec->GetAlphaShift();
867 if (alpha > 2.*TMath::Pi()) alpha -= 2.*TMath::Pi();
868 if (alpha < 0. ) alpha += 2.*TMath::Pi();
869 Int_t ns=Int_t(alpha/fInnerSec->GetAlpha())%fkNIS;
871 alpha=ns*fInnerSec->GetAlpha()+fInnerSec->GetAlphaShift()-t.GetAlpha();
873 if (t.Rotate(alpha)) {
874 if (FollowProlongation(t)) {
875 if (t.GetNumberOfClusters() >= Int_t(0.4*nrows)) {
877 CookLabel(pt,0.1); //For comparison only
878 pt->PropagateTo(fParam->GetInnerRadiusLow());
886 delete fSeeds->RemoveAt(i);
889 Info("Clusters2Tracks","Number of found tracks : %d",found);
893 fSeeds->Clear(); delete fSeeds; fSeeds=0;
898 //_____________________________________________________________________________
899 Int_t AliTPCtracker::RefitInward(AliESD* event) {
901 // The function propagates tracks throught TPC inward
902 // using already associated clusters.
903 // The clusters must be already loaded !
906 Int_t nTracks = event->GetNumberOfTracks();
909 for (Int_t i = 0; i < nTracks; i++) {
910 AliESDtrack* track = event->GetTrack(i);
911 ULong_t status = track->GetStatus();
913 if ( (status & AliESDtrack::kTPCout ) == 0 ) continue;
914 if ( (status & AliESDtrack::kTPCrefit) != 0 ) continue;
916 AliTPCtrack* tpcTrack = new AliTPCtrack(*track);
917 AliTPCseed* seed = new AliTPCseed(*tpcTrack, tpcTrack->GetAlpha());
919 fSectors = fOuterSec;
921 Int_t res = FollowRefitInward(seed, tpcTrack);
923 Int_t nc = seed->GetNumberOfClusters();
925 fSectors = fInnerSec;
927 res = FollowRefitInward(seed, tpcTrack);
928 UseClusters(seed, nc);
931 seed->PropagateTo(fParam->GetInnerRadiusLow());
932 seed->SetLabel(tpcTrack->GetLabel());
933 seed->SetdEdx(tpcTrack->GetdEdx());
934 track->UpdateTrackParams(seed, AliESDtrack::kTPCrefit);
942 Info("RefitInward","Number of refitted tracks : %d",nRefited);
947 //_____________________________________________________________________________
948 Int_t AliTPCtracker::RefitInward(TTree */*in*/, TTree */*out*/) {
950 // The function propagates tracks throught TPC inward
951 // using already associated clusters.
953 Error("RefitInward","This method is not converted to NewIO yet\n");
956 if (!inSeeds->IsOpen()) {
957 cout << "Input File with seeds not open !\n" << endl;
961 if (!inTracks->IsOpen()) {
962 cout << "Input File not open !\n" << endl;
966 if (!outTracks->IsOpen()) {
967 cout << "Output File not open !\n" << endl;
971 TDirectory *savedir = gDirectory;
974 // Connect to input seeds tree
977 sprintf(tname, "seedTRDtoTPC_%d", GetEventNumber());
978 TTree *seedsTree = (TTree*)inSeeds->Get(tname);
979 TBranch *inSeedBranch = seedsTree->GetBranch("tracks");
980 AliTPCtrack *inSeedTrack = 0;
981 inSeedBranch->SetAddress(&inSeedTrack);
983 Int_t nSeeds = (Int_t)seedsTree->GetEntries();
985 // Connect to input tree
987 sprintf(tname,"TreeT_TPC_%d",GetEventNumber());
988 // sprintf(tname,"seedsTPCtoTRD_%d",GetEventNumber());
989 TTree *inputTree = (TTree*)inTracks->Get(tname);
990 TBranch *inBranch = inputTree->GetBranch("tracks");
991 AliTPCtrack *inTrack = 0;
992 inBranch->SetAddress(&inTrack);
994 Int_t nTracks = (Int_t)inputTree->GetEntries();
996 // Create output tree
998 AliTPCtrack *outTrack = new AliTPCtrack();
999 sprintf(tname, "tracksTPC_%d", GetEventNumber());
1000 TTree *outputTree = new TTree(tname,"Refited TPC tracks");
1001 outputTree->Branch("tracks", "AliTPCtrack", &outTrack, 32000, 3);
1005 // create table of track labels
1006 Int_t* inLab = new Int_t[nTracks];
1007 for (Int_t i = 0; i < nTracks; i++) inLab[i] = -1;
1009 // [SR, 01.04.2003] - Barrel tracks
1010 if (IsStoringBarrel()) SetBarrelTree("refit");
1012 for(Int_t t=0; t < nSeeds; t++) {
1014 seedsTree->GetEntry(t);
1015 // find TPC track for seed
1016 Int_t lab = TMath::Abs(inSeedTrack->GetLabel());
1017 for(Int_t i=0; i < nTracks; i++) {
1018 if (inLab[i] == lab) {
1019 inputTree->GetEntry(i);
1021 } else if (inLab[i] < 0) {
1022 inputTree->GetEntry(i);
1023 inLab[i] = TMath::Abs(inTrack->GetLabel());
1024 if (inLab[i] == lab) break;
1027 if (TMath::Abs(inTrack->GetLabel()) != lab) continue;
1028 AliTPCseed *seed = new AliTPCseed(*inSeedTrack, inTrack->GetAlpha());
1029 if (IsStoringBarrel()) StoreBarrelTrack(seed, 4, 2);
1031 fSectors = fOuterSec;
1033 Int_t res = FollowRefitInward(seed, inTrack);
1035 Int_t nc = seed->GetNumberOfClusters();
1037 if (IsStoringBarrel()) StoreBarrelTrack(seed, 3, 2);
1038 if (IsStoringBarrel()) StoreBarrelTrack(seed, 2, 2);
1040 fSectors = fInnerSec;
1042 res = FollowRefitInward(seed, inTrack);
1043 UseClusters(seed, nc);
1045 if (IsStoringBarrel()) StoreBarrelTrack(seed, 1, 2);
1048 seed->PropagateTo(fParam->GetInnerRadiusLow());
1049 outTrack = (AliTPCtrack*)seed;
1050 outTrack->SetLabel(inTrack->GetLabel());
1051 outTrack->SetdEdx(inTrack->GetdEdx());
1056 if (IsStoringBarrel()) fBarrelTree->Fill();
1060 cout << "Refitted " << nRefited << " tracks." << endl;
1063 outputTree->Write();
1068 if (IsStoringBarrel()) {
1070 fBarrelTree->Write();
1071 fBarrelFile->Flush();
1073 delete fBarrelArray;
1077 if (inputTree) delete inputTree;
1078 if (outputTree) delete outputTree;
1087 Int_t AliTPCtracker::PropagateBack(AliESD *event) {
1088 //-----------------------------------------------------------------
1089 // This function propagates tracks back through the TPC.
1090 // The clusters must be already loaded !
1091 //-----------------------------------------------------------------
1092 Int_t nentr=event->GetNumberOfTracks();
1093 Info("PropagateBack", "Number of ESD tracks: %d\n", nentr);
1096 for (Int_t i=0; i<nentr; i++) {
1097 AliESDtrack *esd=event->GetTrack(i);
1098 ULong_t status=esd->GetStatus();
1100 if ( (status & AliESDtrack::kTPCin ) == 0 ) continue;
1101 if ( (status & AliESDtrack::kTPCout) != 0 ) continue;
1103 const AliTPCtrack t(*esd);
1104 AliTPCseed s(t,t.GetAlpha());
1106 if ( (status & AliESDtrack::kITSout) == 0 ) s.ResetCovariance();
1111 Int_t nc=t.GetNumberOfClusters();
1112 s.SetNumber(nc); //set number of the cluster to start with
1115 fSectors=fInnerSec; fN=fkNIS;
1117 Double_t alpha=s.GetAlpha() - fSectors->GetAlphaShift();
1118 if (alpha > 2.*TMath::Pi()) alpha -= 2.*TMath::Pi();
1119 if (alpha < 0. ) alpha += 2.*TMath::Pi();
1120 Int_t ns=Int_t(alpha/fSectors->GetAlpha())%fN;
1121 alpha =ns*fSectors->GetAlpha() + fSectors->GetAlphaShift();
1122 alpha-=s.GetAlpha();
1124 if (!s.Rotate(alpha)) continue;
1125 if (!FollowBackProlongation(s,t)) continue;
1130 fSectors=fOuterSec; fN=fkNOS;
1132 nc=s.GetNumberOfClusters();
1134 alpha=s.GetAlpha() - fSectors->GetAlphaShift();
1135 if (alpha > 2.*TMath::Pi()) alpha -= 2.*TMath::Pi();
1136 if (alpha < 0. ) alpha += 2.*TMath::Pi();
1137 ns=Int_t(alpha/fSectors->GetAlpha())%fN;
1139 alpha =ns*fSectors->GetAlpha() + fSectors->GetAlphaShift();
1140 alpha-=s.GetAlpha();
1142 if (!s.Rotate(alpha)) continue;
1143 if (!FollowBackProlongation(s,t)) continue;
1145 Int_t nrows=fOuterSec->GetNRows()+fInnerSec->GetNRows();
1146 if (s.GetNumberOfClusters() < Int_t(0.4*nrows)) continue;
1148 s.PropagateTo(fParam->GetOuterRadiusUp());
1150 CookLabel(&s,0.1); //For comparison only
1152 esd->UpdateTrackParams(&s,AliESDtrack::kTPCout);
1155 Info("PropagateBack","Number of back propagated tracks: %d",ntrk);
1160 //_____________________________________________________________________________
1161 Int_t AliTPCtracker::PropagateBack(TTree */*in*/, TTree */*out*/) {
1162 //-----------------------------------------------------------------
1163 // This function propagates tracks back through the TPC.
1164 //-----------------------------------------------------------------
1165 Error("PropagateBack","This method is not converted to NewIO yet\n");
1168 fSeeds=new TObjArray(15000);
1170 TFile *in=(TFile*)inp;
1171 TFile *in2=(TFile*)inp2;
1172 TDirectory *savedir=gDirectory;
1174 if (!in->IsOpen()) {
1175 cerr<<"AliTPCtracker::PropagateBack(): ";
1176 cerr<<"file with TPC (or back propagated ITS) tracks is not open !\n";
1180 if (!out->IsOpen()) {
1181 cerr<<"AliTPCtracker::PropagateBack(): ";
1182 cerr<<"file for back propagated TPC tracks is not open !\n";
1190 sprintf(tName,"TreeT_ITSb_%d",GetEventNumber());
1191 TTree *bckTree=(TTree*)in->Get(tName);
1192 if (!bckTree && inp2) bckTree=(TTree*)in2->Get(tName);
1194 cerr<<"AliTPCtracker::PropagateBack() ";
1195 cerr<<"can't get a tree with back propagated ITS tracks !\n";
1198 AliTPCtrack *bckTrack=new AliTPCtrack;
1199 if (bckTree) bckTree->SetBranchAddress("tracks",&bckTrack);
1201 sprintf(tName,"TreeT_TPC_%d",GetEventNumber());
1202 TTree *tpcTree=(TTree*)in->Get(tName);
1204 cerr<<"AliTPCtracker::PropagateBack() ";
1205 cerr<<"can't get a tree with TPC tracks !\n";
1208 AliTPCtrack *tpcTrack=new AliTPCtrack;
1209 tpcTree->SetBranchAddress("tracks",&tpcTrack);
1211 // [SR, 01.04.2003] - Barrel tracks
1212 if (IsStoringBarrel()) SetBarrelTree("back");
1214 // Prepare an array of tracks to be back propagated
1215 Int_t nup=fOuterSec->GetNRows(), nlow=fInnerSec->GetNRows();
1216 Int_t nrows=nlow+nup;
1219 // Match ITS tracks with old TPC tracks
1220 // the tracks do not have to be sorted [SR, GSI, 18.02.2003]
1222 // the algorithm is linear and uses LUT for sorting
1223 // cost of the algorithm: 2 * number of tracks
1226 TObjArray tracks(15000);
1228 Int_t tpcN= (Int_t)tpcTree->GetEntries();
1229 Int_t bckN= (bckTree)? (Int_t)bckTree->GetEntries() : 0;
1232 const Int_t nLab = 200000; // limit on number of primaries (arbitrary)
1234 for(Int_t i=0; i<nLab; i++) lut[i] = -1;
1237 for(Int_t i=0; i<bckN; i++) {
1238 bckTree->GetEvent(i);
1239 Int_t lab = TMath::Abs(bckTrack->GetLabel());
1240 if (lab < nLab) lut[lab] = i;
1242 cerr << "AliTPCtracker: The size of the LUT is too small\n";
1247 for (Int_t i=0; i<tpcN; i++) {
1248 tpcTree->GetEvent(i);
1249 Double_t alpha=tpcTrack->GetAlpha();
1253 // No ITS - use TPC track only
1255 tpcTrack->ResetCovariance();
1256 AliTPCseed * seed = new AliTPCseed(*tpcTrack, tpcTrack->GetAlpha());
1258 fSeeds->AddLast(seed);
1259 tracks.AddLast(new AliTPCtrack(*tpcTrack));
1264 // discard not prolongated tracks (to be discussed)
1266 Int_t lab = TMath::Abs(tpcTrack->GetLabel());
1267 if (lab > nLab || lut[lab] < 0) continue;
1269 bckTree->GetEvent(lut[lab]);
1270 bckTrack->Rotate(alpha-bckTrack->GetAlpha());
1272 fSeeds->AddLast(new AliTPCseed(*bckTrack,bckTrack->GetAlpha()));
1273 tracks.AddLast(new AliTPCtrack(*tpcTrack));
1281 // tree name seedsTPCtoTRD as expected by TRD and as
1282 // discussed and decided in Strasbourg (May 2002)
1283 // [SR, GSI, 18.02.2003]
1285 sprintf(tName,"seedsTPCtoTRD_%d",GetEventNumber());
1286 TTree backTree(tName,"Tree with back propagated TPC tracks");
1288 AliTPCtrack *otrack=0;
1289 backTree.Branch("tracks","AliTPCtrack",&otrack,32000,0);
1294 Int_t nseed=fSeeds->GetEntriesFast();
1296 // loop changed [SR, 01.04.2003]
1298 for (i=0; i<nseed; i++) {
1301 if (IsStoringBarrel()) fBarrelArray->Clear();
1303 AliTPCseed *ps=(AliTPCseed*)fSeeds->UncheckedAt(i), &s=*ps;
1304 const AliTPCtrack *pt=(AliTPCtrack*)tracks.UncheckedAt(i), &t=*pt;
1307 ps->ResetNRotation();
1309 if (IsStoringBarrel()) StoreBarrelTrack(ps, 1, 1);
1311 // Load outer sectors
1315 Int_t nc=t.GetNumberOfClusters();
1316 //s.SetLabel(nc-1); //set number of the cluster to start with
1319 if (FollowBackProlongation(s,t)) UseClusters(&s);
1321 fSeeds->RemoveAt(i);
1325 if (IsStoringBarrel()) StoreBarrelTrack(ps, 2, 1);
1326 if (IsStoringBarrel()) StoreBarrelTrack(ps, 3, 1);
1328 // Load outer sectors
1332 nc=s.GetNumberOfClusters();
1334 Double_t alpha=s.GetAlpha() - fSectors->GetAlphaShift();
1335 if (alpha > 2.*TMath::Pi()) alpha -= 2.*TMath::Pi();
1336 if (alpha < 0. ) alpha += 2.*TMath::Pi();
1337 Int_t ns=Int_t(alpha/fSectors->GetAlpha())%fN;
1339 alpha =ns*fSectors->GetAlpha() + fSectors->GetAlphaShift();
1340 alpha-=s.GetAlpha();
1342 if (s.Rotate(alpha)) {
1343 if (FollowBackProlongation(s,t)) {
1344 if (s.GetNumberOfClusters() >= Int_t(0.4*nrows)) {
1346 s.SetLabel(t.GetLabel());
1349 // Propagate to outer reference plane for comparison reasons
1350 // reason for keeping fParam object [SR, GSI, 18.02.2003]
1352 ps->PropagateTo(fParam->GetOuterRadiusUp());
1356 if (IsStoringBarrel()) StoreBarrelTrack(ps, 4, 1);
1357 cerr<<found++<<'\r';
1362 if (IsStoringBarrel()) fBarrelTree->Fill();
1363 delete fSeeds->RemoveAt(i);
1371 if (IsStoringBarrel()) {
1373 fBarrelTree->Write();
1374 fBarrelFile->Flush();
1376 delete fBarrelArray;
1381 cerr<<"Number of seeds: "<<nseed<<endl;
1382 cerr<<"Number of back propagated ITS tracks: "<<bckN<<endl;
1383 cerr<<"Number of back propagated TPC tracks: "<<found<<endl;
1388 if (bckTree) delete bckTree; //Thanks to Mariana Bondila
1389 delete tpcTree; //Thanks to Mariana Bondila
1397 //_________________________________________________________________________
1398 AliCluster *AliTPCtracker::GetCluster(Int_t index) const {
1399 //--------------------------------------------------------------------
1400 // Return pointer to a given cluster
1401 //--------------------------------------------------------------------
1402 Int_t sec=(index&0xff000000)>>24;
1403 Int_t row=(index&0x00ff0000)>>16;
1404 Int_t ncl=(index&0x0000ffff)>>00;
1406 const AliTPCcluster *cl=0;
1409 cl=fInnerSec[sec][row].GetUnsortedCluster(ncl);
1412 cl=fOuterSec[sec][row].GetUnsortedCluster(ncl);
1415 return (AliCluster*)cl;
1418 //__________________________________________________________________________
1419 void AliTPCtracker::CookLabel(AliKalmanTrack *t, Float_t wrong) const {
1420 //--------------------------------------------------------------------
1421 //This function "cooks" a track label. If label<0, this track is fake.
1422 //--------------------------------------------------------------------
1423 Int_t noc=t->GetNumberOfClusters();
1424 Int_t *lb=new Int_t[noc];
1425 Int_t *mx=new Int_t[noc];
1426 AliCluster **clusters=new AliCluster*[noc];
1429 for (i=0; i<noc; i++) {
1431 Int_t index=t->GetClusterIndex(i);
1432 clusters[i]=GetCluster(index);
1435 Int_t lab=123456789;
1436 for (i=0; i<noc; i++) {
1437 AliCluster *c=clusters[i];
1438 lab=TMath::Abs(c->GetLabel(0));
1440 for (j=0; j<noc; j++) if (lb[j]==lab || mx[j]==0) break;
1446 for (i=0; i<noc; i++) if (mx[i]>max) {max=mx[i]; lab=lb[i];}
1448 for (i=0; i<noc; i++) {
1449 AliCluster *c=clusters[i];
1450 if (TMath::Abs(c->GetLabel(1)) == lab ||
1451 TMath::Abs(c->GetLabel(2)) == lab ) max++;
1454 if ((1.- Float_t(max)/noc) > wrong) lab=-lab;
1457 Int_t tail=Int_t(0.10*noc);
1459 for (i=1; i<=tail; i++) {
1460 AliCluster *c=clusters[noc-i];
1461 if (lab == TMath::Abs(c->GetLabel(0)) ||
1462 lab == TMath::Abs(c->GetLabel(1)) ||
1463 lab == TMath::Abs(c->GetLabel(2))) max++;
1465 if (max < Int_t(0.5*tail)) lab=-lab;
1475 //_________________________________________________________________________
1476 void AliTPCtracker::AliTPCSector::Setup(const AliTPCParam *par, Int_t f) {
1477 //-----------------------------------------------------------------------
1478 // Setup inner sector
1479 //-----------------------------------------------------------------------
1481 fAlpha=par->GetInnerAngle();
1482 fAlphaShift=par->GetInnerAngleShift();
1483 fPadPitchWidth=par->GetInnerPadPitchWidth();
1484 f1PadPitchLength=par->GetInnerPadPitchLength();
1485 f2PadPitchLength=f1PadPitchLength;
1487 fN=par->GetNRowLow();
1488 fRow=new AliTPCRow[fN];
1489 for (Int_t i=0; i<fN; i++) fRow[i].SetX(par->GetPadRowRadiiLow(i));
1491 fAlpha=par->GetOuterAngle();
1492 fAlphaShift=par->GetOuterAngleShift();
1493 fPadPitchWidth=par->GetOuterPadPitchWidth();
1494 f1PadPitchLength=par->GetOuter1PadPitchLength();
1495 f2PadPitchLength=par->GetOuter2PadPitchLength();
1497 fN=par->GetNRowUp();
1498 fRow=new AliTPCRow[fN];
1499 for (Int_t i=0; i<fN; i++){
1500 fRow[i].SetX(par->GetPadRowRadiiUp(i));
1505 //_________________________________________________________________________
1506 void AliTPCtracker::
1507 AliTPCRow::InsertCluster(const AliTPCcluster* c, Int_t sec, Int_t row) {
1508 //-----------------------------------------------------------------------
1509 // Insert a cluster into this pad row in accordence with its y-coordinate
1510 //-----------------------------------------------------------------------
1511 if (fN==kMaxClusterPerRow) {
1512 ::Error("InsertCluster","Too many clusters !");
1516 Int_t index=(((sec<<8)+row)<<16)+fN;
1519 fSize=kMaxClusterPerRow/8;
1520 fClusterArray=new AliTPCcluster[fSize];
1522 fClusterArray[0]=*c; fClusters[fN++]=fClusterArray;
1528 AliTPCcluster *buff=new AliTPCcluster[size];
1529 memcpy(buff,fClusterArray,fSize*sizeof(AliTPCcluster));
1530 for (Int_t i=0; i<fN; i++)
1531 fClusters[i]=buff+(fClusters[i]-fClusterArray);
1532 delete[] fClusterArray;
1537 Int_t i=Find(c->GetY());
1538 memmove(fClusters+i+1 ,fClusters+i,(fN-i)*sizeof(AliTPCcluster*));
1539 memmove(fIndex +i+1 ,fIndex +i,(fN-i)*sizeof(UInt_t));
1541 fClusters[i]=fClusterArray+fN; fClusterArray[fN++]=*c;
1544 //___________________________________________________________________
1545 Int_t AliTPCtracker::AliTPCRow::Find(Double_t y) const {
1546 //-----------------------------------------------------------------------
1547 // Return the index of the nearest cluster
1548 //-----------------------------------------------------------------------
1550 if (y <= fClusters[0]->GetY()) return 0;
1551 if (y > fClusters[fN-1]->GetY()) return fN;
1552 Int_t b=0, e=fN-1, m=(b+e)/2;
1553 for (; b<e; m=(b+e)/2) {
1554 if (y > fClusters[m]->GetY()) b=m+1;
1560 //_____________________________________________________________________________
1561 void AliTPCtracker::AliTPCseed::CookdEdx(Double_t low, Double_t up) {
1562 //-----------------------------------------------------------------
1563 // This funtion calculates dE/dX within the "low" and "up" cuts.
1564 //-----------------------------------------------------------------
1566 Int_t nc=GetNumberOfClusters();
1568 Int_t swap;//stupid sorting
1571 for (i=0; i<nc-1; i++) {
1572 if (fdEdxSample[i]<=fdEdxSample[i+1]) continue;
1573 Float_t tmp=fdEdxSample[i];
1574 fdEdxSample[i]=fdEdxSample[i+1]; fdEdxSample[i+1]=tmp;
1579 Int_t nl=Int_t(low*nc), nu=Int_t(up*nc);
1581 for (i=nl; i<=nu; i++) dedx += fdEdxSample[i];
1586 Double_t p=TMath::Sqrt((1.+ GetTgl()*GetTgl())/(Get1Pt()*Get1Pt()));
1588 Double_t log1=TMath::Log(p+0.45), log2=TMath::Log(p+0.12);
1590 if (dedx < 34 + 30/(p+0.45)/(p+0.45) + 24*log1) {SetMass(0.13957); return;}
1591 if (dedx < 34 + 30/(p+0.12)/(p+0.12) + 24*log2) {SetMass(0.49368); return;}
1592 SetMass(0.93827); return;
1596 if (dedx < 34 + 30/(p+0.12)/(p+0.12) + 24*log2) {SetMass(0.13957); return;}
1597 SetMass(0.93827); return;
1600 SetMass(0.13957); return;