Return the most probable momentum if the magnetic field is too weak and the momentum...
[u/mrichter/AliRoot.git] / STEER / AliESDtrack.cxx
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 //           Implementation of the ESD track class
17 //   ESD = Event Summary Data
18 //   This is the class to deal with during the phisics analysis of data
19 //      Origin: Iouri Belikov, CERN
20 //      e-mail: Jouri.Belikov@cern.ch
21 //-----------------------------------------------------------------
22
23 #include <TMath.h>
24
25 #include "AliESDVertex.h"
26 #include "AliESDtrack.h"
27 #include "AliKalmanTrack.h"
28 #include "AliLog.h"
29 #include "AliTrackPointArray.h"
30
31 ClassImp(AliESDtrack)
32
33 void SetPIDValues(Float_t * dest, const Double_t * src, Int_t n) {
34   // This function copies "n" PID weights from "scr" to "dest"
35   // and normalizes their sum to 1 thus producing conditional probabilities.
36   // The negative weights are set to 0.
37   // In case all the weights are non-positive they are replaced by
38   // uniform probabilities
39
40   if (n<=0) return;
41
42   Float_t uniform = 1./(Float_t)n;
43
44   Float_t sum = 0;
45   for (Int_t i=0; i<n; i++) 
46     if (src[i]>=0) {
47       sum+=src[i];
48       dest[i] = src[i];
49     }
50     else {
51       dest[i] = 0;
52     }
53
54   if(sum>0)
55     for (Int_t i=0; i<n; i++) dest[i] /= sum;
56   else
57     for (Int_t i=0; i<n; i++) dest[i] = uniform;
58 }
59
60 //_______________________________________________________________________
61 AliESDtrack::AliESDtrack() : 
62   AliExternalTrackParam(),
63   fFlags(0),
64   fLabel(0),
65   fID(0),
66   fTrackLength(0),
67   fD(0),fZ(0),
68   fCdd(0),fCdz(0),fCzz(0),
69   fStopVertex(0),
70   fCp(0),
71   fCchi2(1e10),
72   fIp(0),
73   fOp(0),
74   fITSchi2(0),
75   fITSncls(0),
76   fITSsignal(0),
77   fITSLabel(0),
78   fITSFakeRatio(0),
79   fTPCchi2(0),
80   fTPCncls(0),
81   fTPCnclsF(0),
82   fTPCClusterMap(159),//number of padrows
83   fTPCsignal(0),
84   fTPCsignalN(0),
85   fTPCsignalS(0),
86   fTPCLabel(0),
87   fTRDchi2(0),
88   fTRDncls(0),
89   fTRDncls0(0),
90   fTRDsignal(0),
91   fTRDLabel(0),
92   fTRDQuality(0),
93   fTRDBudget(0),
94   fTOFchi2(0),
95   fTOFindex(0),
96   fTOFCalChannel(-1),
97   fTOFsignal(-1),
98   fTOFsignalToT(0),
99   fRICHchi2(1e10),
100   fRICHncls(0),
101   fRICHindex(0),
102   fRICHsignal(-1),
103   fRICHtheta(-1),
104   fRICHphi(-1),
105   fRICHdx(-1),
106   fRICHdy(-1),
107   fRICHmipX(-1),
108   fRICHmipY(-1),
109   fPoints(0),
110   fFriendTrack(new AliESDfriendTrack())
111 {
112   //
113   // The default ESD constructor 
114   //
115   Int_t i, j;
116   for (i=0; i<AliPID::kSPECIES; i++) {
117     fTrackTime[i]=0.;
118     fR[i]=1.;
119     fITSr[i]=1.;
120     fTPCr[i]=1.;
121     fTRDr[i]=1.;
122     fTOFr[i]=1.;
123     fRICHr[i]=1.;
124   }
125   
126   for (i=0; i<3; i++)   { fKinkIndexes[i]=0;}
127   for (i=0; i<3; i++)   { fV0Indexes[i]=-1;}
128   for (i=0;i<kNPlane;i++) {
129     for (j=0;j<kNSlice;j++) {
130       fTRDsignals[i][j]=0.; 
131     }
132     fTRDTimBin[i]=-1;
133   }
134   for (i=0;i<4;i++) {fTPCPoints[i]=-1;}
135   for (i=0;i<3;i++) {fTOFLabel[i]=-1;}
136   for (i=0;i<10;i++) {fTOFInfo[i]=-1;}
137 }
138
139 //_______________________________________________________________________
140 AliESDtrack::AliESDtrack(const AliESDtrack& track):
141   AliExternalTrackParam(track),
142   fFlags(track.fFlags),
143   fLabel(track.fLabel),
144   fID(track.fID),
145   fTrackLength(track.fTrackLength),
146   fD(track.fD),fZ(track.fZ),
147   fCdd(track.fCdd),fCdz(track.fCdz),fCzz(track.fCzz),
148   fStopVertex(track.fStopVertex),
149   fCp(0),
150   fCchi2(track.fCchi2),
151   fIp(0),
152   fOp(0),
153   fITSchi2(track.fITSchi2),
154   fITSncls(track.fITSncls),
155   fITSsignal(track.fITSsignal),
156   fITSLabel(track.fITSLabel),
157   fITSFakeRatio(track.fITSFakeRatio),
158   fTPCchi2(track.fTPCchi2),
159   fTPCncls(track.fTPCncls),
160   fTPCnclsF(track.fTPCnclsF),
161   fTPCClusterMap(track.fTPCClusterMap),
162   fTPCsignal(track.fTPCsignal),
163   fTPCsignalN(track.fTPCsignalN),
164   fTPCsignalS(track.fTPCsignalS),
165   fTPCLabel(track.fTPCLabel),
166   fTRDchi2(track.fTRDchi2),
167   fTRDncls(track.fTRDncls),
168   fTRDncls0(track.fTRDncls0),
169   fTRDsignal(track.fTRDsignal),
170   fTRDLabel(track.fTRDLabel),
171   fTRDQuality(track.fTRDQuality),
172   fTRDBudget(track.fTRDBudget),
173   fTOFchi2(track.fTOFchi2),
174   fTOFindex(track.fTOFindex),
175   fTOFCalChannel(track.fTOFCalChannel),
176   fTOFsignal(track.fTOFsignal),
177   fTOFsignalToT(track.fTOFsignalToT),
178   fRICHchi2(track.fRICHchi2),
179   fRICHncls(track.fRICHncls),
180   fRICHindex(track.fRICHindex),
181   fRICHsignal(track.fRICHsignal),
182   fRICHtheta(track.fRICHtheta),
183   fRICHphi(track.fRICHphi),
184   fRICHdx(track.fRICHdx),
185   fRICHdy(track.fRICHdy),
186   fRICHmipX(track.fRICHmipX),
187   fRICHmipY(track.fRICHmipY),
188   fPoints(0),
189   fFriendTrack(new AliESDfriendTrack(*(track.fFriendTrack)))
190 {
191   //
192   //copy constructor
193   //
194   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTrackTime[i]=track.fTrackTime[i];
195   for (Int_t i=0;i<AliPID::kSPECIES;i++)  fR[i]=track.fR[i];
196   //
197   for (Int_t i=0;i<AliPID::kSPECIES;i++) fITSr[i]=track.fITSr[i]; 
198   //
199   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTPCr[i]=track.fTPCr[i]; 
200   for (Int_t i=0;i<4;i++) {fTPCPoints[i]=track.fTPCPoints[i];}
201   for (Int_t i=0; i<3;i++)   { fKinkIndexes[i]=track.fKinkIndexes[i];}
202   for (Int_t i=0; i<3;i++)   { fV0Indexes[i]=track.fV0Indexes[i];}
203   //
204   for (Int_t i=0;i<kNPlane;i++) {
205     for (Int_t j=0;j<kNSlice;j++) {
206       fTRDsignals[i][j]=track.fTRDsignals[i][j]; 
207     }
208     fTRDTimBin[i]=track.fTRDTimBin[i];
209   }
210   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTRDr[i]=track.fTRDr[i]; 
211   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTOFr[i]=track.fTOFr[i];
212   for (Int_t i=0;i<3;i++) fTOFLabel[i]=track.fTOFLabel[i];
213   for (Int_t i=0;i<10;i++) fTOFInfo[i]=track.fTOFInfo[i];
214   for (Int_t i=0;i<AliPID::kSPECIES;i++) fRICHr[i]=track.fRICHr[i];
215
216   if (track.fCp) fCp=new AliExternalTrackParam(*track.fCp);
217   if (track.fIp) fIp=new AliExternalTrackParam(*track.fIp);
218   if (track.fOp) fOp=new AliExternalTrackParam(*track.fOp);
219
220   if (track.fPoints) fPoints=new AliTrackPointArray(*(track.fPoints));
221 }
222
223 //_______________________________________________________________________
224 AliESDtrack::~AliESDtrack(){ 
225   //
226   // This is destructor according Coding Conventrions 
227   //
228   //printf("Delete track\n");
229   delete fIp; 
230   delete fOp;
231   delete fCp; 
232   delete fFriendTrack;
233   delete fPoints;
234 }
235
236 //_______________________________________________________________________
237 void AliESDtrack::MakeMiniESDtrack(){
238   // Resets everything except
239   // fFlags: Reconstruction status flags 
240   // fLabel: Track label
241   // fID:  Unique ID of the track
242   // fD: Impact parameter in XY-plane
243   // fZ: Impact parameter in Z 
244   // fR[AliPID::kSPECIES]: combined "detector response probability"
245   // Running track parameters
246   // fRalpha: track rotation angle
247   // fRx: X-coordinate of the track reference plane 
248   // fRp[5]: external track parameters  
249   // fRc[15]: external cov. matrix of the track parameters
250   
251   fTrackLength = 0;
252   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTrackTime[i] = 0;
253   fStopVertex = 0;
254
255   // Reset track parameters constrained to the primary vertex
256   fCp = 0;
257   fCchi2 = 0;
258
259   // Reset track parameters at the inner wall of TPC
260   fIp = 0;
261
262   // Reset track parameters at the inner wall of the TRD
263   fOp = 0;
264
265   // Reset ITS track related information
266   fITSchi2 = 0;
267   fITSncls = 0;       
268   fITSsignal = 0;     
269   for (Int_t i=0;i<AliPID::kSPECIES;i++) fITSr[i]=0; 
270   fITSLabel = 0;       
271   fITSFakeRatio = 0;   
272
273   // Reset TPC related track information
274   fTPCchi2 = 0;       
275   fTPCncls = 0;       
276   fTPCnclsF = 0;       
277   fTPCClusterMap = 0;  
278   fTPCsignal= 0;      
279   fTPCsignalS= 0;      
280   fTPCsignalN= 0;      
281   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTPCr[i]=0; 
282   fTPCLabel=0;       
283   for (Int_t i=0;i<4;i++) fTPCPoints[i] = 0;
284   for (Int_t i=0; i<3;i++)   fKinkIndexes[i] = 0;
285   for (Int_t i=0; i<3;i++)   fV0Indexes[i] = 0;
286
287   // Reset TRD related track information
288   fTRDchi2 = 0;        
289   fTRDncls = 0;       
290   fTRDncls0 = 0;       
291   fTRDsignal = 0;      
292   for (Int_t i=0;i<kNPlane;i++) {
293     for (Int_t j=0;j<kNSlice;j++) {
294       fTRDsignals[i][j] = 0; 
295     }
296     fTRDTimBin[i]  = 0;
297   }
298   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTRDr[i] = 0; 
299   fTRDLabel = 0;       
300   fTRDQuality  = 0;
301   fTRDBudget  = 0;
302
303   // Reset TOF related track information
304   fTOFchi2 = 0;        
305   fTOFindex = 0;       
306   fTOFsignal = 0;      
307   fTOFCalChannel = -1;
308   fTOFsignalToT = 0;
309   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTOFr[i] = 0;
310   for (Int_t i=0;i<3;i++) fTOFLabel[i] = 0;
311   for (Int_t i=0;i<10;i++) fTOFInfo[i] = 0;
312
313   // Reset RICH related track information
314   fRICHchi2 = 0;     
315   fRICHncls = 0;     
316   fRICHindex = 0;     
317   fRICHsignal = 0;     
318   for (Int_t i=0;i<AliPID::kSPECIES;i++) fRICHr[i] = 0;
319   fRICHtheta = 0;     
320   fRICHphi = 0;      
321   fRICHdx = 0;     
322   fRICHdy = 0;      
323
324   delete fFriendTrack; fFriendTrack = 0;
325   delete fPoints; fPoints = 0;
326
327 //_______________________________________________________________________
328 Double_t AliESDtrack::GetMass() const {
329   // Returns the mass of the most probable particle type
330   Float_t max=0.;
331   Int_t k=-1;
332   for (Int_t i=0; i<AliPID::kSPECIES; i++) {
333     if (fR[i]>max) {k=i; max=fR[i];}
334   }
335   if (k==0) { // dE/dx "crossing points" in the TPC
336      Double_t p=GetP();
337      if ((p>0.38)&&(p<0.48))
338         if (fR[0]<fR[3]*10.) return AliPID::ParticleMass(AliPID::kKaon);
339      if ((p>0.75)&&(p<0.85))
340         if (fR[0]<fR[4]*10.) return AliPID::ParticleMass(AliPID::kProton);
341      return 0.00051;
342   }
343   if (k==1) return AliPID::ParticleMass(AliPID::kMuon); 
344   if (k==2||k==-1) return AliPID::ParticleMass(AliPID::kPion);
345   if (k==3) return AliPID::ParticleMass(AliPID::kKaon);
346   if (k==4) return AliPID::ParticleMass(AliPID::kProton);
347   AliWarning("Undefined mass !");
348   return AliPID::ParticleMass(AliPID::kPion);
349 }
350
351 //_______________________________________________________________________
352 Bool_t AliESDtrack::UpdateTrackParams(const AliKalmanTrack *t, ULong_t flags){
353   //
354   // This function updates track's running parameters 
355   //
356   Int_t *index=0;
357   Bool_t rc=kTRUE;
358
359   SetStatus(flags);
360   fLabel=t->GetLabel();
361
362   if (t->IsStartedTimeIntegral()) {
363     SetStatus(kTIME);
364     Double_t times[10];t->GetIntegratedTimes(times); SetIntegratedTimes(times);
365     SetIntegratedLength(t->GetIntegratedLength());
366   }
367
368   Set(*t);
369   
370   switch (flags) {
371     
372   case kITSin: case kITSout: case kITSrefit:
373     index=fFriendTrack->GetITSindices(); 
374     for (Int_t i=0;i<AliESDfriendTrack::kMaxITScluster;i++) 
375          index[i]=t->GetClusterIndex(i);
376     fITSncls=t->GetNumberOfClusters();
377     fITSchi2=t->GetChi2();
378     fITSsignal=t->GetPIDsignal();
379     fITSLabel = t->GetLabel();
380     fITSFakeRatio = t->GetFakeRatio();
381     break;
382     
383   case kTPCin: case kTPCrefit:
384     fTPCLabel = t->GetLabel();
385     if (!fIp) fIp=new AliExternalTrackParam(*t);
386     else fIp->Set(*t);
387   case kTPCout:
388     index=fFriendTrack->GetTPCindices(); 
389     if (flags & kTPCout){
390       if (!fOp) fOp=new AliExternalTrackParam(*t);
391       else fOp->Set(*t);
392     }
393     fTPCncls=t->GetNumberOfClusters();    
394     fTPCchi2=t->GetChi2();
395     
396      {//prevrow must be declared in separate namespace, otherwise compiler cries:
397       //"jump to case label crosses initialization of `Int_t prevrow'"
398        Int_t prevrow = -1;
399        //       for (Int_t i=0;i<fTPCncls;i++) 
400        for (Int_t i=0;i<AliESDfriendTrack::kMaxTPCcluster;i++) 
401         {
402           index[i]=t->GetClusterIndex(i);
403           Int_t idx = index[i];
404
405           if (idx<0) continue; 
406
407           // Piotr's Cluster Map for HBT  
408           // ### please change accordingly if cluster array is changing 
409           // to "New TPC Tracking" style (with gaps in array) 
410           Int_t sect = (idx&0xff000000)>>24;
411           Int_t row = (idx&0x00ff0000)>>16;
412           if (sect > 18) row +=63; //if it is outer sector, add number of inner sectors
413
414           fTPCClusterMap.SetBitNumber(row,kTRUE);
415
416           //Fill the gap between previous row and this row with 0 bits
417           //In case  ###  pleas change it as well - just set bit 0 in case there 
418           //is no associated clusters for current "i"
419           if (prevrow < 0) 
420            {
421              prevrow = row;//if previous bit was not assigned yet == this is the first one
422            }
423           else
424            { //we don't know the order (inner to outer or reverse)
425              //just to be save in case it is going to change
426              Int_t n = 0, m = 0;
427              if (prevrow < row)
428               {
429                 n = prevrow;
430                 m = row;
431               }
432              else
433               {
434                 n = row;
435                 m = prevrow;
436               }
437
438              for (Int_t j = n+1; j < m; j++)
439               {
440                 fTPCClusterMap.SetBitNumber(j,kFALSE);
441               }
442              prevrow = row; 
443            }
444           // End Of Piotr's Cluster Map for HBT
445         }
446      }
447     fTPCsignal=t->GetPIDsignal();
448     break;
449
450   case kTRDout: case kTRDin: case kTRDrefit:
451     index=fFriendTrack->GetTRDindices();
452     fTRDLabel = t->GetLabel(); 
453     fTRDncls=t->GetNumberOfClusters();
454     fTRDchi2=t->GetChi2();
455     for (Int_t i=0;i<fTRDncls;i++) index[i]=t->GetClusterIndex(i);
456     fTRDsignal=t->GetPIDsignal();
457     break;
458   case kTRDbackup:
459     if (!fOp) fOp=new AliExternalTrackParam(*t);
460     else fOp->Set(*t);
461     fTRDncls0 = t->GetNumberOfClusters(); 
462     break;
463   case kTOFin: 
464     break;
465   case kTOFout: 
466     break;
467   case kTRDStop:
468     break;
469   default: 
470     AliError("Wrong flag !");
471     return kFALSE;
472   }
473
474   return rc;
475 }
476
477 //_______________________________________________________________________
478 void AliESDtrack::GetExternalParameters(Double_t &x, Double_t p[5]) const {
479   //---------------------------------------------------------------------
480   // This function returns external representation of the track parameters
481   //---------------------------------------------------------------------
482   x=GetX();
483   for (Int_t i=0; i<5; i++) p[i]=GetParameter()[i];
484 }
485
486 //_______________________________________________________________________
487 void AliESDtrack::GetExternalCovariance(Double_t cov[15]) const {
488   //---------------------------------------------------------------------
489   // This function returns external representation of the cov. matrix
490   //---------------------------------------------------------------------
491   for (Int_t i=0; i<15; i++) cov[i]=AliExternalTrackParam::GetCovariance()[i];
492 }
493
494 //_______________________________________________________________________
495 Bool_t AliESDtrack::GetConstrainedExternalParameters
496                  (Double_t &alpha, Double_t &x, Double_t p[5]) const {
497   //---------------------------------------------------------------------
498   // This function returns the constrained external track parameters
499   //---------------------------------------------------------------------
500   if (!fCp) return kFALSE;
501   alpha=fCp->GetAlpha();
502   x=fCp->GetX();
503   for (Int_t i=0; i<5; i++) p[i]=fCp->GetParameter()[i];
504   return kTRUE;
505 }
506
507 //_______________________________________________________________________
508 Bool_t 
509 AliESDtrack::GetConstrainedExternalCovariance(Double_t c[15]) const {
510   //---------------------------------------------------------------------
511   // This function returns the constrained external cov. matrix
512   //---------------------------------------------------------------------
513   if (!fCp) return kFALSE;
514   for (Int_t i=0; i<15; i++) c[i]=fCp->GetCovariance()[i];
515   return kTRUE;
516 }
517
518 Bool_t
519 AliESDtrack::GetInnerExternalParameters
520                  (Double_t &alpha, Double_t &x, Double_t p[5]) const {
521   //---------------------------------------------------------------------
522   // This function returns external representation of the track parameters 
523   // at the inner layer of TPC
524   //---------------------------------------------------------------------
525   if (!fIp) return kFALSE;
526   alpha=fIp->GetAlpha();
527   x=fIp->GetX();
528   for (Int_t i=0; i<5; i++) p[i]=fIp->GetParameter()[i];
529   return kTRUE;
530 }
531
532 Bool_t 
533 AliESDtrack::GetInnerExternalCovariance(Double_t cov[15]) const {
534  //---------------------------------------------------------------------
535  // This function returns external representation of the cov. matrix 
536  // at the inner layer of TPC
537  //---------------------------------------------------------------------
538   if (!fIp) return kFALSE;
539   for (Int_t i=0; i<15; i++) cov[i]=fIp->GetCovariance()[i];
540   return kTRUE;
541 }
542
543 Bool_t 
544 AliESDtrack::GetOuterExternalParameters
545                  (Double_t &alpha, Double_t &x, Double_t p[5]) const {
546   //---------------------------------------------------------------------
547   // This function returns external representation of the track parameters 
548   // at the inner layer of TRD
549   //---------------------------------------------------------------------
550   if (!fOp) return kFALSE;
551   alpha=fOp->GetAlpha();
552   x=fOp->GetX();
553   for (Int_t i=0; i<5; i++) p[i]=fOp->GetParameter()[i];
554   return kTRUE;
555 }
556
557 Bool_t 
558 AliESDtrack::GetOuterExternalCovariance(Double_t cov[15]) const {
559  //---------------------------------------------------------------------
560  // This function returns external representation of the cov. matrix 
561  // at the inner layer of TRD
562  //---------------------------------------------------------------------
563   if (!fOp) return kFALSE;
564   for (Int_t i=0; i<15; i++) cov[i]=fOp->GetCovariance()[i];
565   return kTRUE;
566 }
567
568 Int_t AliESDtrack::GetNcls(Int_t idet) const
569 {
570   // Get number of clusters by subdetector index
571   //
572   Int_t ncls = 0;
573   switch(idet){
574   case 0:
575     ncls = fITSncls;
576     break;
577   case 1:
578     ncls = fTPCncls;
579     break;
580   case 2:
581     ncls = fTRDncls;
582     break;
583   case 3:
584     if (fTOFindex != 0)
585       ncls = 1;
586     break;
587   default:
588     break;
589   }
590   return ncls;
591 }
592
593 Int_t AliESDtrack::GetClusters(Int_t idet, Int_t *idx) const
594 {
595   // Get cluster index array by subdetector index
596   //
597   Int_t ncls = 0;
598   switch(idet){
599   case 0:
600     ncls = GetITSclusters(idx);
601     break;
602   case 1:
603     ncls = GetTPCclusters(idx);
604     break;
605   case 2:
606     ncls = GetTRDclusters(idx);
607     break;
608   case 3:
609     if (fTOFindex != 0) {
610       idx[0] = GetTOFcluster();
611       ncls = 1;
612     }
613     break;
614   default:
615     break;
616   }
617   return ncls;
618 }
619
620 //_______________________________________________________________________
621 void AliESDtrack::GetIntegratedTimes(Double_t *times) const {
622   // Returns the array with integrated times for each particle hypothesis
623   for (Int_t i=0; i<AliPID::kSPECIES; i++) times[i]=fTrackTime[i];
624 }
625
626 //_______________________________________________________________________
627 void AliESDtrack::SetIntegratedTimes(const Double_t *times) {
628   // Sets the array with integrated times for each particle hypotesis
629   for (Int_t i=0; i<AliPID::kSPECIES; i++) fTrackTime[i]=times[i];
630 }
631
632 //_______________________________________________________________________
633 void AliESDtrack::SetITSpid(const Double_t *p) {
634   // Sets values for the probability of each particle type (in ITS)
635   SetPIDValues(fITSr,p,AliPID::kSPECIES);
636   SetStatus(AliESDtrack::kITSpid);
637 }
638
639 //_______________________________________________________________________
640 void AliESDtrack::GetITSpid(Double_t *p) const {
641   // Gets the probability of each particle type (in ITS)
642   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fITSr[i];
643 }
644
645 //_______________________________________________________________________
646 Int_t AliESDtrack::GetITSclusters(Int_t *idx) const {
647   //---------------------------------------------------------------------
648   // This function returns indices of the assgined ITS clusters 
649   //---------------------------------------------------------------------
650   if (idx!=0) {
651      Int_t *index=fFriendTrack->GetITSindices();
652      for (Int_t i=0; i<AliESDfriendTrack::kMaxITScluster; i++) idx[i]=index[i];
653   }
654   return fITSncls;
655 }
656
657 //_______________________________________________________________________
658 Int_t AliESDtrack::GetTPCclusters(Int_t *idx) const {
659   //---------------------------------------------------------------------
660   // This function returns indices of the assgined ITS clusters 
661   //---------------------------------------------------------------------
662   if (idx!=0) {
663     Int_t *index=fFriendTrack->GetTPCindices();
664     for (Int_t i=0; i<AliESDfriendTrack::kMaxTPCcluster; i++) idx[i]=index[i];
665   }
666   return fTPCncls;
667 }
668
669 Float_t AliESDtrack::GetTPCdensity(Int_t row0, Int_t row1) const{
670   //
671   // GetDensity of the clusters on given region between row0 and row1
672   // Dead zone effect takin into acoount
673   //
674   Int_t good  = 0;
675   Int_t found = 0;
676   //  
677   Int_t *index=fFriendTrack->GetTPCindices();
678   for (Int_t i=row0;i<=row1;i++){     
679     Int_t idx = index[i];
680     if (idx!=-1)  good++;             // track outside of dead zone
681     if (idx>0)    found++;
682   }
683   Float_t density=0.5;
684   if (good>(row1-row0)*0.5) density = Float_t(found)/Float_t(good);
685   return density;
686 }
687
688 //_______________________________________________________________________
689 void AliESDtrack::SetTPCpid(const Double_t *p) {  
690   // Sets values for the probability of each particle type (in TPC)
691   SetPIDValues(fTPCr,p,AliPID::kSPECIES);
692   SetStatus(AliESDtrack::kTPCpid);
693 }
694
695 //_______________________________________________________________________
696 void AliESDtrack::GetTPCpid(Double_t *p) const {
697   // Gets the probability of each particle type (in TPC)
698   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fTPCr[i];
699 }
700
701 //_______________________________________________________________________
702 Int_t AliESDtrack::GetTRDclusters(Int_t *idx) const {
703   //---------------------------------------------------------------------
704   // This function returns indices of the assgined TRD clusters 
705   //---------------------------------------------------------------------
706   if (idx!=0) {
707      Int_t *index=fFriendTrack->GetTRDindices();
708      for (Int_t i=0; i<AliESDfriendTrack::kMaxTRDcluster; i++) idx[i]=index[i];
709   }
710   return fTRDncls;
711 }
712
713 //_______________________________________________________________________
714 void AliESDtrack::SetTRDpid(const Double_t *p) {  
715   // Sets values for the probability of each particle type (in TRD)
716   SetPIDValues(fTRDr,p,AliPID::kSPECIES);
717   SetStatus(AliESDtrack::kTRDpid);
718 }
719
720 //_______________________________________________________________________
721 void AliESDtrack::GetTRDpid(Double_t *p) const {
722   // Gets the probability of each particle type (in TRD)
723   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fTRDr[i];
724 }
725
726 //_______________________________________________________________________
727 void    AliESDtrack::SetTRDpid(Int_t iSpecies, Float_t p)
728 {
729   // Sets the probability of particle type iSpecies to p (in TRD)
730   fTRDr[iSpecies] = p;
731 }
732
733 Float_t AliESDtrack::GetTRDpid(Int_t iSpecies) const
734 {
735   // Returns the probability of particle type iSpecies (in TRD)
736   return fTRDr[iSpecies];
737 }
738
739 //_______________________________________________________________________
740 void AliESDtrack::SetTOFpid(const Double_t *p) {  
741   // Sets the probability of each particle type (in TOF)
742   SetPIDValues(fTOFr,p,AliPID::kSPECIES);
743   SetStatus(AliESDtrack::kTOFpid);
744 }
745
746 //_______________________________________________________________________
747 void AliESDtrack::SetTOFLabel(const Int_t *p) {  
748   // Sets  (in TOF)
749   for (Int_t i=0; i<3; i++) fTOFLabel[i]=p[i];
750 }
751
752 //_______________________________________________________________________
753 void AliESDtrack::GetTOFpid(Double_t *p) const {
754   // Gets probabilities of each particle type (in TOF)
755   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fTOFr[i];
756 }
757
758 //_______________________________________________________________________
759 void AliESDtrack::GetTOFLabel(Int_t *p) const {
760   // Gets (in TOF)
761   for (Int_t i=0; i<3; i++) p[i]=fTOFLabel[i];
762 }
763
764 //_______________________________________________________________________
765 void AliESDtrack::GetTOFInfo(Float_t *info) const {
766   // Gets (in TOF)
767   for (Int_t i=0; i<10; i++) info[i]=fTOFInfo[i];
768 }
769
770 //_______________________________________________________________________
771 void AliESDtrack::SetTOFInfo(Float_t*info) {
772   // Gets (in TOF)
773   for (Int_t i=0; i<10; i++) fTOFInfo[i]=info[i];
774 }
775
776
777
778 //_______________________________________________________________________
779 void AliESDtrack::SetRICHpid(const Double_t *p) {  
780   // Sets the probability of each particle type (in RICH)
781   SetPIDValues(fRICHr,p,AliPID::kSPECIES);
782   SetStatus(AliESDtrack::kRICHpid);
783 }
784
785 //_______________________________________________________________________
786 void AliESDtrack::GetRICHpid(Double_t *p) const {
787   // Gets probabilities of each particle type (in RICH)
788   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fRICHr[i];
789 }
790
791
792
793 //_______________________________________________________________________
794 void AliESDtrack::SetESDpid(const Double_t *p) {  
795   // Sets the probability of each particle type for the ESD track
796   SetPIDValues(fR,p,AliPID::kSPECIES);
797   SetStatus(AliESDtrack::kESDpid);
798 }
799
800 //_______________________________________________________________________
801 void AliESDtrack::GetESDpid(Double_t *p) const {
802   // Gets probability of each particle type for the ESD track
803   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fR[i];
804 }
805
806 //_______________________________________________________________________
807 Bool_t AliESDtrack::RelateToVertex
808 (const AliESDVertex *vtx, Double_t b, Double_t maxd) {
809   //
810   // Try to relate this track to the vertex "vtx", 
811   // if the (rough) transverse impact parameter is not bigger then "maxd". 
812   //            Magnetic field is "b" (kG).
813   //
814   // a) The track gets extapolated to the DCA to the vertex.
815   // b) The impact parameters and their covariance matrix are calculated.
816   // c) An attempt to constrain this track to the vertex is done.
817   //
818   //    In the case of success, the returned value is kTRUE
819   //    (otherwise, it's kFALSE)
820   //  
821   Double_t alpha=GetAlpha();
822   Double_t sn=TMath::Sin(alpha), cs=TMath::Cos(alpha);
823   Double_t x=GetX(), y=GetParameter()[0], snp=GetParameter()[2];
824   Double_t xv= vtx->GetXv()*cs + vtx->GetYv()*sn;
825   Double_t yv=-vtx->GetXv()*sn + vtx->GetYv()*cs, zv=vtx->GetZv();
826   x-=xv; y-=yv;
827
828   //Estimate the impact parameter neglecting the track curvature
829   Double_t d=TMath::Abs(x*snp - y*TMath::Sqrt(1.- snp*snp));
830   if (d > maxd) return kFALSE; 
831
832   //Propagate to the DCA
833   Double_t crv=kB2C*b*GetParameter()[4];
834   if (TMath::Abs(b) < kAlmost0Field) crv=0.;
835
836   Double_t tgfv=-(crv*x - snp)/(crv*y + TMath::Sqrt(1.-snp*snp));
837   sn=tgfv/TMath::Sqrt(1.+ tgfv*tgfv);
838   if (TMath::Abs(tgfv)>0.) cs = sn/tgfv;
839   else cs=1.;
840
841   x = xv*cs + yv*sn;
842   yv=-xv*sn + yv*cs; xv=x;
843
844   if (!Propagate(alpha+TMath::ASin(sn),xv,b)) return kFALSE;
845
846   fD = GetParameter()[0] - yv;
847   fZ = GetParameter()[1] - zv;
848   
849   Double_t cov[6]; vtx->GetCovMatrix(cov);
850   fCdd = GetCovariance()[0] + cov[2];      // neglecting non-diagonals
851   fCdz = GetCovariance()[1];               //     in the vertex's    
852   fCzz = GetCovariance()[2] + cov[5];      //    covariance matrix
853
854   {//Try to constrain 
855     Double_t p[2]={yv,zv}, c[3]={cov[2],0.,cov[5]};
856     Double_t chi2=GetPredictedChi2(p,c);
857
858     if (chi2>77.) return kFALSE;
859
860     AliExternalTrackParam tmp(*this);
861     if (!tmp.Update(p,c)) return kFALSE;
862
863     fCchi2=chi2;
864     if (!fCp) fCp=new AliExternalTrackParam();
865     new (fCp) AliExternalTrackParam(tmp);
866   }
867
868   return kTRUE;
869 }
870
871 void AliESDtrack::SetTrackPointArray(AliTrackPointArray *points) { 
872     fPoints=points;
873     //fFriendTrack->SetTrackPointArray(points); 
874 }
875 const AliTrackPointArray *AliESDtrack::GetTrackPointArray() const { 
876     return fPoints;
877   //return fFriendTrack->GetTrackPointArray(); 
878 }
879
880 //_______________________________________________________________________
881 void AliESDtrack::Print(Option_t *) const {
882   // Prints info on the track
883   
884   printf("ESD track info\n") ; 
885   Double_t p[AliPID::kSPECIESN] ; 
886   Int_t index = 0 ; 
887   if( IsOn(kITSpid) ){
888     printf("From ITS: ") ; 
889     GetITSpid(p) ; 
890     for(index = 0 ; index < AliPID::kSPECIES; index++) 
891       printf("%f, ", p[index]) ;
892     printf("\n           signal = %f\n", GetITSsignal()) ;
893   } 
894   if( IsOn(kTPCpid) ){
895     printf("From TPC: ") ; 
896     GetTPCpid(p) ; 
897     for(index = 0 ; index < AliPID::kSPECIES; index++) 
898       printf("%f, ", p[index]) ;
899     printf("\n           signal = %f\n", GetTPCsignal()) ;
900   }
901   if( IsOn(kTRDpid) ){
902     printf("From TRD: ") ; 
903     GetTRDpid(p) ; 
904     for(index = 0 ; index < AliPID::kSPECIES; index++) 
905       printf("%f, ", p[index]) ;
906     printf("\n           signal = %f\n", GetTRDsignal()) ;
907   }
908   if( IsOn(kTOFpid) ){
909     printf("From TOF: ") ; 
910     GetTOFpid(p) ; 
911     for(index = 0 ; index < AliPID::kSPECIES; index++) 
912       printf("%f, ", p[index]) ;
913     printf("\n           signal = %f\n", GetTOFsignal()) ;
914   }
915   if( IsOn(kRICHpid) ){
916     printf("From RICH: ") ; 
917     GetRICHpid(p) ; 
918     for(index = 0 ; index < AliPID::kSPECIES; index++) 
919       printf("%f, ", p[index]) ;
920     printf("\n           signal = %f\n", GetRICHsignal()) ;
921   }
922