0a148c2b643d834a643ba579f19b5015c1f166f1
[u/mrichter/AliRoot.git] / STEER / ESD / 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 //
24 //  What do you need to know before starting analysis
25 //  (by Marian Ivanov: marian.ivanov@cern.ch)
26 //
27 //
28 //   AliESDtrack:
29 //   1.  What is the AliESDtrack
30 //   2.  What informations do we store
31 //   3.  How to use the information for analysis
32 //   
33 //
34 //   1.AliESDtrack is the container of the information about the track/particle
35 //     reconstructed during Barrel Tracking.
36 //     The track information is propagated from one tracking detector to 
37 //     other using the functionality of AliESDtrack - Current parameters.  
38 //
39 //     No global fit model is used.
40 //     Barrel tracking use Kalman filtering technique, it gives optimal local 
41 //     track parameters at given point under certian assumptions.
42 //             
43 //     Kalman filter take into account additional effect which are 
44 //     difficult to handle using global fit.
45 //     Effects:
46 //        a.) Multiple scattering
47 //        b.) Energy loss
48 //        c.) Non homogenous magnetic field
49 //
50 //     In general case, following barrel detectors are contributing to 
51 //     the Kalman track information:
52 //         a. TPC
53 //         b. ITS
54 //         c. TRD
55 //
56 //      In general 3 reconstruction itteration are performed:
57 //         1. Find tracks   - sequence TPC->ITS
58 //         2. PropagateBack - sequence ITS->TPC->TRD -> Outer PID detectors
59 //         3. Refit invward - sequence TRD->TPC->ITS
60 //      The current tracks are updated after each detector (see bellow).
61 //      In specical cases a track  sanpshots are stored.   
62 // 
63 //
64 //      For some type of analysis (+visualization) track local parameters at 
65 //      different position are neccesary. A snapshots during the track 
66 //      propagation are created.
67 //      (See AliExternalTrackParam class for desctiption of variables and 
68 //      functionality)
69 //      Snapshots:
70 //      a. Current parameters - class itself (AliExternalTrackParam)
71 //         Contributors: general case TRD->TPC->ITS
72 //         Preferable usage:  Decission  - primary or secondary track
73 //         NOTICE - By default the track parameters are stored at the DCA point
74 //                  to the primary vertex. optimal for primary tracks, 
75 //                  far from optimal for secondary tracks.
76 //      b. Constrained parameters - Kalman information updated with 
77 //         the Primary vertex information 
78 //         Contributors: general case TRD->TPC->ITS
79 //         Preferable usage: Use only for tracks selected as primary
80 //         NOTICE - not real constrain - taken as additional measurement 
81 //         with corresponding error
82 //         Function:  
83 //       const AliExternalTrackParam *GetConstrainedParam() const {return fCp;}
84 //      c. Inner parameters -  Track parameters at inner wall of the TPC 
85 //         Contributors: general case TRD->TPC
86 //         function:
87 //           const AliExternalTrackParam *GetInnerParam() const { return fIp;}
88 //
89 //      d. TPCinnerparam  - contributors - TPC only
90 //         Contributors:  TPC
91 //         Preferable usage: Requested for HBT study 
92 //                         (smaller correlations as using also ITS information)
93 //         NOTICE - the track parameters are propagated to the DCA to  
94 //         to primary vertex
95 //         Optimal for primary, far from optimal for secondary tracks
96 //         Function:
97 //    const AliExternalTrackParam *GetTPCInnerParam() const {return fTPCInner;}
98 //     
99 //      e. Outer parameters - 
100 //           Contributors-  general case - ITS-> TPC -> TRD
101 //           The last point - Outer parameters radius is determined
102 //           e.a) Local inclination angle bigger than threshold - 
103 //                Low momenta tracks 
104 //           e.a) Catastrofic energy losss in material
105 //           e.b) Not further improvement (no space points)
106 //           Usage:             
107 //              a.) Tracking: Starting parameter for Refit inward 
108 //              b.) Visualization
109 //              c.) QA
110 //         NOTICE: Should be not used for the physic analysis
111 //         Function:
112 //            const AliExternalTrackParam *GetOuterParam() const { return fOp;}
113 //
114 //-----------------------------------------------------------------
115
116 #include <TMath.h>
117 #include <TParticle.h>
118 #include <TDatabasePDG.h>
119 #include <TMatrixD.h>
120
121 #include "AliESDVertex.h"
122 #include "AliESDtrack.h"
123 #include "AliESDEvent.h"
124 #include "AliKalmanTrack.h"
125 #include "AliVTrack.h"
126 #include "AliLog.h"
127 #include "AliTrackPointArray.h"
128 #include "TPolyMarker3D.h"
129 #include "AliTrackerBase.h"
130 #include "AliTPCdEdxInfo.h"
131
132 ClassImp(AliESDtrack)
133
134 void SetPIDValues(Double_t * dest, const Double_t * src, Int_t n) {
135   // This function copies "n" PID weights from "scr" to "dest"
136   // and normalizes their sum to 1 thus producing conditional probabilities.
137   // The negative weights are set to 0.
138   // In case all the weights are non-positive they are replaced by
139   // uniform probabilities
140
141   if (n<=0) return;
142
143   Float_t uniform = 1./(Float_t)n;
144
145   Float_t sum = 0;
146   for (Int_t i=0; i<n; i++) 
147     if (src[i]>=0) {
148       sum+=src[i];
149       dest[i] = src[i];
150     }
151     else {
152       dest[i] = 0;
153     }
154
155   if(sum>0)
156     for (Int_t i=0; i<n; i++) dest[i] /= sum;
157   else
158     for (Int_t i=0; i<n; i++) dest[i] = uniform;
159 }
160
161 //_______________________________________________________________________
162 AliESDtrack::AliESDtrack() : 
163   AliExternalTrackParam(),
164   fCp(0),
165   fIp(0),
166   fTPCInner(0),
167   fOp(0),
168   fHMPIDp(0),  
169   fFriendTrack(NULL),
170   fTPCFitMap(159),//number of padrows
171   fTPCClusterMap(159),//number of padrows
172   fTPCSharedMap(159),//number of padrows
173   fFlags(0),
174   fID(0),
175   fLabel(0),
176   fITSLabel(0),
177   fTPCLabel(0),
178   fTRDLabel(0),
179   fTOFCalChannel(-1),
180   fTOFindex(-1),
181   fHMPIDqn(0),
182   fHMPIDcluIdx(-1),
183   fCaloIndex(kEMCALNoMatch),
184   fHMPIDtrkTheta(0),
185   fHMPIDtrkPhi(0),
186   fHMPIDsignal(0),
187   fTrackLength(0),
188   fdTPC(0),fzTPC(0),
189   fCddTPC(0),fCdzTPC(0),fCzzTPC(0),
190   fCchi2TPC(0),
191   fD(0),fZ(0),
192   fCdd(0),fCdz(0),fCzz(0),
193   fCchi2(0),
194   fITSchi2(0),
195   fTPCchi2(0),
196   fTPCchi2Iter1(0),
197   fTRDchi2(0),
198   fTOFchi2(0),
199   fHMPIDchi2(0),
200   fGlobalChi2(0),
201   fITSsignal(0),
202   fTPCsignal(0),
203   fTPCsignalTuned(0),
204   fTPCsignalS(0),
205   fTPCdEdxInfo(0),
206   fTRDsignal(0),
207   fTRDQuality(0),
208   fTRDBudget(0),
209   fTOFsignal(99999),
210   fTOFsignalToT(99999),
211   fTOFsignalRaw(99999),
212   fTOFsignalDz(999),
213   fTOFsignalDx(999),
214   fTOFdeltaBC(999),
215   fTOFl0l1(999),
216   fCaloDx(0),
217   fCaloDz(0),
218   fHMPIDtrkX(0),
219   fHMPIDtrkY(0),
220   fHMPIDmipX(0),
221   fHMPIDmipY(0),
222   fTPCncls(0),
223   fTPCnclsF(0),
224   fTPCsignalN(0),
225   fTPCnclsIter1(0),
226   fTPCnclsFIter1(0),
227   fITSncls(0),
228   fITSClusterMap(0),
229   fITSSharedMap(0),
230   fTRDncls(0),
231   fTRDncls0(0),
232   fTRDntracklets(0),
233   fTRDNchamberdEdx(0),
234   fTRDNclusterdEdx(0),
235   fTRDnSlices(0),
236   fTRDslices(0x0),
237   fVertexID(-2),// -2 means an orphan track 
238   fESDEvent(0),
239   fCacheNCrossedRows(-10),
240   fCacheChi2TPCConstrainedVsGlobal(-10),
241   fCacheChi2TPCConstrainedVsGlobalVertex(0),
242   fTrackPhiOnEMCal(-999),
243   fTrackEtaOnEMCal(-999)
244 {
245   //
246   // The default ESD constructor 
247   //
248   if (!OnlineMode()) fFriendTrack=new AliESDfriendTrack();
249
250   Int_t i;
251   for (i=kNITSchi2Std;i--;) fITSchi2Std[i] = 0;
252   for (i=0; i<AliPID::kSPECIES; i++) {
253     fTrackTime[i]=0.;
254     fR[i]=0.;
255     fITSr[i]=0.;
256     fTPCr[i]=0.;
257     fTRDr[i]=0.;
258     fTOFr[i]=0.;
259     fHMPIDr[i]=0.;
260   }
261   
262   for (i=0; i<3; i++)   { fKinkIndexes[i]=0;}
263   for (i=0; i<3; i++)   { fV0Indexes[i]=0;}
264   for (i=0;i<kTRDnPlanes;i++) {
265     fTRDTimBin[i]=0;
266   }
267   for (i=0;i<4;i++) {fITSdEdxSamples[i]=0.;}
268   for (i=0;i<4;i++) {fTPCPoints[i]=0;}
269   for (i=0;i<3;i++) {fTOFLabel[i]=-1;}
270   for (i=0;i<10;i++) {fTOFInfo[i]=0;}
271   for (i=0;i<12;i++) {fITSModule[i]=-1;}
272 }
273
274 bool AliESDtrack::fgkOnlineMode=false;
275
276 //_______________________________________________________________________
277 AliESDtrack::AliESDtrack(const AliESDtrack& track):
278   AliExternalTrackParam(track),
279   fCp(0),
280   fIp(0),
281   fTPCInner(0),
282   fOp(0),
283   fHMPIDp(0),  
284   fFriendTrack(0),
285   fTPCFitMap(track.fTPCFitMap),
286   fTPCClusterMap(track.fTPCClusterMap),
287   fTPCSharedMap(track.fTPCSharedMap),
288   fFlags(track.fFlags),
289   fID(track.fID),
290   fLabel(track.fLabel),
291   fITSLabel(track.fITSLabel),
292   fTPCLabel(track.fTPCLabel),
293   fTRDLabel(track.fTRDLabel),
294   fTOFCalChannel(track.fTOFCalChannel),
295   fTOFindex(track.fTOFindex),
296   fHMPIDqn(track.fHMPIDqn),
297   fHMPIDcluIdx(track.fHMPIDcluIdx),
298   fCaloIndex(track.fCaloIndex),
299   fHMPIDtrkTheta(track.fHMPIDtrkTheta),
300   fHMPIDtrkPhi(track.fHMPIDtrkPhi),
301   fHMPIDsignal(track.fHMPIDsignal),
302   fTrackLength(track.fTrackLength),
303   fdTPC(track.fdTPC),fzTPC(track.fzTPC),
304   fCddTPC(track.fCddTPC),fCdzTPC(track.fCdzTPC),fCzzTPC(track.fCzzTPC),
305   fCchi2TPC(track.fCchi2TPC),
306   fD(track.fD),fZ(track.fZ),
307   fCdd(track.fCdd),fCdz(track.fCdz),fCzz(track.fCzz),
308   fCchi2(track.fCchi2),
309   fITSchi2(track.fITSchi2),
310   fTPCchi2(track.fTPCchi2),
311   fTPCchi2Iter1(track.fTPCchi2Iter1),
312   fTRDchi2(track.fTRDchi2),
313   fTOFchi2(track.fTOFchi2),
314   fHMPIDchi2(track.fHMPIDchi2),
315   fGlobalChi2(track.fGlobalChi2),
316   fITSsignal(track.fITSsignal),
317   fTPCsignal(track.fTPCsignal),
318   fTPCsignalTuned(track.fTPCsignalTuned),
319   fTPCsignalS(track.fTPCsignalS),
320   fTPCdEdxInfo(0),
321   fTRDsignal(track.fTRDsignal),
322   fTRDQuality(track.fTRDQuality),
323   fTRDBudget(track.fTRDBudget),
324   fTOFsignal(track.fTOFsignal),
325   fTOFsignalToT(track.fTOFsignalToT),
326   fTOFsignalRaw(track.fTOFsignalRaw),
327   fTOFsignalDz(track.fTOFsignalDz),
328   fTOFsignalDx(track.fTOFsignalDx),
329   fTOFdeltaBC(track.fTOFdeltaBC),
330   fTOFl0l1(track.fTOFl0l1),
331   fCaloDx(track.fCaloDx),
332   fCaloDz(track.fCaloDz),
333   fHMPIDtrkX(track.fHMPIDtrkX),
334   fHMPIDtrkY(track.fHMPIDtrkY),
335   fHMPIDmipX(track.fHMPIDmipX),
336   fHMPIDmipY(track.fHMPIDmipY),
337   fTPCncls(track.fTPCncls),
338   fTPCnclsF(track.fTPCnclsF),
339   fTPCsignalN(track.fTPCsignalN),
340   fTPCnclsIter1(track.fTPCnclsIter1),
341   fTPCnclsFIter1(track.fTPCnclsIter1),
342   fITSncls(track.fITSncls),
343   fITSClusterMap(track.fITSClusterMap),
344   fITSSharedMap(track.fITSSharedMap),
345   fTRDncls(track.fTRDncls),
346   fTRDncls0(track.fTRDncls0),
347   fTRDntracklets(track.fTRDntracklets),
348   fTRDNchamberdEdx(track.fTRDNchamberdEdx),
349   fTRDNclusterdEdx(track.fTRDNclusterdEdx),
350   fTRDnSlices(track.fTRDnSlices),
351   fTRDslices(0x0),
352   fVertexID(track.fVertexID),
353   fESDEvent(track.fESDEvent),
354   fCacheNCrossedRows(track.fCacheNCrossedRows),
355   fCacheChi2TPCConstrainedVsGlobal(track.fCacheChi2TPCConstrainedVsGlobal),
356   fCacheChi2TPCConstrainedVsGlobalVertex(track.fCacheChi2TPCConstrainedVsGlobalVertex),
357   fTrackPhiOnEMCal(track.fTrackPhiOnEMCal),
358   fTrackEtaOnEMCal(track.fTrackEtaOnEMCal)
359 {
360   //
361   //copy constructor
362   //
363   for (Int_t i=kNITSchi2Std;i--;) fITSchi2Std[i] = track.fTrackTime[i];
364   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTrackTime[i]=track.fTrackTime[i];
365   for (Int_t i=0;i<AliPID::kSPECIES;i++)  fR[i]=track.fR[i];
366   //
367   for (Int_t i=0;i<AliPID::kSPECIES;i++) fITSr[i]=track.fITSr[i]; 
368   //
369   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTPCr[i]=track.fTPCr[i]; 
370   for (Int_t i=0;i<4;i++) {fITSdEdxSamples[i]=track.fITSdEdxSamples[i];}
371   for (Int_t i=0;i<4;i++) {fTPCPoints[i]=track.fTPCPoints[i];}
372   for (Int_t i=0; i<3;i++)   { fKinkIndexes[i]=track.fKinkIndexes[i];}
373   for (Int_t i=0; i<3;i++)   { fV0Indexes[i]=track.fV0Indexes[i];}
374   //
375   for (Int_t i=0;i<kTRDnPlanes;i++) {
376     fTRDTimBin[i]=track.fTRDTimBin[i];
377   }
378
379   if (fTRDnSlices) {
380     fTRDslices=new Double32_t[fTRDnSlices];
381     for (Int_t i=0; i<fTRDnSlices; i++) fTRDslices[i]=track.fTRDslices[i];
382   }
383
384
385   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTRDr[i]=track.fTRDr[i]; 
386   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTOFr[i]=track.fTOFr[i];
387   for (Int_t i=0;i<3;i++) fTOFLabel[i]=track.fTOFLabel[i];
388   for (Int_t i=0;i<10;i++) fTOFInfo[i]=track.fTOFInfo[i];
389   for (Int_t i=0;i<12;i++) fITSModule[i]=track.fITSModule[i];
390   for (Int_t i=0;i<AliPID::kSPECIES;i++) fHMPIDr[i]=track.fHMPIDr[i];
391
392   if (track.fCp) fCp=new AliExternalTrackParam(*track.fCp);
393   if (track.fIp) fIp=new AliExternalTrackParam(*track.fIp);
394   if (track.fTPCInner) fTPCInner=new AliExternalTrackParam(*track.fTPCInner);
395   if (track.fOp) fOp=new AliExternalTrackParam(*track.fOp);
396   if (track.fHMPIDp) fHMPIDp=new AliExternalTrackParam(*track.fHMPIDp);
397   if (track.fTPCdEdxInfo) fTPCdEdxInfo = new AliTPCdEdxInfo(*track.fTPCdEdxInfo);
398
399   
400   if (track.fFriendTrack) fFriendTrack=new AliESDfriendTrack(*(track.fFriendTrack));
401 }
402
403 //_______________________________________________________________________
404 AliESDtrack::AliESDtrack(const AliVTrack *track) : 
405   AliExternalTrackParam(track),
406   fCp(0),
407   fIp(0),
408   fTPCInner(0),
409   fOp(0),
410   fHMPIDp(0),  
411   fFriendTrack(0),
412   fTPCFitMap(159),//number of padrows
413   fTPCClusterMap(159),//number of padrows
414   fTPCSharedMap(159),//number of padrows
415   fFlags(0),
416   fID(),
417   fLabel(0),
418   fITSLabel(0),
419   fTPCLabel(0),
420   fTRDLabel(0),
421   fTOFCalChannel(-1),
422   fTOFindex(-1),
423   fHMPIDqn(0),
424   fHMPIDcluIdx(-1),
425   fCaloIndex(kEMCALNoMatch),
426   fHMPIDtrkTheta(0),
427   fHMPIDtrkPhi(0),
428   fHMPIDsignal(0),
429   fTrackLength(0),
430   fdTPC(0),fzTPC(0),
431   fCddTPC(0),fCdzTPC(0),fCzzTPC(0),
432   fCchi2TPC(0),
433   fD(0),fZ(0),
434   fCdd(0),fCdz(0),fCzz(0),
435   fCchi2(0),
436   fITSchi2(0),
437   fTPCchi2(0),
438   fTPCchi2Iter1(0),
439   fTRDchi2(0),
440   fTOFchi2(0),
441   fHMPIDchi2(0),
442   fGlobalChi2(0),
443   fITSsignal(0),
444   fTPCsignal(0),
445   fTPCsignalTuned(0),
446   fTPCsignalS(0),
447   fTPCdEdxInfo(0),
448   fTRDsignal(0),
449   fTRDQuality(0),
450   fTRDBudget(0),
451   fTOFsignal(99999),
452   fTOFsignalToT(99999),
453   fTOFsignalRaw(99999),
454   fTOFsignalDz(999),
455   fTOFsignalDx(999),
456   fTOFdeltaBC(999),
457   fTOFl0l1(999),
458   fCaloDx(0),
459   fCaloDz(0),
460   fHMPIDtrkX(0),
461   fHMPIDtrkY(0),
462   fHMPIDmipX(0),
463   fHMPIDmipY(0),
464   fTPCncls(0),
465   fTPCnclsF(0),
466   fTPCsignalN(0),
467   fTPCnclsIter1(0),
468   fTPCnclsFIter1(0),
469   fITSncls(0),
470   fITSClusterMap(0),
471   fITSSharedMap(0),
472   fTRDncls(0),
473   fTRDncls0(0),
474   fTRDntracklets(0),
475   fTRDNchamberdEdx(0),
476   fTRDNclusterdEdx(0),
477   fTRDnSlices(0),
478   fTRDslices(0x0),
479   fVertexID(-2),  // -2 means an orphan track
480   fESDEvent(0),
481   fCacheNCrossedRows(-10),
482   fCacheChi2TPCConstrainedVsGlobal(-10),
483   fCacheChi2TPCConstrainedVsGlobalVertex(0),
484   fTrackPhiOnEMCal(-999),
485   fTrackEtaOnEMCal(-999)
486 {
487   //
488   // ESD track from AliVTrack.
489   // This is not a copy constructor !
490   //
491
492   if (track->InheritsFrom("AliExternalTrackParam")) {
493      AliError("This is not a copy constructor. Use AliESDtrack(const AliESDtrack &) !");
494      AliWarning("Calling the default constructor...");
495      AliESDtrack();
496      return;
497   }
498
499   // Reset all the arrays
500   Int_t i;
501   for (i=kNITSchi2Std;i--;) fITSchi2Std[i] = 0;
502   for (i=0; i<AliPID::kSPECIES; i++) {
503     fTrackTime[i]=0.;
504     fR[i]=0.;
505     fITSr[i]=0.;
506     fTPCr[i]=0.;
507     fTRDr[i]=0.;
508     fTOFr[i]=0.;
509     fHMPIDr[i]=0.;
510   }
511   
512   for (i=0; i<3; i++)   { fKinkIndexes[i]=0;}
513   for (i=0; i<3; i++)   { fV0Indexes[i]=-1;}
514   for (i=0;i<kTRDnPlanes;i++) {
515     fTRDTimBin[i]=0;
516   }
517   for (i=0;i<4;i++) {fITSdEdxSamples[i]=0.;}
518   for (i=0;i<4;i++) {fTPCPoints[i]=0;}
519   for (i=0;i<3;i++) {fTOFLabel[i]=-1;}
520   for (i=0;i<10;i++) {fTOFInfo[i]=0;}
521   for (i=0;i<12;i++) {fITSModule[i]=-1;}
522
523   // Set the ID
524   SetID(track->GetID());
525
526   // Set ITS cluster map
527   fITSClusterMap=track->GetITSClusterMap();
528   fITSSharedMap=0;
529
530   fITSncls=0;
531   for(i=0; i<6; i++) {
532     if(HasPointOnITSLayer(i)) fITSncls++;
533   }
534
535   // Set TPC ncls 
536   fTPCncls=track->GetTPCNcls();
537
538
539   // Set the combined PID
540   const Double_t *pid = track->PID();
541   if(pid){
542     for (i=0; i<AliPID::kSPECIES; i++) fR[i]=pid[i];
543   }
544   // AliESD track label
545   SetLabel(track->GetLabel());
546   // Set the status
547   SetStatus(track->GetStatus());
548 }
549
550 //_______________________________________________________________________
551 AliESDtrack::AliESDtrack(TParticle * part) : 
552   AliExternalTrackParam(),
553   fCp(0),
554   fIp(0),
555   fTPCInner(0),
556   fOp(0),
557   fHMPIDp(0),  
558   fFriendTrack(0),
559   fTPCFitMap(159),//number of padrows
560   fTPCClusterMap(159),//number of padrows
561   fTPCSharedMap(159),//number of padrows
562   fFlags(0),
563   fID(0),
564   fLabel(0),
565   fITSLabel(0),
566   fTPCLabel(0),
567   fTRDLabel(0),
568   fTOFCalChannel(-1),
569   fTOFindex(-1),
570   fHMPIDqn(0),
571   fHMPIDcluIdx(-1),
572   fCaloIndex(kEMCALNoMatch),
573   fHMPIDtrkTheta(0),
574   fHMPIDtrkPhi(0),
575   fHMPIDsignal(0),
576   fTrackLength(0),
577   fdTPC(0),fzTPC(0),
578   fCddTPC(0),fCdzTPC(0),fCzzTPC(0),
579   fCchi2TPC(0),
580   fD(0),fZ(0),
581   fCdd(0),fCdz(0),fCzz(0),
582   fCchi2(0),
583   fITSchi2(0),
584   fTPCchi2(0),
585   fTPCchi2Iter1(0),  
586   fTRDchi2(0),
587   fTOFchi2(0),
588   fHMPIDchi2(0),
589   fGlobalChi2(0),
590   fITSsignal(0),
591   fTPCsignal(0),
592   fTPCsignalTuned(0),
593   fTPCsignalS(0),
594   fTPCdEdxInfo(0),
595   fTRDsignal(0),
596   fTRDQuality(0),
597   fTRDBudget(0),
598   fTOFsignal(99999),
599   fTOFsignalToT(99999),
600   fTOFsignalRaw(99999),
601   fTOFsignalDz(999),
602   fTOFsignalDx(999),
603   fTOFdeltaBC(999),
604   fTOFl0l1(999),
605   fCaloDx(0),
606   fCaloDz(0),
607   fHMPIDtrkX(0),
608   fHMPIDtrkY(0),
609   fHMPIDmipX(0),
610   fHMPIDmipY(0),
611   fTPCncls(0),
612   fTPCnclsF(0),
613   fTPCsignalN(0),
614   fTPCnclsIter1(0),
615   fTPCnclsFIter1(0),
616   fITSncls(0),
617   fITSClusterMap(0),
618   fITSSharedMap(0),
619   fTRDncls(0),
620   fTRDncls0(0),
621   fTRDntracklets(0),
622   fTRDNchamberdEdx(0),
623   fTRDNclusterdEdx(0),
624   fTRDnSlices(0),
625   fTRDslices(0x0),
626   fVertexID(-2),  // -2 means an orphan track
627   fESDEvent(0),
628   fCacheNCrossedRows(-10),
629   fCacheChi2TPCConstrainedVsGlobal(-10),
630   fCacheChi2TPCConstrainedVsGlobalVertex(0),
631   fTrackPhiOnEMCal(-999),
632   fTrackEtaOnEMCal(-999)
633 {
634   //
635   // ESD track from TParticle
636   //
637
638   // Reset all the arrays
639   Int_t i;
640   for (i=kNITSchi2Std;i--;) fITSchi2Std[i] = 0;
641   for (i=0; i<AliPID::kSPECIES; i++) {
642     fTrackTime[i]=0.;
643     fR[i]=0.;
644     fITSr[i]=0.;
645     fTPCr[i]=0.;
646     fTRDr[i]=0.;
647     fTOFr[i]=0.;
648     fHMPIDr[i]=0.;
649   }
650   
651   for (i=0; i<3; i++)   { fKinkIndexes[i]=0;}
652   for (i=0; i<3; i++)   { fV0Indexes[i]=-1;}
653   for (i=0;i<kTRDnPlanes;i++) {
654     fTRDTimBin[i]=0;
655   }
656   for (i=0;i<4;i++) {fITSdEdxSamples[i]=0.;}
657   for (i=0;i<4;i++) {fTPCPoints[i]=0;}
658   for (i=0;i<3;i++) {fTOFLabel[i]=-1;}
659   for (i=0;i<10;i++) {fTOFInfo[i]=0;}
660   for (i=0;i<12;i++) {fITSModule[i]=-1;}
661
662   // Calculate the AliExternalTrackParam content
663
664   Double_t xref;
665   Double_t alpha;
666   Double_t param[5];
667   Double_t covar[15];
668
669   // Calculate alpha: the rotation angle of the corresponding local system (TPC sector)
670   alpha = part->Phi()*180./TMath::Pi();
671   if (alpha<0) alpha+= 360.;
672   if (alpha>360) alpha -= 360.;
673
674   Int_t sector = (Int_t)(alpha/20.);
675   alpha = 10. + 20.*sector;
676   alpha /= 180;
677   alpha *= TMath::Pi();
678
679   // Covariance matrix: no errors, the parameters are exact
680   for (i=0; i<15; i++) covar[i]=0.;
681
682   // Get the vertex of origin and the momentum
683   TVector3 ver(part->Vx(),part->Vy(),part->Vz());
684   TVector3 mom(part->Px(),part->Py(),part->Pz());
685
686   // Rotate to the local coordinate system (TPC sector)
687   ver.RotateZ(-alpha);
688   mom.RotateZ(-alpha);
689
690   // X of the referense plane
691   xref = ver.X();
692
693   Int_t pdgCode = part->GetPdgCode();
694
695   Double_t charge = 
696     TDatabasePDG::Instance()->GetParticle(pdgCode)->Charge();
697
698   param[0] = ver.Y();
699   param[1] = ver.Z();
700   param[2] = TMath::Sin(mom.Phi());
701   param[3] = mom.Pz()/mom.Pt();
702   param[4] = TMath::Sign(1/mom.Pt(),charge);
703
704   // Set AliExternalTrackParam
705   Set(xref, alpha, param, covar);
706
707   // Set the PID
708   Int_t indexPID = 99;
709
710   switch (TMath::Abs(pdgCode)) {
711
712   case  11: // electron
713     indexPID = 0;
714     break;
715
716   case 13: // muon
717     indexPID = 1;
718     break;
719
720   case 211: // pion
721     indexPID = 2;
722     break;
723
724   case 321: // kaon
725     indexPID = 3;
726     break;
727
728   case 2212: // proton
729     indexPID = 4;
730     break;
731
732   default:
733     break;
734   }
735
736   // If the particle is not e,mu,pi,K or p the PID probabilities are set to 0
737   if (indexPID < AliPID::kSPECIES) {
738     fR[indexPID]=1.;
739     fITSr[indexPID]=1.;
740     fTPCr[indexPID]=1.;
741     fTRDr[indexPID]=1.;
742     fTOFr[indexPID]=1.;
743     fHMPIDr[indexPID]=1.;
744
745   }
746   // AliESD track label
747   SetLabel(part->GetUniqueID());
748
749 }
750
751 //_______________________________________________________________________
752 AliESDtrack::~AliESDtrack(){ 
753   //
754   // This is destructor according Coding Conventrions 
755   //
756   //printf("Delete track\n");
757   delete fIp; 
758   delete fTPCInner; 
759   delete fOp;
760   delete fHMPIDp;
761   delete fCp; 
762   delete fFriendTrack;
763   delete fTPCdEdxInfo;
764   if(fTRDnSlices)
765     delete[] fTRDslices;
766
767   //Reset cached values - needed for TClonesArray in AliESDInputHandler
768   fCacheNCrossedRows = -10.;
769   fCacheChi2TPCConstrainedVsGlobal = -10.;
770   if(fCacheChi2TPCConstrainedVsGlobalVertex) fCacheChi2TPCConstrainedVsGlobalVertex = 0;
771
772 }
773
774 AliESDtrack &AliESDtrack::operator=(const AliESDtrack &source){
775   
776
777   if(&source == this) return *this;
778   AliExternalTrackParam::operator=(source);
779
780   
781   if(source.fCp){
782     // we have the trackparam: assign or copy construct
783     if(fCp)*fCp = *source.fCp;
784     else fCp = new AliExternalTrackParam(*source.fCp);
785   }
786   else{
787     // no track param delete the old one
788     delete fCp;
789     fCp = 0;
790   }
791
792   if(source.fIp){
793     // we have the trackparam: assign or copy construct
794     if(fIp)*fIp = *source.fIp;
795     else fIp = new AliExternalTrackParam(*source.fIp);
796   }
797   else{
798     // no track param delete the old one
799     delete fIp;
800     fIp = 0;
801   }
802
803
804   if(source.fTPCInner){
805     // we have the trackparam: assign or copy construct
806     if(fTPCInner) *fTPCInner = *source.fTPCInner;
807     else fTPCInner = new AliExternalTrackParam(*source.fTPCInner);
808   }
809   else{
810     // no track param delete the old one
811     delete fTPCInner;
812     fTPCInner = 0;
813   }
814
815   if(source.fTPCdEdxInfo) {
816     if(fTPCdEdxInfo) *fTPCdEdxInfo = *source.fTPCdEdxInfo;
817     fTPCdEdxInfo = new AliTPCdEdxInfo(*source.fTPCdEdxInfo);
818   }
819
820   if(source.fOp){
821     // we have the trackparam: assign or copy construct
822     if(fOp) *fOp = *source.fOp;
823     else fOp = new AliExternalTrackParam(*source.fOp);
824   }
825   else{
826     // no track param delete the old one
827     delete fOp;
828     fOp = 0;
829   }
830
831   
832   if(source.fHMPIDp){
833     // we have the trackparam: assign or copy construct
834     if(fHMPIDp) *fHMPIDp = *source.fHMPIDp;
835     else fHMPIDp = new AliExternalTrackParam(*source.fHMPIDp);
836   }
837   else{
838     // no track param delete the old one
839     delete fHMPIDp;
840     fHMPIDp = 0;
841   }
842
843   // copy also the friend track 
844   // use copy constructor
845   if(source.fFriendTrack){
846     // we have the trackparam: assign or copy construct
847     delete fFriendTrack; fFriendTrack=new AliESDfriendTrack(*source.fFriendTrack);
848   }
849   else{
850     // no track param delete the old one
851     delete fFriendTrack; fFriendTrack= 0;
852   }
853
854   fTPCFitMap = source.fTPCFitMap; 
855   fTPCClusterMap = source.fTPCClusterMap; 
856   fTPCSharedMap  = source.fTPCSharedMap;  
857   // the simple stuff
858   fFlags    = source.fFlags; 
859   fID       = source.fID;             
860   fLabel    = source.fLabel;
861   fITSLabel = source.fITSLabel;
862   for(int i = 0; i< 12;++i){
863     fITSModule[i] = source.fITSModule[i];
864   }
865   fTPCLabel = source.fTPCLabel; 
866   fTRDLabel = source.fTRDLabel;
867   for(int i = 0; i< 3;++i){
868     fTOFLabel[i] = source.fTOFLabel[i];    
869   }
870   fTOFCalChannel = source.fTOFCalChannel;
871   fTOFindex      = source.fTOFindex;
872   fHMPIDqn       = source.fHMPIDqn;
873   fHMPIDcluIdx   = source.fHMPIDcluIdx; 
874   fCaloIndex    = source.fCaloIndex;
875   for (int i=kNITSchi2Std;i--;) fITSchi2Std[i] = source.fITSchi2Std[i];
876   for(int i = 0; i< 3;++i){
877     fKinkIndexes[i] = source.fKinkIndexes[i]; 
878     fV0Indexes[i]   = source.fV0Indexes[i]; 
879   }
880
881   for(int i = 0; i< AliPID::kSPECIES;++i){
882     fR[i]     = source.fR[i];
883     fITSr[i]  = source.fITSr[i];
884     fTPCr[i]  = source.fTPCr[i];
885     fTRDr[i]  = source.fTRDr[i];
886     fTOFr[i]  = source.fTOFr[i];
887     fHMPIDr[i] = source.fHMPIDr[i];
888     fTrackTime[i] = source.fTrackTime[i];  
889   }
890
891   fHMPIDtrkTheta = source.fHMPIDtrkTheta;
892   fHMPIDtrkPhi   = source.fHMPIDtrkPhi;
893   fHMPIDsignal   = source.fHMPIDsignal; 
894
895   
896   fTrackLength   = source. fTrackLength;
897   fdTPC  = source.fdTPC; 
898   fzTPC  = source.fzTPC; 
899   fCddTPC = source.fCddTPC;
900   fCdzTPC = source.fCdzTPC;
901   fCzzTPC = source.fCzzTPC;
902   fCchi2TPC = source.fCchi2TPC;
903
904   fD  = source.fD; 
905   fZ  = source.fZ; 
906   fCdd = source.fCdd;
907   fCdz = source.fCdz;
908   fCzz = source.fCzz;
909   fCchi2     = source.fCchi2;
910
911   fITSchi2   = source.fITSchi2;             
912   fTPCchi2   = source.fTPCchi2;            
913   fTPCchi2Iter1   = source.fTPCchi2Iter1;            
914   fTRDchi2   = source.fTRDchi2;      
915   fTOFchi2   = source.fTOFchi2;      
916   fHMPIDchi2 = source.fHMPIDchi2;      
917
918   fGlobalChi2 = source.fGlobalChi2;      
919
920   fITSsignal  = source.fITSsignal;     
921   for (Int_t i=0;i<4;i++) {fITSdEdxSamples[i]=source.fITSdEdxSamples[i];}
922   fTPCsignal  = source.fTPCsignal;     
923   fTPCsignalTuned  = source.fTPCsignalTuned;
924   fTPCsignalS = source.fTPCsignalS;    
925   for(int i = 0; i< 4;++i){
926     fTPCPoints[i] = source.fTPCPoints[i];  
927   }
928   fTRDsignal = source.fTRDsignal;
929   fTRDNchamberdEdx = source.fTRDNchamberdEdx;
930   fTRDNclusterdEdx = source.fTRDNclusterdEdx;
931
932   for(int i = 0;i < kTRDnPlanes;++i){
933     fTRDTimBin[i] = source.fTRDTimBin[i];   
934   }
935
936   if(fTRDnSlices)
937     delete[] fTRDslices;
938   fTRDslices=0;
939   fTRDnSlices=source.fTRDnSlices;
940   if (fTRDnSlices) {
941     fTRDslices=new Double32_t[fTRDnSlices];
942     for(int j = 0;j < fTRDnSlices;++j) fTRDslices[j] = source.fTRDslices[j];
943   }
944
945   fTRDQuality =   source.fTRDQuality;     
946   fTRDBudget  =   source.fTRDBudget;      
947   fTOFsignal  =   source.fTOFsignal;     
948   fTOFsignalToT = source.fTOFsignalToT;   
949   fTOFsignalRaw = source.fTOFsignalRaw;  
950   fTOFsignalDz  = source.fTOFsignalDz;      
951   fTOFsignalDx  = source.fTOFsignalDx;      
952   fTOFdeltaBC   = source.fTOFdeltaBC;
953   fTOFl0l1      = source.fTOFl0l1;
954  
955   for(int i = 0;i<10;++i){
956     fTOFInfo[i] = source.fTOFInfo[i];    
957   }
958
959   fHMPIDtrkX = source.fHMPIDtrkX; 
960   fHMPIDtrkY = source.fHMPIDtrkY; 
961   fHMPIDmipX = source.fHMPIDmipX;
962   fHMPIDmipY = source.fHMPIDmipY; 
963
964   fTPCncls    = source.fTPCncls;      
965   fTPCnclsF   = source.fTPCnclsF;     
966   fTPCsignalN = source.fTPCsignalN;   
967   fTPCnclsIter1    = source.fTPCnclsIter1;      
968   fTPCnclsFIter1   = source.fTPCnclsFIter1;     
969
970   fITSncls = source.fITSncls;       
971   fITSClusterMap = source.fITSClusterMap; 
972   fITSSharedMap = source.fITSSharedMap; 
973   fTRDncls   = source.fTRDncls;       
974   fTRDncls0  = source.fTRDncls0;      
975   fTRDntracklets  = source.fTRDntracklets; 
976   fVertexID = source.fVertexID;
977
978   fCacheNCrossedRows = source.fCacheNCrossedRows;
979   fCacheChi2TPCConstrainedVsGlobal = source.fCacheChi2TPCConstrainedVsGlobal;
980   fCacheChi2TPCConstrainedVsGlobalVertex = source.fCacheChi2TPCConstrainedVsGlobalVertex;
981
982   fTrackPhiOnEMCal= source.fTrackPhiOnEMCal;
983   fTrackEtaOnEMCal= source.fTrackEtaOnEMCal;
984
985   return *this;
986 }
987
988
989
990 void AliESDtrack::Copy(TObject &obj) const {
991   
992   // this overwrites the virtual TOBject::Copy()
993   // to allow run time copying without casting
994   // in AliESDEvent
995
996   if(this==&obj)return;
997   AliESDtrack *robj = dynamic_cast<AliESDtrack*>(&obj);
998   if(!robj)return; // not an AliESDtrack
999   *robj = *this;
1000
1001 }
1002
1003
1004
1005 void AliESDtrack::AddCalibObject(TObject * object){
1006   //
1007   // add calib object to the list
1008   //
1009   if (!fFriendTrack) fFriendTrack  = new AliESDfriendTrack;
1010   if (!fFriendTrack) return;
1011   fFriendTrack->AddCalibObject(object);
1012 }
1013
1014 TObject *  AliESDtrack::GetCalibObject(Int_t index){
1015   //
1016   // return calib objct at given position
1017   //
1018   if (!fFriendTrack) return 0;
1019   return fFriendTrack->GetCalibObject(index);
1020 }
1021
1022
1023 Bool_t AliESDtrack::FillTPCOnlyTrack(AliESDtrack &track){
1024   
1025   // Fills the information of the TPC-only first reconstruction pass
1026   // into the passed ESDtrack object. For consistency fTPCInner is also filled
1027   // again
1028
1029
1030
1031   // For data produced before r26675
1032   // RelateToVertexTPC was not properly called during reco
1033   // so you'll have to call it again, before FillTPCOnlyTrack
1034   //  Float_t p[2],cov[3];
1035   // track->GetImpactParametersTPC(p,cov); 
1036   // if(p[0]==0&&p[1]==0) // <- Default values
1037   //  track->RelateToVertexTPC(esd->GetPrimaryVertexTPC(),esd->GetMagneticField(),kVeryBig);
1038   
1039
1040   if(!fTPCInner)return kFALSE;
1041
1042   // fill the TPC track params to the global track parameters
1043   track.Set(fTPCInner->GetX(),fTPCInner->GetAlpha(),fTPCInner->GetParameter(),fTPCInner->GetCovariance());
1044   track.fD = fdTPC;
1045   track.fZ = fzTPC;
1046   track.fCdd = fCddTPC;
1047   track.fCdz = fCdzTPC;
1048   track.fCzz = fCzzTPC;
1049
1050   // copy the inner params
1051   if(track.fIp) *track.fIp = *fIp;
1052   else track.fIp = new AliExternalTrackParam(*fIp);
1053
1054   // copy the TPCinner parameters
1055   if(track.fTPCInner) *track.fTPCInner = *fTPCInner;
1056   else track.fTPCInner = new AliExternalTrackParam(*fTPCInner);
1057   track.fdTPC   = fdTPC;
1058   track.fzTPC   = fzTPC;
1059   track.fCddTPC = fCddTPC;
1060   track.fCdzTPC = fCdzTPC;
1061   track.fCzzTPC = fCzzTPC;
1062   track.fCchi2TPC = fCchi2TPC;
1063
1064   // copy all other TPC specific parameters
1065
1066   // replace label by TPC label
1067   track.fLabel    = fTPCLabel;
1068   track.fTPCLabel = fTPCLabel;
1069
1070   track.fTPCchi2 = fTPCchi2; 
1071   track.fTPCchi2Iter1 = fTPCchi2Iter1; 
1072   track.fTPCsignal = fTPCsignal;
1073   track.fTPCsignalTuned = fTPCsignalTuned;
1074   track.fTPCsignalS = fTPCsignalS;
1075   for(int i = 0;i<4;++i)track.fTPCPoints[i] = fTPCPoints[i];
1076
1077   track.fTPCncls    = fTPCncls;     
1078   track.fTPCnclsF   = fTPCnclsF;     
1079   track.fTPCsignalN =  fTPCsignalN;
1080   track.fTPCnclsIter1    = fTPCnclsIter1;     
1081   track.fTPCnclsFIter1   = fTPCnclsFIter1;     
1082
1083   // PID 
1084   for(int i=0;i<AliPID::kSPECIES;++i){
1085     track.fTPCr[i] = fTPCr[i];
1086     // combined PID is TPC only!
1087     track.fR[i] = fTPCr[i];
1088   }
1089   track.fTPCFitMap = fTPCFitMap;
1090   track.fTPCClusterMap = fTPCClusterMap;
1091   track.fTPCSharedMap = fTPCSharedMap;
1092
1093
1094   // reset the flags
1095   track.fFlags = kTPCin;
1096   track.fID    = fID;
1097
1098   track.fFlags |= fFlags & kTPCpid; //copy the TPCpid status flag
1099  
1100   for (Int_t i=0;i<3;i++) track.fKinkIndexes[i] = fKinkIndexes[i];
1101   
1102   return kTRUE;
1103     
1104 }
1105
1106 //_______________________________________________________________________
1107 void AliESDtrack::MakeMiniESDtrack(){
1108   // Resets everything except
1109   // fFlags: Reconstruction status flags 
1110   // fLabel: Track label
1111   // fID:  Unique ID of the track
1112   // Impact parameter information
1113   // fR[AliPID::kSPECIES]: combined "detector response probability"
1114   // Running track parameters in the base class (AliExternalTrackParam)
1115   
1116   fTrackLength = 0;
1117
1118   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTrackTime[i] = 0;
1119
1120   // Reset track parameters constrained to the primary vertex
1121   delete fCp;fCp = 0;
1122
1123   // Reset track parameters at the inner wall of TPC
1124   delete fIp;fIp = 0;
1125   delete fTPCInner;fTPCInner=0;
1126   // Reset track parameters at the inner wall of the TRD
1127   delete fOp;fOp = 0;
1128   // Reset track parameters at the HMPID
1129   delete fHMPIDp;fHMPIDp = 0;
1130
1131
1132   // Reset ITS track related information
1133   fITSchi2 = 0;
1134   fITSncls = 0;       
1135   fITSClusterMap=0;
1136   fITSSharedMap=0;
1137   fITSsignal = 0;     
1138   for (Int_t i=0;i<4;i++) fITSdEdxSamples[i] = 0.;
1139   for (Int_t i=0;i<AliPID::kSPECIES;i++) fITSr[i]=0; 
1140   fITSLabel = 0;       
1141
1142   // Reset TPC related track information
1143   fTPCchi2 = 0;       
1144   fTPCchi2Iter1 = 0;       
1145   fTPCncls = 0;       
1146   fTPCnclsF = 0;       
1147   fTPCnclsIter1 = 0;       
1148   fTPCnclsFIter1 = 0;  
1149   fTPCFitMap = 0;       
1150   fTPCClusterMap = 0;  
1151   fTPCSharedMap = 0;  
1152   fTPCsignal= 0;      
1153   fTPCsignalTuned= 0;
1154   fTPCsignalS= 0;      
1155   fTPCsignalN= 0;      
1156   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTPCr[i]=0; 
1157   fTPCLabel=0;       
1158   for (Int_t i=0;i<4;i++) fTPCPoints[i] = 0;
1159   for (Int_t i=0; i<3;i++)   fKinkIndexes[i] = 0;
1160   for (Int_t i=0; i<3;i++)   fV0Indexes[i] = 0;
1161
1162   // Reset TRD related track information
1163   fTRDchi2 = 0;        
1164   fTRDncls = 0;       
1165   fTRDncls0 = 0;       
1166   fTRDsignal = 0;      
1167   fTRDNchamberdEdx = 0;
1168   fTRDNclusterdEdx = 0;
1169
1170   for (Int_t i=0;i<kTRDnPlanes;i++) {
1171     fTRDTimBin[i]  = 0;
1172   }
1173   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTRDr[i] = 0; 
1174   fTRDLabel = 0;       
1175   fTRDQuality  = 0;
1176   fTRDntracklets = 0;
1177   if(fTRDnSlices)
1178     delete[] fTRDslices;
1179   fTRDslices=0x0;
1180   fTRDnSlices=0;
1181   fTRDBudget  = 0;
1182
1183   // Reset TOF related track information
1184   fTOFchi2 = 0;        
1185   fTOFindex = -1;       
1186   fTOFsignal = 99999;      
1187   fTOFCalChannel = -1;
1188   fTOFsignalToT = 99999;
1189   fTOFsignalRaw = 99999;
1190   fTOFsignalDz = 999;
1191   fTOFsignalDx = 999;
1192   fTOFdeltaBC = 999;
1193   fTOFl0l1 = 999;
1194   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTOFr[i] = 0;
1195   for (Int_t i=0;i<3;i++) fTOFLabel[i] = -1;
1196   for (Int_t i=0;i<10;i++) fTOFInfo[i] = 0;
1197
1198   // Reset HMPID related track information
1199   fHMPIDchi2 = 0;     
1200   fHMPIDqn = 0;     
1201   fHMPIDcluIdx = -1;     
1202   fHMPIDsignal = 0;     
1203   for (Int_t i=0;i<AliPID::kSPECIES;i++) fHMPIDr[i] = 0;
1204   fHMPIDtrkTheta = 0;     
1205   fHMPIDtrkPhi = 0;      
1206   fHMPIDtrkX = 0;     
1207   fHMPIDtrkY = 0;      
1208   fHMPIDmipX = 0;
1209   fHMPIDmipY = 0;
1210   fCaloIndex = kEMCALNoMatch;
1211
1212   // reset global track chi2
1213   fGlobalChi2 = 0;
1214
1215   fVertexID = -2; // an orphan track
1216
1217   delete fFriendTrack; fFriendTrack = 0;
1218
1219
1220 //_______________________________________________________________________
1221 Int_t AliESDtrack::GetPID(Bool_t tpcOnly) const 
1222 {
1223   // Returns the particle most probable id
1224   Int_t i;
1225   const Double32_t *prob = 0;
1226   if (tpcOnly) { // check if TPCpid is valid
1227     prob = fTPCr;
1228     for (i=0; i<AliPID::kSPECIES-1; i++) if (prob[i] != prob[i+1]) break;
1229     if (i == AliPID::kSPECIES-1) prob = 0; // not valid, try with combined pid
1230   }
1231   if (!prob) { // either requested TPCpid is not valid or comb.pid is requested 
1232     prob = fR;
1233     for (i=0; i<AliPID::kSPECIES-1; i++) if (prob[i] != prob[i+1]) break;
1234     if (i == AliPID::kSPECIES-1) return AliPID::kPion;  // If all the probabilities are equal, return the pion mass
1235   }
1236   //
1237   Float_t max=0.;
1238   Int_t k=-1;
1239   for (i=0; i<AliPID::kSPECIES; i++) if (prob[i]>max) {k=i; max=prob[i];}
1240   //
1241   if (k==0) { // dE/dx "crossing points" in the TPC
1242     Double_t p=GetP();
1243     if ((p>0.38)&&(p<0.48))
1244       if (prob[0]<prob[3]*10.) return AliPID::kKaon;
1245     if ((p>0.75)&&(p<0.85))
1246       if (prob[0]<prob[4]*10.) return AliPID::kProton;
1247     return AliPID::kElectron;
1248   }
1249   if (k==1) return AliPID::kMuon; 
1250   if (k==2||k==-1) return AliPID::kPion;
1251   if (k==3) return AliPID::kKaon;
1252   if (k==4) return AliPID::kProton;
1253   AliWarning("Undefined PID !");
1254   return AliPID::kPion;
1255 }
1256
1257 //_______________________________________________________________________
1258 Int_t AliESDtrack::GetTOFBunchCrossing(Double_t b, Bool_t pidTPConly) const 
1259 {
1260   // Returns the number of bunch crossings after trigger (assuming 25ns spacing)
1261   const double kSpacing = 25e3; // min interbanch spacing
1262   const double kShift = 0;
1263   Int_t bcid = kTOFBCNA; // defualt one
1264   if (!IsOn(kTOFout) || !IsOn(kESDpid)) return bcid; // no info
1265   //
1266   double tdif = fTOFsignal;
1267   if (IsOn(kTIME)) { // integrated time info is there
1268     int pid = GetPID(pidTPConly);
1269     tdif -= fTrackTime[pid];
1270   }
1271   else { // assume integrated time info from TOF radius and momentum
1272     const double kRTOF = 385.;
1273     const double kCSpeed = 3.e-2; // cm/ps
1274     double p = GetP();
1275     if (p<0.01) return bcid;
1276     double m = GetMass(pidTPConly);
1277     double curv = GetC(b);
1278     double path = TMath::Abs(curv)>kAlmost0 ? // account for curvature
1279       2./curv*TMath::ASin(kRTOF*curv/2.)*TMath::Sqrt(1.+GetTgl()*GetTgl()) : kRTOF;
1280     tdif -= path/kCSpeed*TMath::Sqrt(1.+m*m/(p*p));
1281   }
1282   bcid = TMath::Nint((tdif - kShift)/kSpacing);
1283   return bcid;
1284 }
1285
1286 //______________________________________________________________________________
1287 Double_t AliESDtrack::M() const
1288 {
1289   // Returns the assumed mass
1290   // (the pion mass, if the particle can't be identified properly).
1291   static Bool_t printerr=kTRUE;
1292   if (printerr) {
1293      AliWarning("WARNING !!! ... THIS WILL BE PRINTED JUST ONCE !!!");
1294      printerr = kFALSE;
1295      AliWarning("This is the ESD mass. Use it with care !"); 
1296   }
1297   return GetMass(); 
1298 }
1299   
1300 //______________________________________________________________________________
1301 Double_t AliESDtrack::E() const
1302 {
1303   // Returns the energy of the particle given its assumed mass.
1304   // Assumes the pion mass if the particle can't be identified properly.
1305   
1306   Double_t m = M();
1307   Double_t p = P();
1308   return TMath::Sqrt(p*p + m*m);
1309 }
1310
1311 //______________________________________________________________________________
1312 Double_t AliESDtrack::Y() const
1313 {
1314   // Returns the rapidity of a particle given its assumed mass.
1315   // Assumes the pion mass if the particle can't be identified properly.
1316   
1317   Double_t e = E();
1318   Double_t pz = Pz();
1319   if (e != TMath::Abs(pz)) { // energy was not equal to pz
1320     return 0.5*TMath::Log((e+pz)/(e-pz));
1321   } else { // energy was equal to pz
1322     return -999.;
1323   }
1324 }
1325
1326 //_______________________________________________________________________
1327 Bool_t AliESDtrack::UpdateTrackParams(const AliKalmanTrack *t, ULong_t flags){
1328   //
1329   // This function updates track's running parameters 
1330   //
1331   Bool_t rc=kTRUE;
1332
1333   SetStatus(flags);
1334   fLabel=t->GetLabel();
1335
1336   if (t->IsStartedTimeIntegral()) {
1337     SetStatus(kTIME);
1338     Double_t times[10];t->GetIntegratedTimes(times); SetIntegratedTimes(times);
1339     SetIntegratedLength(t->GetIntegratedLength());
1340   }
1341
1342   Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
1343   if (fFriendTrack) {
1344   if (flags==kITSout) fFriendTrack->SetITSOut(*t);
1345   if (flags==kTPCout) fFriendTrack->SetTPCOut(*t);
1346   if (flags==kTRDrefit) fFriendTrack->SetTRDIn(*t);
1347   }
1348   
1349   switch (flags) {
1350     
1351   case kITSin: 
1352     fITSchi2Std[0] = t->GetChi2();
1353     //
1354   case kITSout: 
1355     fITSchi2Std[1] = t->GetChi2();
1356   case kITSrefit:
1357     {
1358     fITSchi2Std[2] = t->GetChi2();
1359     fITSClusterMap=0;
1360     fITSncls=t->GetNumberOfClusters();
1361     if (fFriendTrack) {
1362     Int_t* indexITS = new Int_t[AliESDfriendTrack::kMaxITScluster];
1363     for (Int_t i=0;i<AliESDfriendTrack::kMaxITScluster;i++) {
1364         indexITS[i]=t->GetClusterIndex(i);
1365
1366         if (i<fITSncls) {
1367           Int_t l=(indexITS[i] & 0xf0000000) >> 28;
1368            SETBIT(fITSClusterMap,l);                 
1369         }
1370     }
1371     fFriendTrack->SetITSIndices(indexITS,AliESDfriendTrack::kMaxITScluster);
1372     delete [] indexITS;
1373     }
1374
1375     fITSchi2=t->GetChi2();
1376     fITSsignal=t->GetPIDsignal();
1377     fITSLabel = t->GetLabel();
1378     // keep in fOp the parameters outside ITS for ITS stand-alone tracks 
1379     if (flags==kITSout) { 
1380       if (!fOp) fOp=new AliExternalTrackParam(*t);
1381       else 
1382         fOp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
1383     }   
1384     }
1385     break;
1386     
1387   case kTPCin: case kTPCrefit:
1388     {
1389     fTPCLabel = t->GetLabel();
1390     if (flags==kTPCin)  {
1391       fTPCInner=new AliExternalTrackParam(*t); 
1392       fTPCnclsIter1=t->GetNumberOfClusters();    
1393       fTPCchi2Iter1=t->GetChi2();
1394     }
1395     if (!fIp) fIp=new AliExternalTrackParam(*t);
1396     else 
1397       fIp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
1398     }
1399     // Intentionally no break statement; need to set general TPC variables as well
1400   case kTPCout:
1401     {
1402     if (flags & kTPCout){
1403       if (!fOp) fOp=new AliExternalTrackParam(*t);
1404       else 
1405         fOp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
1406     }
1407     fTPCncls=t->GetNumberOfClusters();    
1408     fTPCchi2=t->GetChi2();
1409     
1410     if (fFriendTrack) {  // Copy cluster indices
1411       Int_t* indexTPC = new Int_t[AliESDfriendTrack::kMaxTPCcluster];
1412       for (Int_t i=0;i<AliESDfriendTrack::kMaxTPCcluster;i++)         
1413         indexTPC[i]=t->GetClusterIndex(i);
1414       fFriendTrack->SetTPCIndices(indexTPC,AliESDfriendTrack::kMaxTPCcluster);
1415       delete [] indexTPC;
1416     }
1417     fTPCsignal=t->GetPIDsignal();
1418     }
1419     break;
1420
1421   case kTRDin: case kTRDrefit:
1422     break;
1423   case kTRDout:
1424     {
1425     fTRDLabel = t->GetLabel(); 
1426     fTRDchi2  = t->GetChi2();
1427     fTRDncls  = t->GetNumberOfClusters();
1428     if (fFriendTrack) {
1429       Int_t* indexTRD = new Int_t[AliESDfriendTrack::kMaxTRDcluster];
1430       for (Int_t i=0;i<AliESDfriendTrack::kMaxTRDcluster;i++) indexTRD[i]=-2;
1431       for (Int_t i=0;i<6;i++) indexTRD[i]=t->GetTrackletIndex(i);
1432       fFriendTrack->SetTRDIndices(indexTRD,AliESDfriendTrack::kMaxTRDcluster);
1433       delete [] indexTRD;
1434     }    
1435     
1436     //commented out by Xianguo
1437     //fTRDsignal=t->GetPIDsignal();
1438     }
1439     break;
1440   case kTRDbackup:
1441     if (!fOp) fOp=new AliExternalTrackParam(*t);
1442     else 
1443       fOp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
1444     fTRDncls0 = t->GetNumberOfClusters(); 
1445     break;
1446   case kTOFin: 
1447     break;
1448   case kTOFout: 
1449     break;
1450   case kTRDStop:
1451     break;
1452   case kHMPIDout:
1453   if (!fHMPIDp) fHMPIDp=new AliExternalTrackParam(*t);
1454     else 
1455       fHMPIDp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
1456     break;
1457   default: 
1458     AliError("Wrong flag !");
1459     return kFALSE;
1460   }
1461
1462   return rc;
1463 }
1464
1465 //_______________________________________________________________________
1466 void AliESDtrack::GetExternalParameters(Double_t &x, Double_t p[5]) const {
1467   //---------------------------------------------------------------------
1468   // This function returns external representation of the track parameters
1469   //---------------------------------------------------------------------
1470   x=GetX();
1471   for (Int_t i=0; i<5; i++) p[i]=GetParameter()[i];
1472 }
1473
1474 //_______________________________________________________________________
1475 void AliESDtrack::GetExternalCovariance(Double_t cov[15]) const {
1476   //---------------------------------------------------------------------
1477   // This function returns external representation of the cov. matrix
1478   //---------------------------------------------------------------------
1479   for (Int_t i=0; i<15; i++) cov[i]=AliExternalTrackParam::GetCovariance()[i];
1480 }
1481
1482 //_______________________________________________________________________
1483 Bool_t AliESDtrack::GetConstrainedExternalParameters
1484                  (Double_t &alpha, Double_t &x, Double_t p[5]) const {
1485   //---------------------------------------------------------------------
1486   // This function returns the constrained external track parameters
1487   //---------------------------------------------------------------------
1488   if (!fCp) return kFALSE;
1489   alpha=fCp->GetAlpha();
1490   x=fCp->GetX();
1491   for (Int_t i=0; i<5; i++) p[i]=fCp->GetParameter()[i];
1492   return kTRUE;
1493 }
1494
1495 //_______________________________________________________________________
1496 Bool_t 
1497 AliESDtrack::GetConstrainedExternalCovariance(Double_t c[15]) const {
1498   //---------------------------------------------------------------------
1499   // This function returns the constrained external cov. matrix
1500   //---------------------------------------------------------------------
1501   if (!fCp) return kFALSE;
1502   for (Int_t i=0; i<15; i++) c[i]=fCp->GetCovariance()[i];
1503   return kTRUE;
1504 }
1505
1506 Bool_t
1507 AliESDtrack::GetInnerExternalParameters
1508                  (Double_t &alpha, Double_t &x, Double_t p[5]) const {
1509   //---------------------------------------------------------------------
1510   // This function returns external representation of the track parameters 
1511   // at the inner layer of TPC
1512   //---------------------------------------------------------------------
1513   if (!fIp) return kFALSE;
1514   alpha=fIp->GetAlpha();
1515   x=fIp->GetX();
1516   for (Int_t i=0; i<5; i++) p[i]=fIp->GetParameter()[i];
1517   return kTRUE;
1518 }
1519
1520 Bool_t 
1521 AliESDtrack::GetInnerExternalCovariance(Double_t cov[15]) const {
1522  //---------------------------------------------------------------------
1523  // This function returns external representation of the cov. matrix 
1524  // at the inner layer of TPC
1525  //---------------------------------------------------------------------
1526   if (!fIp) return kFALSE;
1527   for (Int_t i=0; i<15; i++) cov[i]=fIp->GetCovariance()[i];
1528   return kTRUE;
1529 }
1530
1531 void 
1532 AliESDtrack::SetOuterParam(const AliExternalTrackParam *p, ULong_t flags) {
1533   //
1534   // This is a direct setter for the outer track parameters
1535   //
1536   SetStatus(flags);
1537   if (fOp) delete fOp;
1538   fOp=new AliExternalTrackParam(*p);
1539 }
1540
1541 void 
1542 AliESDtrack::SetOuterHmpParam(const AliExternalTrackParam *p, ULong_t flags) {
1543   //
1544   // This is a direct setter for the outer track parameters
1545   //
1546   SetStatus(flags);
1547   if (fHMPIDp) delete fHMPIDp;
1548   fHMPIDp=new AliExternalTrackParam(*p);
1549 }
1550
1551 Bool_t 
1552 AliESDtrack::GetOuterExternalParameters
1553                  (Double_t &alpha, Double_t &x, Double_t p[5]) const {
1554   //---------------------------------------------------------------------
1555   // This function returns external representation of the track parameters 
1556   // at the inner layer of TRD
1557   //---------------------------------------------------------------------
1558   if (!fOp) return kFALSE;
1559   alpha=fOp->GetAlpha();
1560   x=fOp->GetX();
1561   for (Int_t i=0; i<5; i++) p[i]=fOp->GetParameter()[i];
1562   return kTRUE;
1563 }
1564
1565 Bool_t 
1566 AliESDtrack::GetOuterHmpExternalParameters
1567                  (Double_t &alpha, Double_t &x, Double_t p[5]) const {
1568   //---------------------------------------------------------------------
1569   // This function returns external representation of the track parameters 
1570   // at the inner layer of TRD
1571   //---------------------------------------------------------------------
1572   if (!fHMPIDp) return kFALSE;
1573   alpha=fHMPIDp->GetAlpha();
1574   x=fHMPIDp->GetX();
1575   for (Int_t i=0; i<5; i++) p[i]=fHMPIDp->GetParameter()[i];
1576   return kTRUE;
1577 }
1578
1579 Bool_t 
1580 AliESDtrack::GetOuterExternalCovariance(Double_t cov[15]) const {
1581  //---------------------------------------------------------------------
1582  // This function returns external representation of the cov. matrix 
1583  // at the inner layer of TRD
1584  //---------------------------------------------------------------------
1585   if (!fOp) return kFALSE;
1586   for (Int_t i=0; i<15; i++) cov[i]=fOp->GetCovariance()[i];
1587   return kTRUE;
1588 }
1589
1590 Bool_t 
1591 AliESDtrack::GetOuterHmpExternalCovariance(Double_t cov[15]) const {
1592  //---------------------------------------------------------------------
1593  // This function returns external representation of the cov. matrix 
1594  // at the inner layer of TRD
1595  //---------------------------------------------------------------------
1596   if (!fHMPIDp) return kFALSE;
1597   for (Int_t i=0; i<15; i++) cov[i]=fHMPIDp->GetCovariance()[i];
1598   return kTRUE;
1599 }
1600
1601 Int_t AliESDtrack::GetNcls(Int_t idet) const
1602 {
1603   // Get number of clusters by subdetector index
1604   //
1605   Int_t ncls = 0;
1606   switch(idet){
1607   case 0:
1608     ncls = fITSncls;
1609     break;
1610   case 1:
1611     ncls = fTPCncls;
1612     break;
1613   case 2:
1614     ncls = fTRDncls;
1615     break;
1616   case 3:
1617     if (fTOFindex != -1)
1618       ncls = 1;
1619     break;
1620   case 4: //PHOS
1621     break;
1622   case 5: //HMPID
1623     if ((fHMPIDcluIdx >= 0) && (fHMPIDcluIdx < 7000000)) {
1624       if ((fHMPIDcluIdx%1000000 != 9999) && (fHMPIDcluIdx%1000000 != 99999)) {
1625         ncls = 1;
1626       }
1627     }    
1628     break;
1629   default:
1630     break;
1631   }
1632   return ncls;
1633 }
1634
1635 Int_t AliESDtrack::GetClusters(Int_t idet, Int_t *idx) const
1636 {
1637   // Get cluster index array by subdetector index
1638   //
1639   Int_t ncls = 0;
1640   switch(idet){
1641   case 0:
1642     ncls = GetITSclusters(idx);
1643     break;
1644   case 1:
1645     ncls = GetTPCclusters(idx);
1646     break;
1647   case 2:
1648     ncls = GetTRDclusters(idx);
1649     break;
1650   case 3:
1651     if (fTOFindex != -1) {
1652       idx[0] = fTOFindex;
1653       ncls = 1;
1654     }
1655     break;
1656   case 4: //PHOS
1657     break;
1658   case 5:
1659     if ((fHMPIDcluIdx >= 0) && (fHMPIDcluIdx < 7000000)) {
1660       if ((fHMPIDcluIdx%1000000 != 9999) && (fHMPIDcluIdx%1000000 != 99999)) {
1661         idx[0] = GetHMPIDcluIdx();
1662         ncls = 1;
1663       }
1664     }    
1665     break;
1666   case 6: //EMCAL
1667     break;
1668   default:
1669     break;
1670   }
1671   return ncls;
1672 }
1673
1674 //_______________________________________________________________________
1675 void AliESDtrack::GetIntegratedTimes(Double_t *times) const {
1676   // Returns the array with integrated times for each particle hypothesis
1677   for (Int_t i=0; i<AliPID::kSPECIES; i++) times[i]=fTrackTime[i];
1678 }
1679
1680 //_______________________________________________________________________
1681 void AliESDtrack::SetIntegratedTimes(const Double_t *times) {
1682   // Sets the array with integrated times for each particle hypotesis
1683   for (Int_t i=0; i<AliPID::kSPECIES; i++) fTrackTime[i]=times[i];
1684 }
1685
1686 //_______________________________________________________________________
1687 void AliESDtrack::SetITSpid(const Double_t *p) {
1688   // Sets values for the probability of each particle type (in ITS)
1689   SetPIDValues(fITSr,p,AliPID::kSPECIES);
1690   SetStatus(AliESDtrack::kITSpid);
1691 }
1692
1693 //_______________________________________________________________________
1694 void AliESDtrack::GetITSpid(Double_t *p) const {
1695   // Gets the probability of each particle type (in ITS)
1696   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fITSr[i];
1697 }
1698
1699 //_______________________________________________________________________
1700 Char_t AliESDtrack::GetITSclusters(Int_t *idx) const {
1701   //---------------------------------------------------------------------
1702   // This function returns indices of the assgined ITS clusters 
1703   //---------------------------------------------------------------------
1704   if (idx && fFriendTrack) {
1705     Int_t *index=fFriendTrack->GetITSindices();
1706     for (Int_t i=0; i<AliESDfriendTrack::kMaxITScluster; i++) {
1707       if ( (i>=fITSncls) && (i<6) ) idx[i]=-1;
1708       else {
1709         if (index) {
1710           idx[i]=index[i];
1711         }
1712         else idx[i]= -2;
1713       }
1714     }
1715   }
1716   return fITSncls;
1717 }
1718
1719 //_______________________________________________________________________
1720 Bool_t AliESDtrack::GetITSModuleIndexInfo(Int_t ilayer,Int_t &idet,Int_t &status,
1721                                          Float_t &xloc,Float_t &zloc) const {
1722   //----------------------------------------------------------------------
1723   // This function encodes in the module number also the status of cluster association
1724   // "status" can have the following values: 
1725   // 1 "found" (cluster is associated), 
1726   // 2 "dead" (module is dead from OCDB), 
1727   // 3 "skipped" (module or layer forced to be skipped),
1728   // 4 "outinz" (track out of z acceptance), 
1729   // 5 "nocls" (no clusters in the road), 
1730   // 6 "norefit" (cluster rejected during refit), 
1731   // 7 "deadzspd" (holes in z in SPD)
1732   // Also given are the coordinates of the crossing point of track and module
1733   // (in the local module ref. system)
1734   // WARNING: THIS METHOD HAS TO BE SYNCHRONIZED WITH AliITStrackV2::GetModuleIndexInfo()!
1735   //----------------------------------------------------------------------
1736
1737   if(fITSModule[ilayer]==-1) {
1738     idet = -1;
1739     status=0;
1740     xloc=-99.; zloc=-99.;
1741     return kFALSE;
1742   }
1743
1744   Int_t module = fITSModule[ilayer];
1745
1746   idet = Int_t(module/1000000);
1747
1748   module -= idet*1000000;
1749
1750   status = Int_t(module/100000);
1751
1752   module -= status*100000;
1753
1754   Int_t signs = Int_t(module/10000);
1755
1756   module-=signs*10000;
1757
1758   Int_t xInt = Int_t(module/100);
1759   module -= xInt*100;
1760
1761   Int_t zInt = module;
1762
1763   if(signs==1) { xInt*=1; zInt*=1; }
1764   if(signs==2) { xInt*=1; zInt*=-1; }
1765   if(signs==3) { xInt*=-1; zInt*=1; }
1766   if(signs==4) { xInt*=-1; zInt*=-1; }
1767
1768   xloc = 0.1*(Float_t)xInt;
1769   zloc = 0.1*(Float_t)zInt;
1770
1771   if(status==4) idet = -1;
1772
1773   return kTRUE;
1774 }
1775
1776 //_______________________________________________________________________
1777 UShort_t AliESDtrack::GetTPCclusters(Int_t *idx) const {
1778   //---------------------------------------------------------------------
1779   // This function returns indices of the assgined ITS clusters 
1780   //---------------------------------------------------------------------
1781   if (idx && fFriendTrack) {
1782     Int_t *index=fFriendTrack->GetTPCindices();
1783
1784     if (index){
1785       for (Int_t i=0; i<AliESDfriendTrack::kMaxTPCcluster; i++) idx[i]=index[i];
1786     }
1787     else {
1788       for (Int_t i=0; i<AliESDfriendTrack::kMaxTPCcluster; i++) idx[i]=-2;
1789     }
1790   }
1791   return fTPCncls;
1792 }
1793
1794 //_______________________________________________________________________
1795 Float_t AliESDtrack::GetTPCCrossedRows() const
1796 {
1797   // This function calls GetTPCClusterInfo with some default parameters which are used in the track selection and caches the outcome
1798   // because GetTPCClusterInfo is quite time-consuming
1799   
1800   if (fCacheNCrossedRows > -1)
1801     return fCacheNCrossedRows;
1802   
1803   fCacheNCrossedRows = GetTPCClusterInfo(2, 1);
1804   return fCacheNCrossedRows;
1805 }
1806
1807 //_______________________________________________________________________
1808 Float_t AliESDtrack::GetTPCClusterInfo(Int_t nNeighbours/*=3*/, Int_t type/*=0*/, Int_t row0, Int_t row1, Int_t bitType ) const
1809 {
1810   //
1811   // TPC cluster information
1812   // type 0: get fraction of found/findable clusters with neighbourhood definition
1813   //      1: findable clusters with neighbourhood definition
1814   //      2: found clusters
1815   // bitType:
1816   //      0 - all cluster used
1817   //      1 - clusters  used for the kalman update
1818   // definition of findable clusters:
1819   //            a cluster is defined as findable if there is another cluster
1820   //           within +- nNeighbours pad rows. The idea is to overcome threshold
1821   //           effects with a very simple algorithm.
1822   //
1823
1824   
1825   Int_t found=0;
1826   Int_t findable=0;
1827   Int_t last=-nNeighbours;
1828   const TBits & clusterMap = (bitType%2==0) ? fTPCClusterMap : fTPCFitMap;
1829   
1830   Int_t upperBound=clusterMap.GetNbits();
1831   if (upperBound>row1) upperBound=row1;
1832   for (Int_t i=row0; i<upperBound; ++i){
1833     //look to current row
1834     if (clusterMap[i]) {
1835       last=i;
1836       ++found;
1837       ++findable;
1838       continue;
1839     }
1840     //look to nNeighbours before
1841     if ((i-last)<=nNeighbours) {
1842       ++findable;
1843       continue;
1844     }
1845     //look to nNeighbours after
1846     for (Int_t j=i+1; j<i+1+nNeighbours; ++j){
1847       if (clusterMap[j]){
1848         ++findable;
1849         break;
1850       }
1851     }
1852   }
1853   if (type==2) return found;
1854   if (type==1) return findable;
1855   
1856   if (type==0){
1857     Float_t fraction=0;
1858     if (findable>0) 
1859       fraction=(Float_t)found/(Float_t)findable;
1860     else 
1861       fraction=0;
1862     return fraction;
1863   }  
1864   return 0;  // undefined type - default value
1865 }
1866
1867 //_______________________________________________________________________
1868 Float_t AliESDtrack::GetTPCClusterDensity(Int_t nNeighbours/*=3*/, Int_t type/*=0*/, Int_t row0, Int_t row1, Int_t bitType ) const
1869 {
1870   //
1871   // TPC cluster density -  only rows where signal before and after given row are used
1872   //                     -  slower function
1873   // type 0: get fraction of found/findable clusters with neighbourhood definition
1874   //      1: findable clusters with neighbourhood definition
1875   //      2: found clusters
1876   // bitType:
1877   //      0 - all cluster used
1878   //      1 - clusters  used for the kalman update
1879   // definition of findable clusters:
1880   //            a cluster is defined as findable if there is another cluster
1881   //           within +- nNeighbours pad rows. The idea is to overcome threshold
1882   //           effects with a very simple algorithm.
1883   //  
1884   Int_t found=0;
1885   Int_t findable=0;
1886   //  Int_t last=-nNeighbours;
1887   const TBits & clusterMap = (bitType%2==0) ? fTPCClusterMap : fTPCFitMap;
1888   Int_t upperBound=clusterMap.GetNbits();
1889   if (upperBound>row1) upperBound=row1;
1890   for (Int_t i=row0; i<upperBound; ++i){
1891     Bool_t isUp=kFALSE;
1892     Bool_t isDown=kFALSE;
1893     for (Int_t idelta=1; idelta<=nNeighbours; idelta++){
1894       if (i-idelta>=0 && clusterMap[i-idelta]) isDown=kTRUE;
1895       if (i+idelta<upperBound && clusterMap[i+idelta]) isUp=kTRUE;
1896     }
1897     if (isUp&&isDown){
1898       ++findable;
1899       if (clusterMap[i]) ++found;
1900     }
1901   }
1902   if (type==2) return found;
1903   if (type==1) return findable;
1904   
1905   if (type==0){
1906     Float_t fraction=0;
1907     if (findable>0) 
1908       fraction=(Float_t)found/(Float_t)findable;
1909     else 
1910       fraction=0;
1911     return fraction;
1912   }  
1913   return 0;  // undefined type - default value
1914 }
1915
1916
1917
1918
1919 //_______________________________________________________________________
1920 Double_t AliESDtrack::GetTPCdensity(Int_t row0, Int_t row1) const{
1921   //
1922   // GetDensity of the clusters on given region between row0 and row1
1923   // Dead zone effect takin into acoount
1924   //
1925   if (!fFriendTrack) return 0.0;
1926   Int_t good  = 0;
1927   Int_t found = 0;
1928   //  
1929   Int_t *index=fFriendTrack->GetTPCindices();
1930   for (Int_t i=row0;i<=row1;i++){     
1931     Int_t idx = index[i];
1932     if (idx!=-1)  good++;             // track outside of dead zone
1933     if (idx>0)    found++;
1934   }
1935   Float_t density=0.5;
1936   if (good>TMath::Max((row1-row0)*0.5,0.0)) density = Float_t(found)/Float_t(good);
1937   return density;
1938 }
1939
1940 //_______________________________________________________________________
1941 void AliESDtrack::SetTPCpid(const Double_t *p) {  
1942   // Sets values for the probability of each particle type (in TPC)
1943   SetPIDValues(fTPCr,p,AliPID::kSPECIES);
1944   SetStatus(AliESDtrack::kTPCpid);
1945 }
1946
1947 //_______________________________________________________________________
1948 void AliESDtrack::GetTPCpid(Double_t *p) const {
1949   // Gets the probability of each particle type (in TPC)
1950   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fTPCr[i];
1951 }
1952
1953 //_______________________________________________________________________
1954 UChar_t AliESDtrack::GetTRDclusters(Int_t *idx) const {
1955   //---------------------------------------------------------------------
1956   // This function returns indices of the assgined TRD clusters 
1957   //---------------------------------------------------------------------
1958   if (idx && fFriendTrack) {
1959     Int_t *index=fFriendTrack->GetTRDindices();
1960
1961     if (index) {
1962       for (Int_t i=0; i<AliESDfriendTrack::kMaxTRDcluster; i++) idx[i]=index[i];
1963     }
1964     else {
1965       for (Int_t i=0; i<AliESDfriendTrack::kMaxTRDcluster; i++) idx[i]=-2;
1966     }
1967   }
1968   return fTRDncls;
1969 }
1970
1971 //_______________________________________________________________________
1972 UChar_t AliESDtrack::GetTRDtracklets(Int_t *idx) const {
1973 //
1974 // This function returns the number of TRD tracklets used in tracking
1975 // and it fills the indices of these tracklets in the array "idx" as they 
1976 // are registered in the TRD track list. 
1977 // 
1978 // Caution :
1979 //   1. The idx array has to be allocated with a size >= AliESDtrack::kTRDnPlanes
1980 //   2. The idx array store not only the index but also the layer of the tracklet. 
1981 //      Therefore tracks with TRD gaps contain default values for indices [-1] 
1982
1983   if (!fFriendTrack) return 0;
1984   if (!idx) return GetTRDntracklets();
1985   Int_t *index=fFriendTrack->GetTRDindices();
1986   Int_t n = 0;
1987   for (Int_t i=0; i<kTRDnPlanes; i++){ 
1988     if (index){
1989       if(index[i]>=0) n++;
1990       idx[i]=index[i];
1991     }
1992     else idx[i] = -2;
1993   }
1994   return n;
1995 }
1996
1997 //_______________________________________________________________________
1998 void AliESDtrack::SetTRDpid(const Double_t *p) {  
1999   // Sets values for the probability of each particle type (in TRD)
2000   SetPIDValues(fTRDr,p,AliPID::kSPECIES);
2001   SetStatus(AliESDtrack::kTRDpid);
2002 }
2003
2004 //_______________________________________________________________________
2005 void AliESDtrack::GetTRDpid(Double_t *p) const {
2006   // Gets the probability of each particle type (in TRD)
2007   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fTRDr[i];
2008 }
2009
2010 //_______________________________________________________________________
2011 void    AliESDtrack::SetTRDpid(Int_t iSpecies, Float_t p)
2012 {
2013   // Sets the probability of particle type iSpecies to p (in TRD)
2014   fTRDr[iSpecies] = p;
2015 }
2016
2017 Double_t AliESDtrack::GetTRDpid(Int_t iSpecies) const
2018 {
2019   // Returns the probability of particle type iSpecies (in TRD)
2020   return fTRDr[iSpecies];
2021 }
2022
2023 //____________________________________________________
2024 Int_t AliESDtrack::GetNumberOfTRDslices() const 
2025 {
2026   // built in backward compatibility
2027   Int_t idx = fTRDnSlices - (kTRDnPlanes<<1);
2028   return idx<18 ? fTRDnSlices/kTRDnPlanes : idx/kTRDnPlanes;
2029 }
2030
2031 //____________________________________________________
2032 Double_t AliESDtrack::GetTRDmomentum(Int_t plane, Double_t *sp) const
2033 {
2034 //Returns momentum estimation and optional its error (sp)
2035 // in TRD layer "plane".
2036
2037   if (!fTRDnSlices) {
2038     AliDebug(2, "No TRD info allocated for this track.");
2039     return -1.;
2040   }
2041   if ((plane<0) || (plane>=kTRDnPlanes)) {
2042     AliWarning(Form("Request for TRD plane[%d] outside range.", plane)); 
2043     return -1.;
2044   }
2045
2046   Int_t idx = fTRDnSlices-(kTRDnPlanes<<1)+plane;
2047   // Protection for backward compatibility
2048   if(idx<(GetNumberOfTRDslices()*kTRDnPlanes)) return -1.;
2049
2050   if(sp) (*sp) = fTRDslices[idx+kTRDnPlanes];
2051   return fTRDslices[idx];
2052 }
2053
2054 //____________________________________________________
2055 Double_t  AliESDtrack::GetTRDslice(Int_t plane, Int_t slice) const {
2056   //Gets the charge from the slice of the plane
2057
2058   if(!fTRDslices) {
2059     //AliError("No TRD slices allocated for this track !");
2060     return -1.;
2061   }
2062   if ((plane<0) || (plane>=kTRDnPlanes)) {
2063     AliError("Info for TRD plane not available !");
2064     return -1.;
2065   }
2066   Int_t ns=GetNumberOfTRDslices();
2067   if ((slice<-1) || (slice>=ns)) {
2068     //AliError("Wrong TRD slice !");  
2069     return -1.;
2070   }
2071
2072   if(slice>=0) return fTRDslices[plane*ns + slice];
2073
2074   // return average of the dEdx measurements
2075   Double_t q=0.; Double32_t *s = &fTRDslices[plane*ns];
2076   for (Int_t i=0; i<ns; i++, s++) if((*s)>0.) q+=(*s);
2077   return q/ns;
2078 }
2079
2080 //____________________________________________________
2081 void  AliESDtrack::SetNumberOfTRDslices(Int_t n) {
2082   //Sets the number of slices used for PID 
2083   if (fTRDnSlices) return;
2084
2085   fTRDnSlices=n;
2086   fTRDslices=new Double32_t[fTRDnSlices];
2087   
2088   // set-up correctly the allocated memory
2089   memset(fTRDslices, 0, n*sizeof(Double32_t));
2090   for (Int_t i=GetNumberOfTRDslices(); i--;) fTRDslices[i]=-1.;
2091 }
2092
2093 //____________________________________________________
2094 void  AliESDtrack::SetTRDslice(Double_t q, Int_t plane, Int_t slice) {
2095   //Sets the charge q in the slice of the plane
2096   if(!fTRDslices) {
2097     AliError("No TRD slices allocated for this track !");
2098     return;
2099   }
2100   if ((plane<0) || (plane>=kTRDnPlanes)) {
2101     AliError("Info for TRD plane not allocated !");
2102     return;
2103   }
2104   Int_t ns=GetNumberOfTRDslices();
2105   if ((slice<0) || (slice>=ns)) {
2106     AliError("Wrong TRD slice !");
2107     return;
2108   }
2109   Int_t n=plane*ns + slice;
2110   fTRDslices[n]=q;
2111 }
2112
2113
2114 //____________________________________________________
2115 void AliESDtrack::SetTRDmomentum(Double_t p, Int_t plane, Double_t *sp)
2116 {
2117   if(!fTRDslices) {
2118     AliError("No TRD slices allocated for this track !");
2119     return;
2120   }
2121   if ((plane<0) || (plane>=kTRDnPlanes)) {
2122     AliError("Info for TRD plane not allocated !");
2123     return;
2124   }
2125
2126   Int_t idx = fTRDnSlices-(kTRDnPlanes<<1)+plane;
2127   // Protection for backward compatibility
2128   if(idx<GetNumberOfTRDslices()*kTRDnPlanes) return;
2129
2130   if(sp) fTRDslices[idx+kTRDnPlanes] = (*sp);
2131   fTRDslices[idx] = p;
2132 }
2133
2134
2135 //_______________________________________________________________________
2136 void AliESDtrack::SetTOFpid(const Double_t *p) {  
2137   // Sets the probability of each particle type (in TOF)
2138   SetPIDValues(fTOFr,p,AliPID::kSPECIES);
2139   SetStatus(AliESDtrack::kTOFpid);
2140 }
2141
2142 //_______________________________________________________________________
2143 void AliESDtrack::SetTOFLabel(const Int_t *p) {  
2144   // Sets  (in TOF)
2145   for (Int_t i=0; i<3; i++) fTOFLabel[i]=p[i];
2146 }
2147
2148 //_______________________________________________________________________
2149 void AliESDtrack::GetTOFpid(Double_t *p) const {
2150   // Gets probabilities of each particle type (in TOF)
2151   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fTOFr[i];
2152 }
2153
2154 //_______________________________________________________________________
2155 void AliESDtrack::GetTOFLabel(Int_t *p) const {
2156   // Gets (in TOF)
2157   for (Int_t i=0; i<3; i++) p[i]=fTOFLabel[i];
2158 }
2159
2160 //_______________________________________________________________________
2161 void AliESDtrack::GetTOFInfo(Float_t *info) const {
2162   // Gets (in TOF)
2163   for (Int_t i=0; i<10; i++) info[i]=fTOFInfo[i];
2164 }
2165
2166 //_______________________________________________________________________
2167 void AliESDtrack::SetTOFInfo(Float_t*info) {
2168   // Gets (in TOF)
2169   for (Int_t i=0; i<10; i++) fTOFInfo[i]=info[i];
2170 }
2171
2172
2173
2174 //_______________________________________________________________________
2175 void AliESDtrack::SetHMPIDpid(const Double_t *p) {  
2176   // Sets the probability of each particle type (in HMPID)
2177   SetPIDValues(fHMPIDr,p,AliPID::kSPECIES);
2178   SetStatus(AliESDtrack::kHMPIDpid);
2179 }
2180
2181 //_______________________________________________________________________
2182 void  AliESDtrack::SetTPCdEdxInfo(AliTPCdEdxInfo * dEdxInfo){ 
2183   if(fTPCdEdxInfo) delete fTPCdEdxInfo;
2184   fTPCdEdxInfo = dEdxInfo; 
2185 }
2186
2187 //_______________________________________________________________________
2188 void AliESDtrack::GetHMPIDpid(Double_t *p) const {
2189   // Gets probabilities of each particle type (in HMPID)
2190   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fHMPIDr[i];
2191 }
2192
2193
2194
2195 //_______________________________________________________________________
2196 void AliESDtrack::SetESDpid(const Double_t *p) {  
2197   // Sets the probability of each particle type for the ESD track
2198   SetPIDValues(fR,p,AliPID::kSPECIES);
2199   SetStatus(AliESDtrack::kESDpid);
2200 }
2201
2202 //_______________________________________________________________________
2203 void AliESDtrack::GetESDpid(Double_t *p) const {
2204   // Gets probability of each particle type for the ESD track
2205   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fR[i];
2206 }
2207
2208 //_______________________________________________________________________
2209 Bool_t AliESDtrack::RelateToVertexTPC(const AliESDVertex *vtx, 
2210 Double_t b, Double_t maxd, AliExternalTrackParam *cParam) {
2211   //
2212   // Try to relate the TPC-only track parameters to the vertex "vtx", 
2213   // if the (rough) transverse impact parameter is not bigger then "maxd". 
2214   //            Magnetic field is "b" (kG).
2215   //
2216   // a) The TPC-only paramters are extapolated to the DCA to the vertex.
2217   // b) The impact parameters and their covariance matrix are calculated.
2218   // c) An attempt to constrain the TPC-only params to the vertex is done.
2219   //    The constrained params are returned via "cParam".
2220   //
2221   // In the case of success, the returned value is kTRUE
2222   // otherwise, it's kFALSE)
2223   // 
2224
2225   if (!fTPCInner) return kFALSE;
2226   if (!vtx) return kFALSE;
2227
2228   Double_t dz[2],cov[3];
2229   if (!fTPCInner->PropagateToDCA(vtx, b, maxd, dz, cov)) return kFALSE;
2230
2231   fdTPC = dz[0];
2232   fzTPC = dz[1];  
2233   fCddTPC = cov[0];
2234   fCdzTPC = cov[1];
2235   fCzzTPC = cov[2];
2236   
2237   Double_t covar[6]; vtx->GetCovMatrix(covar);
2238   Double_t p[2]={GetParameter()[0]-dz[0],GetParameter()[1]-dz[1]};
2239   Double_t c[3]={covar[2],0.,covar[5]};
2240
2241   Double_t chi2=GetPredictedChi2(p,c);
2242   if (chi2>kVeryBig) return kFALSE;
2243
2244   fCchi2TPC=chi2;
2245
2246   if (!cParam) return kTRUE;
2247
2248   *cParam = *fTPCInner;
2249   if (!cParam->Update(p,c)) return kFALSE;
2250
2251   return kTRUE;
2252 }
2253
2254 //_______________________________________________________________________
2255 Bool_t AliESDtrack::RelateToVertexTPCBxByBz(const AliESDVertex *vtx, 
2256 Double_t b[3], Double_t maxd, AliExternalTrackParam *cParam) {
2257   //
2258   // Try to relate the TPC-only track parameters to the vertex "vtx", 
2259   // if the (rough) transverse impact parameter is not bigger then "maxd". 
2260   //
2261   // All three components of the magnetic field ,"b[3]" (kG), 
2262   // are taken into account.
2263   //
2264   // a) The TPC-only paramters are extapolated to the DCA to the vertex.
2265   // b) The impact parameters and their covariance matrix are calculated.
2266   // c) An attempt to constrain the TPC-only params to the vertex is done.
2267   //    The constrained params are returned via "cParam".
2268   //
2269   // In the case of success, the returned value is kTRUE
2270   // otherwise, it's kFALSE)
2271   // 
2272
2273   if (!fTPCInner) return kFALSE;
2274   if (!vtx) return kFALSE;
2275
2276   Double_t dz[2],cov[3];
2277   if (!fTPCInner->PropagateToDCABxByBz(vtx, b, maxd, dz, cov)) return kFALSE;
2278
2279   fdTPC = dz[0];
2280   fzTPC = dz[1];  
2281   fCddTPC = cov[0];
2282   fCdzTPC = cov[1];
2283   fCzzTPC = cov[2];
2284   
2285   Double_t covar[6]; vtx->GetCovMatrix(covar);
2286   Double_t p[2]={GetParameter()[0]-dz[0],GetParameter()[1]-dz[1]};
2287   Double_t c[3]={covar[2],0.,covar[5]};
2288
2289   Double_t chi2=GetPredictedChi2(p,c);
2290   if (chi2>kVeryBig) return kFALSE;
2291
2292   fCchi2TPC=chi2;
2293
2294   if (!cParam) return kTRUE;
2295
2296   *cParam = *fTPCInner;
2297   if (!cParam->Update(p,c)) return kFALSE;
2298
2299   return kTRUE;
2300 }
2301
2302 //_______________________________________________________________________
2303 Bool_t AliESDtrack::RelateToVertex(const AliESDVertex *vtx, 
2304 Double_t b, Double_t maxd, AliExternalTrackParam *cParam) {
2305   //
2306   // Try to relate this track to the vertex "vtx", 
2307   // if the (rough) transverse impact parameter is not bigger then "maxd". 
2308   //            Magnetic field is "b" (kG).
2309   //
2310   // a) The track gets extapolated to the DCA to the vertex.
2311   // b) The impact parameters and their covariance matrix are calculated.
2312   // c) An attempt to constrain this track to the vertex is done.
2313   //    The constrained params are returned via "cParam".
2314   //
2315   // In the case of success, the returned value is kTRUE
2316   // (otherwise, it's kFALSE)
2317   //  
2318
2319   if (!vtx) return kFALSE;
2320
2321   Double_t dz[2],cov[3];
2322   if (!PropagateToDCA(vtx, b, maxd, dz, cov)) return kFALSE;
2323
2324   fD = dz[0];
2325   fZ = dz[1];  
2326   fCdd = cov[0];
2327   fCdz = cov[1];
2328   fCzz = cov[2];
2329   
2330   Double_t covar[6]; vtx->GetCovMatrix(covar);
2331   Double_t p[2]={GetParameter()[0]-dz[0],GetParameter()[1]-dz[1]};
2332   Double_t c[3]={covar[2],0.,covar[5]};
2333
2334   Double_t chi2=GetPredictedChi2(p,c);
2335   if (chi2>kVeryBig) return kFALSE;
2336
2337   fCchi2=chi2;
2338
2339
2340   //--- Could now these lines be removed ? ---
2341   delete fCp;
2342   fCp=new AliExternalTrackParam(*this);  
2343
2344   if (!fCp->Update(p,c)) {delete fCp; fCp=0; return kFALSE;}
2345   //----------------------------------------
2346
2347   fVertexID = vtx->GetID();
2348
2349   if (!cParam) return kTRUE;
2350
2351   *cParam = *this;
2352   if (!cParam->Update(p,c)) return kFALSE; 
2353
2354   return kTRUE;
2355 }
2356
2357 //_______________________________________________________________________
2358 Bool_t AliESDtrack::RelateToVertexBxByBz(const AliESDVertex *vtx, 
2359 Double_t b[3], Double_t maxd, AliExternalTrackParam *cParam) {
2360   //
2361   // Try to relate this track to the vertex "vtx", 
2362   // if the (rough) transverse impact parameter is not bigger then "maxd". 
2363   //            Magnetic field is "b" (kG).
2364   //
2365   // a) The track gets extapolated to the DCA to the vertex.
2366   // b) The impact parameters and their covariance matrix are calculated.
2367   // c) An attempt to constrain this track to the vertex is done.
2368   //    The constrained params are returned via "cParam".
2369   //
2370   // In the case of success, the returned value is kTRUE
2371   // (otherwise, it's kFALSE)
2372   //  
2373
2374   if (!vtx) return kFALSE;
2375
2376   Double_t dz[2],cov[3];
2377   if (!PropagateToDCABxByBz(vtx, b, maxd, dz, cov)) return kFALSE;
2378
2379   fD = dz[0];
2380   fZ = dz[1];  
2381   fCdd = cov[0];
2382   fCdz = cov[1];
2383   fCzz = cov[2];
2384   
2385   Double_t covar[6]; vtx->GetCovMatrix(covar);
2386   Double_t p[2]={GetParameter()[0]-dz[0],GetParameter()[1]-dz[1]};
2387   Double_t c[3]={covar[2],0.,covar[5]};
2388
2389   Double_t chi2=GetPredictedChi2(p,c);
2390   if (chi2>kVeryBig) return kFALSE;
2391
2392   fCchi2=chi2;
2393
2394
2395   //--- Could now these lines be removed ? ---
2396   delete fCp;
2397   fCp=new AliExternalTrackParam(*this);  
2398
2399   if (!fCp->Update(p,c)) {delete fCp; fCp=0; return kFALSE;}
2400   //----------------------------------------
2401
2402   fVertexID = vtx->GetID();
2403
2404   if (!cParam) return kTRUE;
2405
2406   *cParam = *this;
2407   if (!cParam->Update(p,c)) return kFALSE; 
2408
2409   return kTRUE;
2410 }
2411
2412 //_______________________________________________________________________
2413 void AliESDtrack::Print(Option_t *) const {
2414   // Prints info on the track
2415   AliExternalTrackParam::Print();
2416   printf("ESD track info\n") ; 
2417   Double_t p[AliPID::kSPECIESN] ; 
2418   Int_t index = 0 ; 
2419   if( IsOn(kITSpid) ){
2420     printf("From ITS: ") ; 
2421     GetITSpid(p) ; 
2422     for(index = 0 ; index < AliPID::kSPECIES; index++) 
2423       printf("%f, ", p[index]) ;
2424     printf("\n           signal = %f\n", GetITSsignal()) ;
2425   } 
2426   if( IsOn(kTPCpid) ){
2427     printf("From TPC: ") ; 
2428     GetTPCpid(p) ; 
2429     for(index = 0 ; index < AliPID::kSPECIES; index++) 
2430       printf("%f, ", p[index]) ;
2431     printf("\n           signal = %f\n", GetTPCsignal()) ;
2432   }
2433   if( IsOn(kTRDpid) ){
2434     printf("From TRD: ") ; 
2435     GetTRDpid(p) ; 
2436     for(index = 0 ; index < AliPID::kSPECIES; index++) 
2437       printf("%f, ", p[index]) ;
2438       printf("\n           signal = %f\n", GetTRDsignal()) ;
2439       printf("\n           NchamberdEdx = %d\n", GetTRDNchamberdEdx()) ;
2440       printf("\n           NclusterdEdx = %d\n", GetTRDNclusterdEdx()) ;
2441   }
2442   if( IsOn(kTOFpid) ){
2443     printf("From TOF: ") ; 
2444     GetTOFpid(p) ; 
2445     for(index = 0 ; index < AliPID::kSPECIES; index++) 
2446       printf("%f, ", p[index]) ;
2447     printf("\n           signal = %f\n", GetTOFsignal()) ;
2448   }
2449   if( IsOn(kHMPIDpid) ){
2450     printf("From HMPID: ") ; 
2451     GetHMPIDpid(p) ; 
2452     for(index = 0 ; index < AliPID::kSPECIES; index++) 
2453       printf("%f, ", p[index]) ;
2454     printf("\n           signal = %f\n", GetHMPIDsignal()) ;
2455   }
2456
2457
2458
2459 //
2460 // Draw functionality
2461 // Origin: Marian Ivanov, Marian.Ivanov@cern.ch
2462 //
2463 void AliESDtrack::FillPolymarker(TPolyMarker3D *pol, Float_t magF, Float_t minR, Float_t maxR, Float_t stepR){
2464   //
2465   // Fill points in the polymarker
2466   //
2467   TObjArray arrayRef;
2468   arrayRef.AddLast(new AliExternalTrackParam(*this));
2469   if (fIp) arrayRef.AddLast(new AliExternalTrackParam(*fIp));
2470   if (fOp) arrayRef.AddLast(new AliExternalTrackParam(*fOp));
2471   if (fHMPIDp) arrayRef.AddLast(new AliExternalTrackParam(*fHMPIDp));
2472   //
2473   Double_t mpos[3]={0,0,0};
2474   Int_t entries=arrayRef.GetEntries();
2475   for (Int_t i=0;i<entries;i++){
2476     Double_t pos[3];
2477     ((AliExternalTrackParam*)arrayRef.At(i))->GetXYZ(pos);
2478     mpos[0]+=pos[0]/entries;
2479     mpos[1]+=pos[1]/entries;
2480     mpos[2]+=pos[2]/entries;    
2481   }
2482   // Rotate to the mean position
2483   //
2484   Float_t fi= TMath::ATan2(mpos[1],mpos[0]);
2485   for (Int_t i=0;i<entries;i++){
2486     Bool_t res = ((AliExternalTrackParam*)arrayRef.At(i))->Rotate(fi);
2487     if (!res) delete arrayRef.RemoveAt(i);
2488   }
2489   Int_t counter=0;
2490   for (Double_t r=minR; r<maxR; r+=stepR){
2491     Double_t sweight=0;
2492     Double_t mlpos[3]={0,0,0};
2493     for (Int_t i=0;i<entries;i++){
2494       Double_t point[3]={0,0,0};
2495       AliExternalTrackParam *param = ((AliExternalTrackParam*)arrayRef.At(i));
2496       if (!param) continue;
2497       if (param->GetXYZAt(r,magF,point)){
2498         Double_t weight = 1./(10.+(r-param->GetX())*(r-param->GetX()));
2499         sweight+=weight;
2500         mlpos[0]+=point[0]*weight;
2501         mlpos[1]+=point[1]*weight;
2502         mlpos[2]+=point[2]*weight;
2503       }
2504     }
2505     if (sweight>0){
2506       mlpos[0]/=sweight;
2507       mlpos[1]/=sweight;
2508       mlpos[2]/=sweight;      
2509       pol->SetPoint(counter,mlpos[0],mlpos[1], mlpos[2]);
2510       //      printf("xyz\t%f\t%f\t%f\n",mlpos[0], mlpos[1],mlpos[2]);
2511       counter++;
2512     }
2513   }
2514 }
2515
2516 //_______________________________________________________________________
2517 void AliESDtrack::SetITSdEdxSamples(const Double_t s[4]) {
2518   //
2519   // Store the dE/dx samples measured by the two SSD and two SDD layers.
2520   // These samples are corrected for the track segment length. 
2521   //
2522   for (Int_t i=0; i<4; i++) fITSdEdxSamples[i]=s[i];
2523 }
2524
2525 //_______________________________________________________________________
2526 void AliESDtrack::GetITSdEdxSamples(Double_t *s) const {
2527   //
2528   // Get the dE/dx samples measured by the two SSD and two SDD layers.  
2529   // These samples are corrected for the track segment length.
2530   //
2531   for (Int_t i=0; i<4; i++) s[i]=fITSdEdxSamples[i];
2532 }
2533
2534
2535 UShort_t   AliESDtrack::GetTPCnclsS(Int_t i0,Int_t i1) const{
2536   //
2537   // get number of shared TPC clusters
2538   //
2539   return  fTPCSharedMap.CountBits(i0)-fTPCSharedMap.CountBits(i1);
2540 }
2541
2542 UShort_t   AliESDtrack::GetTPCncls(Int_t i0,Int_t i1) const{
2543   //
2544   // get number of TPC clusters
2545   //
2546   return  fTPCClusterMap.CountBits(i0)-fTPCClusterMap.CountBits(i1);
2547 }
2548
2549 //____________________________________________________________________
2550 Double_t AliESDtrack::GetChi2TPCConstrainedVsGlobal(const AliESDVertex* vtx) const
2551 {
2552   // Calculates the chi2 between the TPC track (TPCinner) constrained to the primary vertex and the global track
2553   //
2554   // Returns -1 in case the calculation failed
2555   //
2556   // Value is cached as a non-persistent member.
2557   //
2558   // Code adapted from original code by GSI group (Jacek, Marian, Michael)
2559   
2560   // cache, ignoring that a different vertex might be passed
2561   if (fCacheChi2TPCConstrainedVsGlobalVertex == vtx)
2562     return fCacheChi2TPCConstrainedVsGlobal;
2563   
2564   fCacheChi2TPCConstrainedVsGlobal = -1;
2565   fCacheChi2TPCConstrainedVsGlobalVertex = vtx;
2566   
2567   Double_t x[3];
2568   GetXYZ(x);
2569   Double_t b[3];
2570   AliTrackerBase::GetBxByBz(x,b);
2571
2572   if (!fTPCInner)  { 
2573     AliWarning("Could not get TPC Inner Param.");
2574     return fCacheChi2TPCConstrainedVsGlobal;
2575   }
2576   
2577   // clone for constraining
2578   AliExternalTrackParam* tpcInnerC = new AliExternalTrackParam(*fTPCInner);
2579   if (!tpcInnerC) { 
2580     AliWarning("Clone of TPCInnerParam failed.");
2581     return fCacheChi2TPCConstrainedVsGlobal;  
2582   }
2583   
2584   // transform to the track reference frame 
2585   Bool_t isOK = tpcInnerC->Rotate(GetAlpha());
2586   isOK &= tpcInnerC->PropagateTo(GetX(), b[2]);
2587   if (!isOK) { 
2588     delete tpcInnerC;
2589     tpcInnerC = 0; 
2590     AliWarning("Rotation/Propagation of track failed.") ; 
2591     return fCacheChi2TPCConstrainedVsGlobal;    
2592   }  
2593
2594   // constrain TPCinner 
2595   isOK = tpcInnerC->ConstrainToVertex(vtx, b);
2596   
2597   // transform to the track reference frame 
2598   isOK &= tpcInnerC->Rotate(GetAlpha());
2599   isOK &= tpcInnerC->PropagateTo(GetX(), b[2]);
2600
2601   if (!isOK) {
2602     AliWarning("ConstrainTPCInner failed.") ;
2603     delete tpcInnerC;
2604     tpcInnerC = 0; 
2605     return fCacheChi2TPCConstrainedVsGlobal;  
2606   }
2607   
2608   // calculate chi2 between vi and vj vectors
2609   // with covi and covj covariance matrices
2610   // chi2ij = (vi-vj)^(T)*(covi+covj)^(-1)*(vi-vj)
2611   TMatrixD deltaT(5,1);
2612   TMatrixD delta(1,5);
2613   TMatrixD covarM(5,5);
2614
2615   for (Int_t ipar=0; ipar<5; ipar++) {
2616     deltaT(ipar,0) = tpcInnerC->GetParameter()[ipar] - GetParameter()[ipar];
2617     delta(0,ipar) = tpcInnerC->GetParameter()[ipar] - GetParameter()[ipar];
2618
2619     for (Int_t jpar=0; jpar<5; jpar++) {
2620       Int_t index = GetIndex(ipar,jpar);
2621       covarM(ipar,jpar) = GetCovariance()[index]+tpcInnerC->GetCovariance()[index];
2622     }
2623   }
2624   // chi2 distance TPC constrained and TPC+ITS
2625   TMatrixD covarMInv = covarM.Invert();
2626   TMatrixD mat2 = covarMInv*deltaT;
2627   TMatrixD chi2 = delta*mat2; 
2628   
2629   delete tpcInnerC; 
2630   tpcInnerC = 0;
2631   
2632   fCacheChi2TPCConstrainedVsGlobal = chi2(0,0);
2633   return fCacheChi2TPCConstrainedVsGlobal;
2634 }