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