STEER/ESD/AliESDtrack.cxx

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