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