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