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