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