dded ITS chi2 on all 3 steps of tracking
[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
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: 
1266     fITSchi2Std[0] = t->GetChi2();
1267     //
1268   case kITSout: 
1269     fITSchi2Std[1] = t->GetChi2();
1270   case kITSrefit:
1271     {
1272     fITSchi2Std[2] = t->GetChi2();
1273     fITSClusterMap=0;
1274     fITSncls=t->GetNumberOfClusters();
1275     if (fFriendTrack) {
1276     Int_t* indexITS = new Int_t[AliESDfriendTrack::kMaxITScluster];
1277     for (Int_t i=0;i<AliESDfriendTrack::kMaxITScluster;i++) {
1278         indexITS[i]=t->GetClusterIndex(i);
1279
1280         if (i<fITSncls) {
1281           Int_t l=(indexITS[i] & 0xf0000000) >> 28;
1282            SETBIT(fITSClusterMap,l);                 
1283         }
1284     }
1285     fFriendTrack->SetITSIndices(indexITS,AliESDfriendTrack::kMaxITScluster);
1286     delete [] indexITS;
1287     }
1288
1289     fITSchi2=t->GetChi2();
1290     fITSsignal=t->GetPIDsignal();
1291     fITSLabel = t->GetLabel();
1292     // keep in fOp the parameters outside ITS for ITS stand-alone tracks 
1293     if (flags==kITSout) { 
1294       if (!fOp) fOp=new AliExternalTrackParam(*t);
1295       else 
1296         fOp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
1297     }   
1298     }
1299     break;
1300     
1301   case kTPCin: case kTPCrefit:
1302     {
1303     fTPCLabel = t->GetLabel();
1304     if (flags==kTPCin)  {
1305       fTPCInner=new AliExternalTrackParam(*t); 
1306       fTPCnclsIter1=t->GetNumberOfClusters();    
1307       fTPCchi2Iter1=t->GetChi2();
1308     }
1309     if (!fIp) fIp=new AliExternalTrackParam(*t);
1310     else 
1311       fIp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
1312     }
1313     // Intentionally no break statement; need to set general TPC variables as well
1314   case kTPCout:
1315     {
1316     if (flags & kTPCout){
1317       if (!fOp) fOp=new AliExternalTrackParam(*t);
1318       else 
1319         fOp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
1320     }
1321     fTPCncls=t->GetNumberOfClusters();    
1322     fTPCchi2=t->GetChi2();
1323     
1324     if (fFriendTrack) {  // Copy cluster indices
1325       Int_t* indexTPC = new Int_t[AliESDfriendTrack::kMaxTPCcluster];
1326       for (Int_t i=0;i<AliESDfriendTrack::kMaxTPCcluster;i++)         
1327         indexTPC[i]=t->GetClusterIndex(i);
1328       fFriendTrack->SetTPCIndices(indexTPC,AliESDfriendTrack::kMaxTPCcluster);
1329       delete [] indexTPC;
1330     }
1331     fTPCsignal=t->GetPIDsignal();
1332     }
1333     break;
1334
1335   case kTRDin: case kTRDrefit:
1336     break;
1337   case kTRDout:
1338     {
1339     fTRDLabel = t->GetLabel(); 
1340     fTRDchi2  = t->GetChi2();
1341     fTRDncls  = t->GetNumberOfClusters();
1342     if (fFriendTrack) {
1343       Int_t* indexTRD = new Int_t[AliESDfriendTrack::kMaxTRDcluster];
1344       for (Int_t i=0;i<AliESDfriendTrack::kMaxTRDcluster;i++) indexTRD[i]=-2;
1345       for (Int_t i=0;i<6;i++) indexTRD[i]=t->GetTrackletIndex(i);
1346       fFriendTrack->SetTRDIndices(indexTRD,AliESDfriendTrack::kMaxTRDcluster);
1347       delete [] indexTRD;
1348     }    
1349     
1350     fTRDsignal=t->GetPIDsignal();
1351     }
1352     break;
1353   case kTRDbackup:
1354     if (!fOp) fOp=new AliExternalTrackParam(*t);
1355     else 
1356       fOp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
1357     fTRDncls0 = t->GetNumberOfClusters(); 
1358     break;
1359   case kTOFin: 
1360     break;
1361   case kTOFout: 
1362     break;
1363   case kTRDStop:
1364     break;
1365   case kHMPIDout:
1366   if (!fHMPIDp) fHMPIDp=new AliExternalTrackParam(*t);
1367     else 
1368       fHMPIDp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
1369     break;
1370   default: 
1371     AliError("Wrong flag !");
1372     return kFALSE;
1373   }
1374
1375   return rc;
1376 }
1377
1378 //_______________________________________________________________________
1379 void AliESDtrack::GetExternalParameters(Double_t &x, Double_t p[5]) const {
1380   //---------------------------------------------------------------------
1381   // This function returns external representation of the track parameters
1382   //---------------------------------------------------------------------
1383   x=GetX();
1384   for (Int_t i=0; i<5; i++) p[i]=GetParameter()[i];
1385 }
1386
1387 //_______________________________________________________________________
1388 void AliESDtrack::GetExternalCovariance(Double_t cov[15]) const {
1389   //---------------------------------------------------------------------
1390   // This function returns external representation of the cov. matrix
1391   //---------------------------------------------------------------------
1392   for (Int_t i=0; i<15; i++) cov[i]=AliExternalTrackParam::GetCovariance()[i];
1393 }
1394
1395 //_______________________________________________________________________
1396 Bool_t AliESDtrack::GetConstrainedExternalParameters
1397                  (Double_t &alpha, Double_t &x, Double_t p[5]) const {
1398   //---------------------------------------------------------------------
1399   // This function returns the constrained external track parameters
1400   //---------------------------------------------------------------------
1401   if (!fCp) return kFALSE;
1402   alpha=fCp->GetAlpha();
1403   x=fCp->GetX();
1404   for (Int_t i=0; i<5; i++) p[i]=fCp->GetParameter()[i];
1405   return kTRUE;
1406 }
1407
1408 //_______________________________________________________________________
1409 Bool_t 
1410 AliESDtrack::GetConstrainedExternalCovariance(Double_t c[15]) const {
1411   //---------------------------------------------------------------------
1412   // This function returns the constrained external cov. matrix
1413   //---------------------------------------------------------------------
1414   if (!fCp) return kFALSE;
1415   for (Int_t i=0; i<15; i++) c[i]=fCp->GetCovariance()[i];
1416   return kTRUE;
1417 }
1418
1419 Bool_t
1420 AliESDtrack::GetInnerExternalParameters
1421                  (Double_t &alpha, Double_t &x, Double_t p[5]) const {
1422   //---------------------------------------------------------------------
1423   // This function returns external representation of the track parameters 
1424   // at the inner layer of TPC
1425   //---------------------------------------------------------------------
1426   if (!fIp) return kFALSE;
1427   alpha=fIp->GetAlpha();
1428   x=fIp->GetX();
1429   for (Int_t i=0; i<5; i++) p[i]=fIp->GetParameter()[i];
1430   return kTRUE;
1431 }
1432
1433 Bool_t 
1434 AliESDtrack::GetInnerExternalCovariance(Double_t cov[15]) const {
1435  //---------------------------------------------------------------------
1436  // This function returns external representation of the cov. matrix 
1437  // at the inner layer of TPC
1438  //---------------------------------------------------------------------
1439   if (!fIp) return kFALSE;
1440   for (Int_t i=0; i<15; i++) cov[i]=fIp->GetCovariance()[i];
1441   return kTRUE;
1442 }
1443
1444 void 
1445 AliESDtrack::SetOuterParam(const AliExternalTrackParam *p, ULong_t flags) {
1446   //
1447   // This is a direct setter for the outer track parameters
1448   //
1449   SetStatus(flags);
1450   if (fOp) delete fOp;
1451   fOp=new AliExternalTrackParam(*p);
1452 }
1453
1454 void 
1455 AliESDtrack::SetOuterHmpParam(const AliExternalTrackParam *p, ULong_t flags) {
1456   //
1457   // This is a direct setter for the outer track parameters
1458   //
1459   SetStatus(flags);
1460   if (fHMPIDp) delete fHMPIDp;
1461   fHMPIDp=new AliExternalTrackParam(*p);
1462 }
1463
1464 Bool_t 
1465 AliESDtrack::GetOuterExternalParameters
1466                  (Double_t &alpha, Double_t &x, Double_t p[5]) const {
1467   //---------------------------------------------------------------------
1468   // This function returns external representation of the track parameters 
1469   // at the inner layer of TRD
1470   //---------------------------------------------------------------------
1471   if (!fOp) return kFALSE;
1472   alpha=fOp->GetAlpha();
1473   x=fOp->GetX();
1474   for (Int_t i=0; i<5; i++) p[i]=fOp->GetParameter()[i];
1475   return kTRUE;
1476 }
1477
1478 Bool_t 
1479 AliESDtrack::GetOuterHmpExternalParameters
1480                  (Double_t &alpha, Double_t &x, Double_t p[5]) const {
1481   //---------------------------------------------------------------------
1482   // This function returns external representation of the track parameters 
1483   // at the inner layer of TRD
1484   //---------------------------------------------------------------------
1485   if (!fHMPIDp) return kFALSE;
1486   alpha=fHMPIDp->GetAlpha();
1487   x=fHMPIDp->GetX();
1488   for (Int_t i=0; i<5; i++) p[i]=fHMPIDp->GetParameter()[i];
1489   return kTRUE;
1490 }
1491
1492 Bool_t 
1493 AliESDtrack::GetOuterExternalCovariance(Double_t cov[15]) const {
1494  //---------------------------------------------------------------------
1495  // This function returns external representation of the cov. matrix 
1496  // at the inner layer of TRD
1497  //---------------------------------------------------------------------
1498   if (!fOp) return kFALSE;
1499   for (Int_t i=0; i<15; i++) cov[i]=fOp->GetCovariance()[i];
1500   return kTRUE;
1501 }
1502
1503 Bool_t 
1504 AliESDtrack::GetOuterHmpExternalCovariance(Double_t cov[15]) const {
1505  //---------------------------------------------------------------------
1506  // This function returns external representation of the cov. matrix 
1507  // at the inner layer of TRD
1508  //---------------------------------------------------------------------
1509   if (!fHMPIDp) return kFALSE;
1510   for (Int_t i=0; i<15; i++) cov[i]=fHMPIDp->GetCovariance()[i];
1511   return kTRUE;
1512 }
1513
1514 Int_t AliESDtrack::GetNcls(Int_t idet) const
1515 {
1516   // Get number of clusters by subdetector index
1517   //
1518   Int_t ncls = 0;
1519   switch(idet){
1520   case 0:
1521     ncls = fITSncls;
1522     break;
1523   case 1:
1524     ncls = fTPCncls;
1525     break;
1526   case 2:
1527     ncls = fTRDncls;
1528     break;
1529   case 3:
1530     if (fTOFindex != -1)
1531       ncls = 1;
1532     break;
1533   case 4: //PHOS
1534     break;
1535   case 5: //HMPID
1536     if ((fHMPIDcluIdx >= 0) && (fHMPIDcluIdx < 7000000)) {
1537       if ((fHMPIDcluIdx%1000000 != 9999) && (fHMPIDcluIdx%1000000 != 99999)) {
1538         ncls = 1;
1539       }
1540     }    
1541     break;
1542   default:
1543     break;
1544   }
1545   return ncls;
1546 }
1547
1548 Int_t AliESDtrack::GetClusters(Int_t idet, Int_t *idx) const
1549 {
1550   // Get cluster index array by subdetector index
1551   //
1552   Int_t ncls = 0;
1553   switch(idet){
1554   case 0:
1555     ncls = GetITSclusters(idx);
1556     break;
1557   case 1:
1558     ncls = GetTPCclusters(idx);
1559     break;
1560   case 2:
1561     ncls = GetTRDclusters(idx);
1562     break;
1563   case 3:
1564     if (fTOFindex != -1) {
1565       idx[0] = fTOFindex;
1566       ncls = 1;
1567     }
1568     break;
1569   case 4: //PHOS
1570     break;
1571   case 5:
1572     if ((fHMPIDcluIdx >= 0) && (fHMPIDcluIdx < 7000000)) {
1573       if ((fHMPIDcluIdx%1000000 != 9999) && (fHMPIDcluIdx%1000000 != 99999)) {
1574         idx[0] = GetHMPIDcluIdx();
1575         ncls = 1;
1576       }
1577     }    
1578     break;
1579   case 6: //EMCAL
1580     break;
1581   default:
1582     break;
1583   }
1584   return ncls;
1585 }
1586
1587 //_______________________________________________________________________
1588 void AliESDtrack::GetIntegratedTimes(Double_t *times) const {
1589   // Returns the array with integrated times for each particle hypothesis
1590   for (Int_t i=0; i<AliPID::kSPECIES; i++) times[i]=fTrackTime[i];
1591 }
1592
1593 //_______________________________________________________________________
1594 void AliESDtrack::SetIntegratedTimes(const Double_t *times) {
1595   // Sets the array with integrated times for each particle hypotesis
1596   for (Int_t i=0; i<AliPID::kSPECIES; i++) fTrackTime[i]=times[i];
1597 }
1598
1599 //_______________________________________________________________________
1600 void AliESDtrack::SetITSpid(const Double_t *p) {
1601   // Sets values for the probability of each particle type (in ITS)
1602   SetPIDValues(fITSr,p,AliPID::kSPECIES);
1603   SetStatus(AliESDtrack::kITSpid);
1604 }
1605
1606 //_______________________________________________________________________
1607 void AliESDtrack::GetITSpid(Double_t *p) const {
1608   // Gets the probability of each particle type (in ITS)
1609   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fITSr[i];
1610 }
1611
1612 //_______________________________________________________________________
1613 Char_t AliESDtrack::GetITSclusters(Int_t *idx) const {
1614   //---------------------------------------------------------------------
1615   // This function returns indices of the assgined ITS clusters 
1616   //---------------------------------------------------------------------
1617   if (idx && fFriendTrack) {
1618     Int_t *index=fFriendTrack->GetITSindices();
1619     for (Int_t i=0; i<AliESDfriendTrack::kMaxITScluster; i++) {
1620       if ( (i>=fITSncls) && (i<6) ) idx[i]=-1;
1621       else {
1622         if (index) {
1623           idx[i]=index[i];
1624         }
1625         else idx[i]= -2;
1626       }
1627     }
1628   }
1629   return fITSncls;
1630 }
1631
1632 //_______________________________________________________________________
1633 Bool_t AliESDtrack::GetITSModuleIndexInfo(Int_t ilayer,Int_t &idet,Int_t &status,
1634                                          Float_t &xloc,Float_t &zloc) const {
1635   //----------------------------------------------------------------------
1636   // This function encodes in the module number also the status of cluster association
1637   // "status" can have the following values: 
1638   // 1 "found" (cluster is associated), 
1639   // 2 "dead" (module is dead from OCDB), 
1640   // 3 "skipped" (module or layer forced to be skipped),
1641   // 4 "outinz" (track out of z acceptance), 
1642   // 5 "nocls" (no clusters in the road), 
1643   // 6 "norefit" (cluster rejected during refit), 
1644   // 7 "deadzspd" (holes in z in SPD)
1645   // Also given are the coordinates of the crossing point of track and module
1646   // (in the local module ref. system)
1647   // WARNING: THIS METHOD HAS TO BE SYNCHRONIZED WITH AliITStrackV2::GetModuleIndexInfo()!
1648   //----------------------------------------------------------------------
1649
1650   if(fITSModule[ilayer]==-1) {
1651     idet = -1;
1652     status=0;
1653     xloc=-99.; zloc=-99.;
1654     return kFALSE;
1655   }
1656
1657   Int_t module = fITSModule[ilayer];
1658
1659   idet = Int_t(module/1000000);
1660
1661   module -= idet*1000000;
1662
1663   status = Int_t(module/100000);
1664
1665   module -= status*100000;
1666
1667   Int_t signs = Int_t(module/10000);
1668
1669   module-=signs*10000;
1670
1671   Int_t xInt = Int_t(module/100);
1672   module -= xInt*100;
1673
1674   Int_t zInt = module;
1675
1676   if(signs==1) { xInt*=1; zInt*=1; }
1677   if(signs==2) { xInt*=1; zInt*=-1; }
1678   if(signs==3) { xInt*=-1; zInt*=1; }
1679   if(signs==4) { xInt*=-1; zInt*=-1; }
1680
1681   xloc = 0.1*(Float_t)xInt;
1682   zloc = 0.1*(Float_t)zInt;
1683
1684   if(status==4) idet = -1;
1685
1686   return kTRUE;
1687 }
1688
1689 //_______________________________________________________________________
1690 UShort_t AliESDtrack::GetTPCclusters(Int_t *idx) const {
1691   //---------------------------------------------------------------------
1692   // This function returns indices of the assgined ITS clusters 
1693   //---------------------------------------------------------------------
1694   if (idx && fFriendTrack) {
1695     Int_t *index=fFriendTrack->GetTPCindices();
1696
1697     if (index){
1698       for (Int_t i=0; i<AliESDfriendTrack::kMaxTPCcluster; i++) idx[i]=index[i];
1699     }
1700     else {
1701       for (Int_t i=0; i<AliESDfriendTrack::kMaxTPCcluster; i++) idx[i]=-2;
1702     }
1703   }
1704   return fTPCncls;
1705 }
1706
1707 //_______________________________________________________________________
1708 Float_t AliESDtrack::GetTPCClusterInfo(Int_t nNeighbours/*=3*/, Int_t type/*=0*/, Int_t row0, Int_t row1) const
1709 {
1710   //
1711   // TPC cluster information
1712   // type 0: get fraction of found/findable clusters with neighbourhood definition
1713   //      1: findable clusters with neighbourhood definition
1714   //      2: found clusters
1715   //
1716   // definition of findable clusters:
1717   //            a cluster is defined as findable if there is another cluster
1718   //           within +- nNeighbours pad rows. The idea is to overcome threshold
1719   //           effects with a very simple algorithm.
1720   //
1721
1722   if (type==2) return fTPCClusterMap.CountBits();
1723   
1724   Int_t found=0;
1725   Int_t findable=0;
1726   Int_t last=-nNeighbours;
1727   
1728   Int_t upperBound=fTPCClusterMap.GetNbits();
1729   if (upperBound>row1) upperBound=row1;
1730   for (Int_t i=row0; i<upperBound; ++i){
1731     //look to current row
1732     if (fTPCClusterMap[i]) {
1733       last=i;
1734       ++found;
1735       ++findable;
1736       continue;
1737     }
1738     //look to nNeighbours before
1739     if ((i-last)<=nNeighbours) {
1740       ++findable;
1741       continue;
1742     }
1743     //look to nNeighbours after
1744     for (Int_t j=i+1; j<i+1+nNeighbours; ++j){
1745       if (fTPCClusterMap[j]){
1746         ++findable;
1747         break;
1748       }
1749     }
1750   }
1751   if (type==1) return findable;
1752   
1753   if (type==0){
1754     Float_t fraction=0;
1755     if (findable>0) 
1756       fraction=(Float_t)found/(Float_t)findable;
1757     else 
1758       fraction=0;
1759     return fraction;
1760   }  
1761   return 0;  // undefined type - default value
1762 }
1763
1764 //_______________________________________________________________________
1765 Double_t AliESDtrack::GetTPCdensity(Int_t row0, Int_t row1) const{
1766   //
1767   // GetDensity of the clusters on given region between row0 and row1
1768   // Dead zone effect takin into acoount
1769   //
1770   if (!fFriendTrack) return 0.0;
1771   Int_t good  = 0;
1772   Int_t found = 0;
1773   //  
1774   Int_t *index=fFriendTrack->GetTPCindices();
1775   for (Int_t i=row0;i<=row1;i++){     
1776     Int_t idx = index[i];
1777     if (idx!=-1)  good++;             // track outside of dead zone
1778     if (idx>0)    found++;
1779   }
1780   Float_t density=0.5;
1781   if (good>(row1-row0)*0.5) density = Float_t(found)/Float_t(good);
1782   return density;
1783 }
1784
1785 //_______________________________________________________________________
1786 void AliESDtrack::SetTPCpid(const Double_t *p) {  
1787   // Sets values for the probability of each particle type (in TPC)
1788   SetPIDValues(fTPCr,p,AliPID::kSPECIES);
1789   SetStatus(AliESDtrack::kTPCpid);
1790 }
1791
1792 //_______________________________________________________________________
1793 void AliESDtrack::GetTPCpid(Double_t *p) const {
1794   // Gets the probability of each particle type (in TPC)
1795   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fTPCr[i];
1796 }
1797
1798 //_______________________________________________________________________
1799 UChar_t AliESDtrack::GetTRDclusters(Int_t *idx) const {
1800   //---------------------------------------------------------------------
1801   // This function returns indices of the assgined TRD clusters 
1802   //---------------------------------------------------------------------
1803   if (idx && fFriendTrack) {
1804     Int_t *index=fFriendTrack->GetTRDindices();
1805
1806     if (index) {
1807       for (Int_t i=0; i<AliESDfriendTrack::kMaxTRDcluster; i++) idx[i]=index[i];
1808     }
1809     else {
1810       for (Int_t i=0; i<AliESDfriendTrack::kMaxTRDcluster; i++) idx[i]=-2;
1811     }
1812   }
1813   return fTRDncls;
1814 }
1815
1816 //_______________________________________________________________________
1817 UChar_t AliESDtrack::GetTRDtracklets(Int_t *idx) const {
1818 //
1819 // This function returns the number of TRD tracklets used in tracking
1820 // and it fills the indices of these tracklets in the array "idx" as they 
1821 // are registered in the TRD track list. 
1822 // 
1823 // Caution :
1824 //   1. The idx array has to be allocated with a size >= AliESDtrack::kTRDnPlanes
1825 //   2. The idx array store not only the index but also the layer of the tracklet. 
1826 //      Therefore tracks with TRD gaps contain default values for indices [-1] 
1827
1828   if (!fFriendTrack) return 0;
1829   if (!idx) return GetTRDntracklets();
1830   Int_t *index=fFriendTrack->GetTRDindices();
1831   Int_t n = 0;
1832   for (Int_t i=0; i<kTRDnPlanes; i++){ 
1833     if (index){
1834       if(index[i]>=0) n++;
1835       idx[i]=index[i];
1836     }
1837     else idx[i] = -2;
1838   }
1839   return n;
1840 }
1841
1842 //_______________________________________________________________________
1843 void AliESDtrack::SetTRDpid(const Double_t *p) {  
1844   // Sets values for the probability of each particle type (in TRD)
1845   SetPIDValues(fTRDr,p,AliPID::kSPECIES);
1846   SetStatus(AliESDtrack::kTRDpid);
1847 }
1848
1849 //_______________________________________________________________________
1850 void AliESDtrack::GetTRDpid(Double_t *p) const {
1851   // Gets the probability of each particle type (in TRD)
1852   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fTRDr[i];
1853 }
1854
1855 //_______________________________________________________________________
1856 void    AliESDtrack::SetTRDpid(Int_t iSpecies, Float_t p)
1857 {
1858   // Sets the probability of particle type iSpecies to p (in TRD)
1859   fTRDr[iSpecies] = p;
1860 }
1861
1862 Double_t AliESDtrack::GetTRDpid(Int_t iSpecies) const
1863 {
1864   // Returns the probability of particle type iSpecies (in TRD)
1865   return fTRDr[iSpecies];
1866 }
1867
1868 //____________________________________________________
1869 Int_t AliESDtrack::GetNumberOfTRDslices() const 
1870 {
1871   // built in backward compatibility
1872   Int_t idx = fTRDnSlices - (kTRDnPlanes<<1);
1873   return idx<18 ? fTRDnSlices/kTRDnPlanes : idx/kTRDnPlanes;
1874 }
1875
1876 //____________________________________________________
1877 Double_t AliESDtrack::GetTRDmomentum(Int_t plane, Double_t *sp) const
1878 {
1879 //Returns momentum estimation and optional its error (sp)
1880 // in TRD layer "plane".
1881
1882   if (!fTRDnSlices) {
1883     AliDebug(2, "No TRD info allocated for this track.");
1884     return -1.;
1885   }
1886   if ((plane<0) || (plane>=kTRDnPlanes)) {
1887     AliWarning(Form("Request for TRD plane[%d] outside range.", plane)); 
1888     return -1.;
1889   }
1890
1891   Int_t idx = fTRDnSlices-(kTRDnPlanes<<1)+plane;
1892   // Protection for backward compatibility
1893   if(idx<(GetNumberOfTRDslices()*kTRDnPlanes)) return -1.;
1894
1895   if(sp) (*sp) = fTRDslices[idx+kTRDnPlanes];
1896   return fTRDslices[idx];
1897 }
1898
1899 //____________________________________________________
1900 Double_t  AliESDtrack::GetTRDslice(Int_t plane, Int_t slice) const {
1901   //Gets the charge from the slice of the plane
1902
1903   if(!fTRDslices) {
1904     //AliError("No TRD slices allocated for this track !");
1905     return -1.;
1906   }
1907   if ((plane<0) || (plane>=kTRDnPlanes)) {
1908     AliError("Info for TRD plane not available !");
1909     return -1.;
1910   }
1911   Int_t ns=GetNumberOfTRDslices();
1912   if ((slice<-1) || (slice>=ns)) {
1913     //AliError("Wrong TRD slice !");  
1914     return -1.;
1915   }
1916
1917   if(slice>=0) return fTRDslices[plane*ns + slice];
1918
1919   // return average of the dEdx measurements
1920   Double_t q=0.; Double32_t *s = &fTRDslices[plane*ns];
1921   for (Int_t i=0; i<ns; i++, s++) if((*s)>0.) q+=(*s);
1922   return q/ns;
1923 }
1924
1925 //____________________________________________________
1926 void  AliESDtrack::SetNumberOfTRDslices(Int_t n) {
1927   //Sets the number of slices used for PID 
1928   if (fTRDnSlices) return;
1929
1930   fTRDnSlices=n;
1931   fTRDslices=new Double32_t[fTRDnSlices];
1932   
1933   // set-up correctly the allocated memory
1934   memset(fTRDslices, 0, n*sizeof(Double32_t));
1935   for (Int_t i=GetNumberOfTRDslices(); i--;) fTRDslices[i]=-1.;
1936 }
1937
1938 //____________________________________________________
1939 void  AliESDtrack::SetTRDslice(Double_t q, Int_t plane, Int_t slice) {
1940   //Sets the charge q in the slice of the plane
1941   if(!fTRDslices) {
1942     AliError("No TRD slices allocated for this track !");
1943     return;
1944   }
1945   if ((plane<0) || (plane>=kTRDnPlanes)) {
1946     AliError("Info for TRD plane not allocated !");
1947     return;
1948   }
1949   Int_t ns=GetNumberOfTRDslices();
1950   if ((slice<0) || (slice>=ns)) {
1951     AliError("Wrong TRD slice !");
1952     return;
1953   }
1954   Int_t n=plane*ns + slice;
1955   fTRDslices[n]=q;
1956 }
1957
1958
1959 //____________________________________________________
1960 void AliESDtrack::SetTRDmomentum(Double_t p, Int_t plane, Double_t *sp)
1961 {
1962   if(!fTRDslices) {
1963     AliError("No TRD slices allocated for this track !");
1964     return;
1965   }
1966   if ((plane<0) || (plane>=kTRDnPlanes)) {
1967     AliError("Info for TRD plane not allocated !");
1968     return;
1969   }
1970
1971   Int_t idx = fTRDnSlices-(kTRDnPlanes<<1)+plane;
1972   // Protection for backward compatibility
1973   if(idx<GetNumberOfTRDslices()*kTRDnPlanes) return;
1974
1975   if(sp) fTRDslices[idx+kTRDnPlanes] = (*sp);
1976   fTRDslices[idx] = p;
1977 }
1978
1979
1980 //_______________________________________________________________________
1981 void AliESDtrack::SetTOFpid(const Double_t *p) {  
1982   // Sets the probability of each particle type (in TOF)
1983   SetPIDValues(fTOFr,p,AliPID::kSPECIES);
1984   SetStatus(AliESDtrack::kTOFpid);
1985 }
1986
1987 //_______________________________________________________________________
1988 void AliESDtrack::SetTOFLabel(const Int_t *p) {  
1989   // Sets  (in TOF)
1990   for (Int_t i=0; i<3; i++) fTOFLabel[i]=p[i];
1991 }
1992
1993 //_______________________________________________________________________
1994 void AliESDtrack::GetTOFpid(Double_t *p) const {
1995   // Gets probabilities of each particle type (in TOF)
1996   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fTOFr[i];
1997 }
1998
1999 //_______________________________________________________________________
2000 void AliESDtrack::GetTOFLabel(Int_t *p) const {
2001   // Gets (in TOF)
2002   for (Int_t i=0; i<3; i++) p[i]=fTOFLabel[i];
2003 }
2004
2005 //_______________________________________________________________________
2006 void AliESDtrack::GetTOFInfo(Float_t *info) const {
2007   // Gets (in TOF)
2008   for (Int_t i=0; i<10; i++) info[i]=fTOFInfo[i];
2009 }
2010
2011 //_______________________________________________________________________
2012 void AliESDtrack::SetTOFInfo(Float_t*info) {
2013   // Gets (in TOF)
2014   for (Int_t i=0; i<10; i++) fTOFInfo[i]=info[i];
2015 }
2016
2017
2018
2019 //_______________________________________________________________________
2020 void AliESDtrack::SetHMPIDpid(const Double_t *p) {  
2021   // Sets the probability of each particle type (in HMPID)
2022   SetPIDValues(fHMPIDr,p,AliPID::kSPECIES);
2023   SetStatus(AliESDtrack::kHMPIDpid);
2024 }
2025
2026 //_______________________________________________________________________
2027 void AliESDtrack::GetHMPIDpid(Double_t *p) const {
2028   // Gets probabilities of each particle type (in HMPID)
2029   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fHMPIDr[i];
2030 }
2031
2032
2033
2034 //_______________________________________________________________________
2035 void AliESDtrack::SetESDpid(const Double_t *p) {  
2036   // Sets the probability of each particle type for the ESD track
2037   SetPIDValues(fR,p,AliPID::kSPECIES);
2038   SetStatus(AliESDtrack::kESDpid);
2039 }
2040
2041 //_______________________________________________________________________
2042 void AliESDtrack::GetESDpid(Double_t *p) const {
2043   // Gets probability of each particle type for the ESD track
2044   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fR[i];
2045 }
2046
2047 //_______________________________________________________________________
2048 Bool_t AliESDtrack::RelateToVertexTPC(const AliESDVertex *vtx, 
2049 Double_t b, Double_t maxd, AliExternalTrackParam *cParam) {
2050   //
2051   // Try to relate the TPC-only track parameters to the vertex "vtx", 
2052   // if the (rough) transverse impact parameter is not bigger then "maxd". 
2053   //            Magnetic field is "b" (kG).
2054   //
2055   // a) The TPC-only paramters are extapolated to the DCA to the vertex.
2056   // b) The impact parameters and their covariance matrix are calculated.
2057   // c) An attempt to constrain the TPC-only params to the vertex is done.
2058   //    The constrained params are returned via "cParam".
2059   //
2060   // In the case of success, the returned value is kTRUE
2061   // otherwise, it's kFALSE)
2062   // 
2063
2064   if (!fTPCInner) return kFALSE;
2065   if (!vtx) return kFALSE;
2066
2067   Double_t dz[2],cov[3];
2068   if (!fTPCInner->PropagateToDCA(vtx, b, maxd, dz, cov)) return kFALSE;
2069
2070   fdTPC = dz[0];
2071   fzTPC = dz[1];  
2072   fCddTPC = cov[0];
2073   fCdzTPC = cov[1];
2074   fCzzTPC = cov[2];
2075   
2076   Double_t covar[6]; vtx->GetCovMatrix(covar);
2077   Double_t p[2]={GetParameter()[0]-dz[0],GetParameter()[1]-dz[1]};
2078   Double_t c[3]={covar[2],0.,covar[5]};
2079
2080   Double_t chi2=GetPredictedChi2(p,c);
2081   if (chi2>kVeryBig) return kFALSE;
2082
2083   fCchi2TPC=chi2;
2084
2085   if (!cParam) return kTRUE;
2086
2087   *cParam = *fTPCInner;
2088   if (!cParam->Update(p,c)) return kFALSE;
2089
2090   return kTRUE;
2091 }
2092
2093 //_______________________________________________________________________
2094 Bool_t AliESDtrack::RelateToVertexTPCBxByBz(const AliESDVertex *vtx, 
2095 Double_t b[3], Double_t maxd, AliExternalTrackParam *cParam) {
2096   //
2097   // Try to relate the TPC-only track parameters to the vertex "vtx", 
2098   // if the (rough) transverse impact parameter is not bigger then "maxd". 
2099   //
2100   // All three components of the magnetic field ,"b[3]" (kG), 
2101   // are taken into account.
2102   //
2103   // a) The TPC-only paramters are extapolated to the DCA to the vertex.
2104   // b) The impact parameters and their covariance matrix are calculated.
2105   // c) An attempt to constrain the TPC-only params to the vertex is done.
2106   //    The constrained params are returned via "cParam".
2107   //
2108   // In the case of success, the returned value is kTRUE
2109   // otherwise, it's kFALSE)
2110   // 
2111
2112   if (!fTPCInner) return kFALSE;
2113   if (!vtx) return kFALSE;
2114
2115   Double_t dz[2],cov[3];
2116   if (!fTPCInner->PropagateToDCABxByBz(vtx, b, maxd, dz, cov)) return kFALSE;
2117
2118   fdTPC = dz[0];
2119   fzTPC = dz[1];  
2120   fCddTPC = cov[0];
2121   fCdzTPC = cov[1];
2122   fCzzTPC = cov[2];
2123   
2124   Double_t covar[6]; vtx->GetCovMatrix(covar);
2125   Double_t p[2]={GetParameter()[0]-dz[0],GetParameter()[1]-dz[1]};
2126   Double_t c[3]={covar[2],0.,covar[5]};
2127
2128   Double_t chi2=GetPredictedChi2(p,c);
2129   if (chi2>kVeryBig) return kFALSE;
2130
2131   fCchi2TPC=chi2;
2132
2133   if (!cParam) return kTRUE;
2134
2135   *cParam = *fTPCInner;
2136   if (!cParam->Update(p,c)) return kFALSE;
2137
2138   return kTRUE;
2139 }
2140
2141 //_______________________________________________________________________
2142 Bool_t AliESDtrack::RelateToVertex(const AliESDVertex *vtx, 
2143 Double_t b, Double_t maxd, AliExternalTrackParam *cParam) {
2144   //
2145   // Try to relate this track to the vertex "vtx", 
2146   // if the (rough) transverse impact parameter is not bigger then "maxd". 
2147   //            Magnetic field is "b" (kG).
2148   //
2149   // a) The track gets extapolated to the DCA to the vertex.
2150   // b) The impact parameters and their covariance matrix are calculated.
2151   // c) An attempt to constrain this track to the vertex is done.
2152   //    The constrained params are returned via "cParam".
2153   //
2154   // In the case of success, the returned value is kTRUE
2155   // (otherwise, it's kFALSE)
2156   //  
2157
2158   if (!vtx) return kFALSE;
2159
2160   Double_t dz[2],cov[3];
2161   if (!PropagateToDCA(vtx, b, maxd, dz, cov)) return kFALSE;
2162
2163   fD = dz[0];
2164   fZ = dz[1];  
2165   fCdd = cov[0];
2166   fCdz = cov[1];
2167   fCzz = cov[2];
2168   
2169   Double_t covar[6]; vtx->GetCovMatrix(covar);
2170   Double_t p[2]={GetParameter()[0]-dz[0],GetParameter()[1]-dz[1]};
2171   Double_t c[3]={covar[2],0.,covar[5]};
2172
2173   Double_t chi2=GetPredictedChi2(p,c);
2174   if (chi2>kVeryBig) return kFALSE;
2175
2176   fCchi2=chi2;
2177
2178
2179   //--- Could now these lines be removed ? ---
2180   delete fCp;
2181   fCp=new AliExternalTrackParam(*this);  
2182
2183   if (!fCp->Update(p,c)) {delete fCp; fCp=0; return kFALSE;}
2184   //----------------------------------------
2185
2186   fVertexID = vtx->GetID();
2187
2188   if (!cParam) return kTRUE;
2189
2190   *cParam = *this;
2191   if (!cParam->Update(p,c)) return kFALSE; 
2192
2193   return kTRUE;
2194 }
2195
2196 //_______________________________________________________________________
2197 Bool_t AliESDtrack::RelateToVertexBxByBz(const AliESDVertex *vtx, 
2198 Double_t b[3], Double_t maxd, AliExternalTrackParam *cParam) {
2199   //
2200   // Try to relate this track to the vertex "vtx", 
2201   // if the (rough) transverse impact parameter is not bigger then "maxd". 
2202   //            Magnetic field is "b" (kG).
2203   //
2204   // a) The track gets extapolated to the DCA to the vertex.
2205   // b) The impact parameters and their covariance matrix are calculated.
2206   // c) An attempt to constrain this track to the vertex is done.
2207   //    The constrained params are returned via "cParam".
2208   //
2209   // In the case of success, the returned value is kTRUE
2210   // (otherwise, it's kFALSE)
2211   //  
2212
2213   if (!vtx) return kFALSE;
2214
2215   Double_t dz[2],cov[3];
2216   if (!PropagateToDCABxByBz(vtx, b, maxd, dz, cov)) return kFALSE;
2217
2218   fD = dz[0];
2219   fZ = dz[1];  
2220   fCdd = cov[0];
2221   fCdz = cov[1];
2222   fCzz = cov[2];
2223   
2224   Double_t covar[6]; vtx->GetCovMatrix(covar);
2225   Double_t p[2]={GetParameter()[0]-dz[0],GetParameter()[1]-dz[1]};
2226   Double_t c[3]={covar[2],0.,covar[5]};
2227
2228   Double_t chi2=GetPredictedChi2(p,c);
2229   if (chi2>kVeryBig) return kFALSE;
2230
2231   fCchi2=chi2;
2232
2233
2234   //--- Could now these lines be removed ? ---
2235   delete fCp;
2236   fCp=new AliExternalTrackParam(*this);  
2237
2238   if (!fCp->Update(p,c)) {delete fCp; fCp=0; return kFALSE;}
2239   //----------------------------------------
2240
2241   fVertexID = vtx->GetID();
2242
2243   if (!cParam) return kTRUE;
2244
2245   *cParam = *this;
2246   if (!cParam->Update(p,c)) return kFALSE; 
2247
2248   return kTRUE;
2249 }
2250
2251 //_______________________________________________________________________
2252 void AliESDtrack::Print(Option_t *) const {
2253   // Prints info on the track
2254   AliExternalTrackParam::Print();
2255   printf("ESD track info\n") ; 
2256   Double_t p[AliPID::kSPECIESN] ; 
2257   Int_t index = 0 ; 
2258   if( IsOn(kITSpid) ){
2259     printf("From ITS: ") ; 
2260     GetITSpid(p) ; 
2261     for(index = 0 ; index < AliPID::kSPECIES; index++) 
2262       printf("%f, ", p[index]) ;
2263     printf("\n           signal = %f\n", GetITSsignal()) ;
2264   } 
2265   if( IsOn(kTPCpid) ){
2266     printf("From TPC: ") ; 
2267     GetTPCpid(p) ; 
2268     for(index = 0 ; index < AliPID::kSPECIES; index++) 
2269       printf("%f, ", p[index]) ;
2270     printf("\n           signal = %f\n", GetTPCsignal()) ;
2271   }
2272   if( IsOn(kTRDpid) ){
2273     printf("From TRD: ") ; 
2274     GetTRDpid(p) ; 
2275     for(index = 0 ; index < AliPID::kSPECIES; index++) 
2276       printf("%f, ", p[index]) ;
2277       printf("\n           signal = %f\n", GetTRDsignal()) ;
2278   }
2279   if( IsOn(kTOFpid) ){
2280     printf("From TOF: ") ; 
2281     GetTOFpid(p) ; 
2282     for(index = 0 ; index < AliPID::kSPECIES; index++) 
2283       printf("%f, ", p[index]) ;
2284     printf("\n           signal = %f\n", GetTOFsignal()) ;
2285   }
2286   if( IsOn(kHMPIDpid) ){
2287     printf("From HMPID: ") ; 
2288     GetHMPIDpid(p) ; 
2289     for(index = 0 ; index < AliPID::kSPECIES; index++) 
2290       printf("%f, ", p[index]) ;
2291     printf("\n           signal = %f\n", GetHMPIDsignal()) ;
2292   }
2293
2294
2295
2296 //
2297 // Draw functionality
2298 // Origin: Marian Ivanov, Marian.Ivanov@cern.ch
2299 //
2300 void AliESDtrack::FillPolymarker(TPolyMarker3D *pol, Float_t magF, Float_t minR, Float_t maxR, Float_t stepR){
2301   //
2302   // Fill points in the polymarker
2303   //
2304   TObjArray arrayRef;
2305   arrayRef.AddLast(new AliExternalTrackParam(*this));
2306   if (fIp) arrayRef.AddLast(new AliExternalTrackParam(*fIp));
2307   if (fOp) arrayRef.AddLast(new AliExternalTrackParam(*fOp));
2308   if (fHMPIDp) arrayRef.AddLast(new AliExternalTrackParam(*fHMPIDp));
2309   //
2310   Double_t mpos[3]={0,0,0};
2311   Int_t entries=arrayRef.GetEntries();
2312   for (Int_t i=0;i<entries;i++){
2313     Double_t pos[3];
2314     ((AliExternalTrackParam*)arrayRef.At(i))->GetXYZ(pos);
2315     mpos[0]+=pos[0]/entries;
2316     mpos[1]+=pos[1]/entries;
2317     mpos[2]+=pos[2]/entries;    
2318   }
2319   // Rotate to the mean position
2320   //
2321   Float_t fi= TMath::ATan2(mpos[1],mpos[0]);
2322   for (Int_t i=0;i<entries;i++){
2323     Bool_t res = ((AliExternalTrackParam*)arrayRef.At(i))->Rotate(fi);
2324     if (!res) delete arrayRef.RemoveAt(i);
2325   }
2326   Int_t counter=0;
2327   for (Double_t r=minR; r<maxR; r+=stepR){
2328     Double_t sweight=0;
2329     Double_t mlpos[3]={0,0,0};
2330     for (Int_t i=0;i<entries;i++){
2331       Double_t point[3]={0,0,0};
2332       AliExternalTrackParam *param = ((AliExternalTrackParam*)arrayRef.At(i));
2333       if (!param) continue;
2334       if (param->GetXYZAt(r,magF,point)){
2335         Double_t weight = 1./(10.+(r-param->GetX())*(r-param->GetX()));
2336         sweight+=weight;
2337         mlpos[0]+=point[0]*weight;
2338         mlpos[1]+=point[1]*weight;
2339         mlpos[2]+=point[2]*weight;
2340       }
2341     }
2342     if (sweight>0){
2343       mlpos[0]/=sweight;
2344       mlpos[1]/=sweight;
2345       mlpos[2]/=sweight;      
2346       pol->SetPoint(counter,mlpos[0],mlpos[1], mlpos[2]);
2347       //      printf("xyz\t%f\t%f\t%f\n",mlpos[0], mlpos[1],mlpos[2]);
2348       counter++;
2349     }
2350   }
2351 }
2352
2353 //_______________________________________________________________________
2354 void AliESDtrack::SetITSdEdxSamples(const Double_t s[4]) {
2355   //
2356   // Store the dE/dx samples measured by the two SSD and two SDD layers.
2357   // These samples are corrected for the track segment length. 
2358   //
2359   for (Int_t i=0; i<4; i++) fITSdEdxSamples[i]=s[i];
2360 }
2361
2362 //_______________________________________________________________________
2363 void AliESDtrack::GetITSdEdxSamples(Double_t *s) const {
2364   //
2365   // Get the dE/dx samples measured by the two SSD and two SDD layers.  
2366   // These samples are corrected for the track segment length.
2367   //
2368   for (Int_t i=0; i<4; i++) s[i]=fITSdEdxSamples[i];
2369 }
2370
2371
2372 UShort_t   AliESDtrack::GetTPCnclsS(Int_t i0,Int_t i1) const{
2373   //
2374   // get number of shared TPC clusters
2375   //
2376   return  fTPCSharedMap.CountBits(i0)-fTPCSharedMap.CountBits(i1);
2377 }
2378
2379 UShort_t   AliESDtrack::GetTPCncls(Int_t i0,Int_t i1) const{
2380   //
2381   // get number of TPC clusters
2382   //
2383   return  fTPCClusterMap.CountBits(i0)-fTPCClusterMap.CountBits(i1);
2384 }