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 #include "AliDetectorPID.h"
 132 #include "TTreeStream.h"
 133 #include "TObjArray.h"
 134
 135 ClassImp(AliESDtrack)
 136
 137 void SetPIDValues(Double_t * dest, const Double_t * src, Int_t n) {
 138   // This function copies "n" PID weights from "scr" to "dest"
 139   // and normalizes their sum to 1 thus producing conditional probabilities.
 140   // The negative weights are set to 0.
 141   // In case all the weights are non-positive they are replaced by
 142   // uniform probabilities
 143
 144   if (n<=0) return;
 145
 146   Float_t uniform = 1./(Float_t)n;
 147
 148   Float_t sum = 0;
 149   for (Int_t i=0; i<n; i++)
 150     if (src[i]>=0) {
 151       sum+=src[i];
 152       dest[i] = src[i];
 153     }
 154     else {
 155       dest[i] = 0;
 156     }
 157
 158   if(sum>0)
 159     for (Int_t i=0; i<n; i++) dest[i] /= sum;
 160   else
 161     for (Int_t i=0; i<n; i++) dest[i] = uniform;
 162 }
 163
 164 //_______________________________________________________________________
 165 AliESDtrack::AliESDtrack() :
 166   AliExternalTrackParam(),
 167   fCp(0),
 168   fIp(0),
 169   fTPCInner(0),
 170   fOp(0),
 171   fHMPIDp(0),
 172   fFriendTrack(NULL),
 173   fTPCFitMap(159),//number of padrows
 174   fTPCClusterMap(159),//number of padrows
 175   fTPCSharedMap(159),//number of padrows
 176   fFlags(0),
 177   fID(0),
 178   fLabel(0),
 179   fITSLabel(0),
 180   fTPCLabel(0),
 181   fTRDLabel(0),
 182   fTOFLabel(NULL),
 183   fTOFCalChannel(-1),
 184   fTOFindex(-1),
 185   fHMPIDqn(0),
 186   fHMPIDcluIdx(-1),
 187   fCaloIndex(kEMCALNoMatch),
 188   fMassForTracking(0.13957),
 189   fHMPIDtrkTheta(0),
 190   fHMPIDtrkPhi(0),
 191   fHMPIDsignal(0),
 192   fTrackTime(0),
 193   fTrackLength(0),
 194   fdTPC(0),fzTPC(0),
 195   fCddTPC(0),fCdzTPC(0),fCzzTPC(0),
 196   fCchi2TPC(0),
 197   fD(0),fZ(0),
 198   fCdd(0),fCdz(0),fCzz(0),
 199   fCchi2(0),
 200   fITSchi2(0),
 201   fTPCchi2(0),
 202   fTPCchi2Iter1(0),
 203   fTRDchi2(0),
 204   fTOFchi2(0),
 205   fHMPIDchi2(0),
 206   fGlobalChi2(0),
 207   fITSsignal(0),
 208   fTPCsignal(0),
 209   fTPCsignalTuned(0),
 210   fTPCsignalS(0),
 211   fTPCdEdxInfo(0),
 212   fTRDsignal(0),
 213   fTRDQuality(0),
 214   fTRDBudget(0),
 215   fTOFsignal(99999),
 216   fTOFsignalTuned(99999),
 217   fTOFsignalToT(99999),
 218   fTOFsignalRaw(99999),
 219   fTOFsignalDz(999),
 220   fTOFsignalDx(999),
 221   fTOFdeltaBC(999),
 222   fTOFl0l1(999),
 223   fCaloDx(0),
 224   fCaloDz(0),
 225   fHMPIDtrkX(0),
 226   fHMPIDtrkY(0),
 227   fHMPIDmipX(0),
 228   fHMPIDmipY(0),
 229   fTPCncls(0),
 230   fTPCnclsF(0),
 231   fTPCsignalN(0),
 232   fTPCnclsIter1(0),
 233   fTPCnclsFIter1(0),
 234   fITSncls(0),
 235   fITSClusterMap(0),
 236   fITSSharedMap(0),
 237   fTRDncls(0),
 238   fTRDncls0(0),
 239   fTRDntracklets(0),
 240   fTRDNchamberdEdx(0),
 241   fTRDNclusterdEdx(0),
 242   fTRDnSlices(0),
 243   fTRDslices(0x0),
 244   fVertexID(-2),// -2 means an orphan track
 245   fESDEvent(0),
 246   fCacheNCrossedRows(-10),
 247   fCacheChi2TPCConstrainedVsGlobal(-10),
 248   fCacheChi2TPCConstrainedVsGlobalVertex(0),
 249   fDetectorPID(0x0),
 250   fTrackPhiOnEMCal(-999),
 251   fTrackEtaOnEMCal(-999),
 252   fTrackPtOnEMCal(-999),
 253   fNtofClusters(0),
 254   fTOFcluster(NULL)
 255 {
 256   //
 257   // The default ESD constructor
 258   //
 259   if (!OnlineMode()) fFriendTrack=new AliESDfriendTrack();
 260
 261   Int_t i;
 262   for (i=kNITSchi2Std;i--;) fITSchi2Std[i] = 0;
 263   for (i=0; i<AliPID::kSPECIES; i++) fR[i]=fITSr[i]=fTPCr[i]=fTRDr[i]=fTOFr[i]=fHMPIDr[i]=0.;
 264
 265   for (i=0; i<3; i++)   { fKinkIndexes[i]=0;}
 266   for (i=0; i<3; i++)   { fV0Indexes[i]=0;}
 267   for (i=0;i<kTRDnPlanes;i++) {
 268     fTRDTimBin[i]=0;
 269   }
 270   for (i=0;i<4;i++) {fITSdEdxSamples[i]=0.;}
 271   for (i=0;i<4;i++) {fTPCPoints[i]=0;}
 272   for (i=0;i<10;i++) {fTOFInfo[i]=0;}
 273   for (i=0;i<12;i++) {fITSModule[i]=-1;}
 274 }
 275
 276 bool AliESDtrack::fgkOnlineMode=false;
 277
 278 //_______________________________________________________________________
 279 AliESDtrack::AliESDtrack(const AliESDtrack& track):
 280   AliExternalTrackParam(track),
 281   fCp(0),
 282   fIp(0),
 283   fTPCInner(0),
 284   fOp(0),
 285   fHMPIDp(0),
 286   fFriendTrack(0),
 287   fTPCFitMap(track.fTPCFitMap),
 288   fTPCClusterMap(track.fTPCClusterMap),
 289   fTPCSharedMap(track.fTPCSharedMap),
 290   fFlags(track.fFlags),
 291   fID(track.fID),
 292   fLabel(track.fLabel),
 293   fITSLabel(track.fITSLabel),
 294   fTPCLabel(track.fTPCLabel),
 295   fTRDLabel(track.fTRDLabel),
 296   fTOFLabel(NULL),
 297   fTOFCalChannel(track.fTOFCalChannel),
 298   fTOFindex(track.fTOFindex),
 299   fHMPIDqn(track.fHMPIDqn),
 300   fHMPIDcluIdx(track.fHMPIDcluIdx),
 301   fCaloIndex(track.fCaloIndex),
 302   fMassForTracking(track.fMassForTracking),
 303   fHMPIDtrkTheta(track.fHMPIDtrkTheta),
 304   fHMPIDtrkPhi(track.fHMPIDtrkPhi),
 305   fHMPIDsignal(track.fHMPIDsignal),
 306   fTrackTime(NULL),
 307   fTrackLength(track.fTrackLength),
 308   fdTPC(track.fdTPC),fzTPC(track.fzTPC),
 309   fCddTPC(track.fCddTPC),fCdzTPC(track.fCdzTPC),fCzzTPC(track.fCzzTPC),
 310   fCchi2TPC(track.fCchi2TPC),
 311   fD(track.fD),fZ(track.fZ),
 312   fCdd(track.fCdd),fCdz(track.fCdz),fCzz(track.fCzz),
 313   fCchi2(track.fCchi2),
 314   fITSchi2(track.fITSchi2),
 315   fTPCchi2(track.fTPCchi2),
 316   fTPCchi2Iter1(track.fTPCchi2Iter1),
 317   fTRDchi2(track.fTRDchi2),
 318   fTOFchi2(track.fTOFchi2),
 319   fHMPIDchi2(track.fHMPIDchi2),
 320   fGlobalChi2(track.fGlobalChi2),
 321   fITSsignal(track.fITSsignal),
 322   fTPCsignal(track.fTPCsignal),
 323   fTPCsignalTuned(track.fTPCsignalTuned),
 324   fTPCsignalS(track.fTPCsignalS),
 325   fTPCdEdxInfo(0),
 326   fTRDsignal(track.fTRDsignal),
 327   fTRDQuality(track.fTRDQuality),
 328   fTRDBudget(track.fTRDBudget),
 329   fTOFsignal(track.fTOFsignal),
 330   fTOFsignalTuned(track.fTOFsignalTuned),
 331   fTOFsignalToT(track.fTOFsignalToT),
 332   fTOFsignalRaw(track.fTOFsignalRaw),
 333   fTOFsignalDz(track.fTOFsignalDz),
 334   fTOFsignalDx(track.fTOFsignalDx),
 335   fTOFdeltaBC(track.fTOFdeltaBC),
 336   fTOFl0l1(track.fTOFl0l1),
 337   fCaloDx(track.fCaloDx),
 338   fCaloDz(track.fCaloDz),
 339   fHMPIDtrkX(track.fHMPIDtrkX),
 340   fHMPIDtrkY(track.fHMPIDtrkY),
 341   fHMPIDmipX(track.fHMPIDmipX),
 342   fHMPIDmipY(track.fHMPIDmipY),
 343   fTPCncls(track.fTPCncls),
 344   fTPCnclsF(track.fTPCnclsF),
 345   fTPCsignalN(track.fTPCsignalN),
 346   fTPCnclsIter1(track.fTPCnclsIter1),
 347   fTPCnclsFIter1(track.fTPCnclsIter1),
 348   fITSncls(track.fITSncls),
 349   fITSClusterMap(track.fITSClusterMap),
 350   fITSSharedMap(track.fITSSharedMap),
 351   fTRDncls(track.fTRDncls),
 352   fTRDncls0(track.fTRDncls0),
 353   fTRDntracklets(track.fTRDntracklets),
 354   fTRDNchamberdEdx(track.fTRDNchamberdEdx),
 355   fTRDNclusterdEdx(track.fTRDNclusterdEdx),
 356   fTRDnSlices(track.fTRDnSlices),
 357   fTRDslices(0x0),
 358   fVertexID(track.fVertexID),
 359   fESDEvent(track.fESDEvent),
 360   fCacheNCrossedRows(track.fCacheNCrossedRows),
 361   fCacheChi2TPCConstrainedVsGlobal(track.fCacheChi2TPCConstrainedVsGlobal),
 362   fCacheChi2TPCConstrainedVsGlobalVertex(track.fCacheChi2TPCConstrainedVsGlobalVertex),
 363   fDetectorPID(0x0),
 364   fTrackPhiOnEMCal(track.fTrackPhiOnEMCal),
 365   fTrackEtaOnEMCal(track.fTrackEtaOnEMCal),
 366   fTrackPtOnEMCal(track.fTrackPtOnEMCal),
 367   fNtofClusters(track.fNtofClusters),
 368   fTOFcluster(NULL)
 369 {
 370   //
 371   //copy constructor
 372   //
 373   for (Int_t i=kNITSchi2Std;i--;) fITSchi2Std[i] = track.fITSchi2Std[i];
 374
 375   if(track.fTrackTime){
 376     fTrackTime = new Double32_t[AliPID::kSPECIESC];
 377     for (Int_t i=0;i<AliPID::kSPECIESC;i++) fTrackTime[i]=track.fTrackTime[i];
 378   }
 379
 380   for (Int_t i=0;i<AliPID::kSPECIES;i++)  fR[i]=track.fR[i];
 381   //
 382   for (Int_t i=0;i<AliPID::kSPECIES;i++) fITSr[i]=track.fITSr[i];
 383   //
 384   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTPCr[i]=track.fTPCr[i];
 385   for (Int_t i=0;i<4;i++) {fITSdEdxSamples[i]=track.fITSdEdxSamples[i];}
 386   for (Int_t i=0;i<4;i++) {fTPCPoints[i]=track.fTPCPoints[i];}
 387   for (Int_t i=0; i<3;i++)   { fKinkIndexes[i]=track.fKinkIndexes[i];}
 388   for (Int_t i=0; i<3;i++)   { fV0Indexes[i]=track.fV0Indexes[i];}
 389   //
 390   for (Int_t i=0;i<kTRDnPlanes;i++) {
 391     fTRDTimBin[i]=track.fTRDTimBin[i];
 392   }
 393
 394   if (fTRDnSlices) {
 395     fTRDslices=new Double32_t[fTRDnSlices];
 396     for (Int_t i=0; i<fTRDnSlices; i++) fTRDslices[i]=track.fTRDslices[i];
 397   }
 398
 399   if (track.fDetectorPID) fDetectorPID = new AliDetectorPID(*track.fDetectorPID);
 400
 401   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTRDr[i]=track.fTRDr[i];
 402   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTOFr[i]=track.fTOFr[i];
 403   if(track.fTOFLabel){
 404     if(!fTOFLabel) fTOFLabel = new Int_t[3];
 405     for (Int_t i=0;i<3;i++) fTOFLabel[i]=track.fTOFLabel[i];
 406   }
 407
 408   for (Int_t i=0;i<10;i++) fTOFInfo[i]=track.fTOFInfo[i];
 409   for (Int_t i=0;i<12;i++) fITSModule[i]=track.fITSModule[i];
 410   for (Int_t i=0;i<AliPID::kSPECIES;i++) fHMPIDr[i]=track.fHMPIDr[i];
 411
 412   if (track.fCp) fCp=new AliExternalTrackParam(*track.fCp);
 413   if (track.fIp) fIp=new AliExternalTrackParam(*track.fIp);
 414   if (track.fTPCInner) fTPCInner=new AliExternalTrackParam(*track.fTPCInner);
 415   if (track.fOp) fOp=new AliExternalTrackParam(*track.fOp);
 416   if (track.fHMPIDp) fHMPIDp=new AliExternalTrackParam(*track.fHMPIDp);
 417   if (track.fTPCdEdxInfo) fTPCdEdxInfo = new AliTPCdEdxInfo(*track.fTPCdEdxInfo);
 418
 419
 420   if (track.fFriendTrack) fFriendTrack=new AliESDfriendTrack(*(track.fFriendTrack));
 421
 422   if(fNtofClusters > 0){
 423     fTOFcluster = new Int_t[fNtofClusters];
 424         for(Int_t i=0;i < fNtofClusters;i++) fTOFcluster[i] = track.fTOFcluster[i];
 425   }
 426 }
 427
 428 //_______________________________________________________________________
 429 AliESDtrack::AliESDtrack(const AliVTrack *track) :
 430   AliExternalTrackParam(track),
 431   fCp(0),
 432   fIp(0),
 433   fTPCInner(0),
 434   fOp(0),
 435   fHMPIDp(0),
 436   fFriendTrack(0),
 437   fTPCFitMap(159),//number of padrows
 438   fTPCClusterMap(159),//number of padrows
 439   fTPCSharedMap(159),//number of padrows
 440   fFlags(0),
 441   fID(),
 442   fLabel(0),
 443   fITSLabel(0),
 444   fTPCLabel(0),
 445   fTRDLabel(0),
 446   fTOFLabel(NULL),
 447   fTOFCalChannel(-1),
 448   fTOFindex(-1),
 449   fHMPIDqn(0),
 450   fHMPIDcluIdx(-1),
 451   fCaloIndex(kEMCALNoMatch),
 452   fMassForTracking(0.13957),
 453   fHMPIDtrkTheta(0),
 454   fHMPIDtrkPhi(0),
 455   fHMPIDsignal(0),
 456   fTrackTime(NULL),
 457   fTrackLength(0),
 458   fdTPC(0),fzTPC(0),
 459   fCddTPC(0),fCdzTPC(0),fCzzTPC(0),
 460   fCchi2TPC(0),
 461   fD(0),fZ(0),
 462   fCdd(0),fCdz(0),fCzz(0),
 463   fCchi2(0),
 464   fITSchi2(0),
 465   fTPCchi2(0),
 466   fTPCchi2Iter1(0),
 467   fTRDchi2(0),
 468   fTOFchi2(0),
 469   fHMPIDchi2(0),
 470   fGlobalChi2(0),
 471   fITSsignal(0),
 472   fTPCsignal(0),
 473   fTPCsignalTuned(0),
 474   fTPCsignalS(0),
 475   fTPCdEdxInfo(0),
 476   fTRDsignal(0),
 477   fTRDQuality(0),
 478   fTRDBudget(0),
 479   fTOFsignal(99999),
 480   fTOFsignalTuned(99999),
 481   fTOFsignalToT(99999),
 482   fTOFsignalRaw(99999),
 483   fTOFsignalDz(999),
 484   fTOFsignalDx(999),
 485   fTOFdeltaBC(999),
 486   fTOFl0l1(999),
 487   fCaloDx(0),
 488   fCaloDz(0),
 489   fHMPIDtrkX(0),
 490   fHMPIDtrkY(0),
 491   fHMPIDmipX(0),
 492   fHMPIDmipY(0),
 493   fTPCncls(0),
 494   fTPCnclsF(0),
 495   fTPCsignalN(0),
 496   fTPCnclsIter1(0),
 497   fTPCnclsFIter1(0),
 498   fITSncls(0),
 499   fITSClusterMap(0),
 500   fITSSharedMap(0),
 501   fTRDncls(0),
 502   fTRDncls0(0),
 503   fTRDntracklets(0),
 504   fTRDNchamberdEdx(0),
 505   fTRDNclusterdEdx(0),
 506   fTRDnSlices(0),
 507   fTRDslices(0x0),
 508   fVertexID(-2),  // -2 means an orphan track
 509   fESDEvent(0),
 510   fCacheNCrossedRows(-10),
 511   fCacheChi2TPCConstrainedVsGlobal(-10),
 512   fCacheChi2TPCConstrainedVsGlobalVertex(0),
 513   fDetectorPID(0x0),
 514   fTrackPhiOnEMCal(-999),
 515   fTrackEtaOnEMCal(-999),
 516   fTrackPtOnEMCal(-999),
 517   fNtofClusters(0),
 518   fTOFcluster(NULL)
 519 {
 520   //
 521   // ESD track from AliVTrack.
 522   // This is not a copy constructor !
 523   //
 524
 525   if (track->InheritsFrom("AliExternalTrackParam")) {
 526      AliError("This is not a copy constructor. Use AliESDtrack(const AliESDtrack &) !");
 527      AliWarning("Calling the default constructor...");
 528      AliESDtrack();
 529      return;
 530   }
 531
 532   // Reset all the arrays
 533   Int_t i;
 534   for (i=kNITSchi2Std;i--;) fITSchi2Std[i] = 0;
 535   for (i=0; i<AliPID::kSPECIES; i++) fR[i]=fITSr[i]=fTPCr[i]=fTRDr[i]=fTOFr[i]=fHMPIDr[i]=0.;
 536
 537   for (i=0; i<3; i++)   { fKinkIndexes[i]=0;}
 538   for (i=0; i<3; i++)   { fV0Indexes[i]=-1;}
 539   for (i=0;i<kTRDnPlanes;i++) {
 540     fTRDTimBin[i]=0;
 541   }
 542   for (i=0;i<4;i++) {fITSdEdxSamples[i]=0.;}
 543   for (i=0;i<4;i++) {fTPCPoints[i]=0;}
 544   for (i=0;i<10;i++) {fTOFInfo[i]=0;}
 545   for (i=0;i<12;i++) {fITSModule[i]=-1;}
 546
 547   // Set the ID
 548   SetID(track->GetID());
 549
 550   // Set ITS cluster map
 551   fITSClusterMap=track->GetITSClusterMap();
 552   fITSSharedMap=0;
 553
 554   fITSncls=0;
 555   for(i=0; i<6; i++) {
 556     if(HasPointOnITSLayer(i)) fITSncls++;
 557   }
 558
 559   // Set TPC ncls
 560   fTPCncls=track->GetTPCNcls();
 561   fTPCnclsF=track->GetTPCNclsF();
 562   // TPC cluster maps
 563   const TBits* bmap = track->GetTPCClusterMapPtr();
 564   if (bmap) SetTPCClusterMap(*bmap);
 565   bmap = GetTPCFitMapPtr();
 566   if (bmap) SetTPCFitMap(*bmap);
 567   bmap = GetTPCSharedMapPtr();
 568   if (bmap) SetTPCSharedMap(*bmap);
 569   //
 570   // Set the combined PID
 571   const Double_t *pid = track->PID();
 572   if(pid) for (i=0; i<AliPID::kSPECIES; i++) fR[i]=pid[i];
 573   //
 574   // calo matched cluster id
 575   SetEMCALcluster(track->GetEMCALcluster());
 576   // AliESD track label
 577   //
 578   // PID info
 579   fITSsignal = track->GetITSsignal();
 580   double itsdEdx[4];
 581   track->GetITSdEdxSamples(itsdEdx);
 582   SetITSdEdxSamples(itsdEdx);
 583   //
 584   SetTPCsignal(track->GetTPCsignal(),fTPCsignalS,track->GetTPCsignalN()); // No signalS in AODPi
 585   AliTPCdEdxInfo * dEdxInfo = track->GetTPCdEdxInfo();
 586   if (dEdxInfo) SetTPCdEdxInfo(new AliTPCdEdxInfo(*dEdxInfo));
 587   //
 588   SetTRDsignal(track->GetTRDsignal());
 589   int ntrdsl = track->GetNumberOfTRDslices();
 590   if (ntrdsl>0) {
 591     SetNumberOfTRDslices((ntrdsl+2)*kTRDnPlanes);
 592     for (int ipl=kTRDnPlanes;ipl--;){
 593       for (int isl=ntrdsl;isl--;) SetTRDslice(track->GetTRDslice(ipl,isl),ipl,isl);
 594       Double_t sp, p = track->GetTRDmomentum(ipl, &sp);
 595       SetTRDmomentum(p, ipl, &sp);
 596     }
 597   }
 598   //
 599   fTRDncls = track->GetTRDncls();
 600   fTRDntracklets &= 0xff & track->GetTRDntrackletsPID();
 601   fTRDchi2 = track->GetTRDchi2();
 602   //
 603   SetTOFsignal(track->GetTOFsignal());
 604   Double_t expt[AliPID::kSPECIESC];
 605   track->GetIntegratedTimes(expt);
 606   SetIntegratedTimes(expt);
 607   //
 608   SetTrackPhiEtaPtOnEMCal(track->GetTrackPhiOnEMCal(),track->GetTrackEtaOnEMCal(),track->GetTrackPtOnEMCal());
 609   //
 610   SetLabel(track->GetLabel());
 611   // Set the status
 612   SetStatus(track->GetStatus());
 613 }
 614
 615 //_______________________________________________________________________
 616 AliESDtrack::AliESDtrack(TParticle * part) :
 617   AliExternalTrackParam(),
 618   fCp(0),
 619   fIp(0),
 620   fTPCInner(0),
 621   fOp(0),
 622   fHMPIDp(0),
 623   fFriendTrack(0),
 624   fTPCFitMap(159),//number of padrows
 625   fTPCClusterMap(159),//number of padrows
 626   fTPCSharedMap(159),//number of padrows
 627   fFlags(0),
 628   fID(0),
 629   fLabel(0),
 630   fITSLabel(0),
 631   fTPCLabel(0),
 632   fTRDLabel(0),
 633   fTOFLabel(NULL),
 634   fTOFCalChannel(-1),
 635   fTOFindex(-1),
 636   fHMPIDqn(0),
 637   fHMPIDcluIdx(-1),
 638   fCaloIndex(kEMCALNoMatch),
 639   fMassForTracking(0.13957),
 640   fHMPIDtrkTheta(0),
 641   fHMPIDtrkPhi(0),
 642   fHMPIDsignal(0),
 643   fTrackTime(NULL),
 644   fTrackLength(0),
 645   fdTPC(0),fzTPC(0),
 646   fCddTPC(0),fCdzTPC(0),fCzzTPC(0),
 647   fCchi2TPC(0),
 648   fD(0),fZ(0),
 649   fCdd(0),fCdz(0),fCzz(0),
 650   fCchi2(0),
 651   fITSchi2(0),
 652   fTPCchi2(0),
 653   fTPCchi2Iter1(0),
 654   fTRDchi2(0),
 655   fTOFchi2(0),
 656   fHMPIDchi2(0),
 657   fGlobalChi2(0),
 658   fITSsignal(0),
 659   fTPCsignal(0),
 660   fTPCsignalTuned(0),
 661   fTPCsignalS(0),
 662   fTPCdEdxInfo(0),
 663   fTRDsignal(0),
 664   fTRDQuality(0),
 665   fTRDBudget(0),
 666   fTOFsignal(99999),
 667   fTOFsignalTuned(99999),
 668   fTOFsignalToT(99999),
 669   fTOFsignalRaw(99999),
 670   fTOFsignalDz(999),
 671   fTOFsignalDx(999),
 672   fTOFdeltaBC(999),
 673   fTOFl0l1(999),
 674   fCaloDx(0),
 675   fCaloDz(0),
 676   fHMPIDtrkX(0),
 677   fHMPIDtrkY(0),
 678   fHMPIDmipX(0),
 679   fHMPIDmipY(0),
 680   fTPCncls(0),
 681   fTPCnclsF(0),
 682   fTPCsignalN(0),
 683   fTPCnclsIter1(0),
 684   fTPCnclsFIter1(0),
 685   fITSncls(0),
 686   fITSClusterMap(0),
 687   fITSSharedMap(0),
 688   fTRDncls(0),
 689   fTRDncls0(0),
 690   fTRDntracklets(0),
 691   fTRDNchamberdEdx(0),
 692   fTRDNclusterdEdx(0),
 693   fTRDnSlices(0),
 694   fTRDslices(0x0),
 695   fVertexID(-2),  // -2 means an orphan track
 696   fESDEvent(0),
 697   fCacheNCrossedRows(-10),
 698   fCacheChi2TPCConstrainedVsGlobal(-10),
 699   fCacheChi2TPCConstrainedVsGlobalVertex(0),
 700   fDetectorPID(0x0),
 701   fTrackPhiOnEMCal(-999),
 702   fTrackEtaOnEMCal(-999),
 703   fTrackPtOnEMCal(-999),
 704   fNtofClusters(0),
 705   fTOFcluster(NULL)
 706 {
 707   //
 708   // ESD track from TParticle
 709   //
 710
 711   // Reset all the arrays
 712   Int_t i;
 713   for (i=kNITSchi2Std;i--;) fITSchi2Std[i] = 0;
 714   for (i=0; i<AliPID::kSPECIES; i++) fR[i]=fITSr[i]=fTPCr[i]=fTRDr[i]=fTOFr[i]=fHMPIDr[i]=0.;
 715
 716   for (i=0; i<3; i++)   { fKinkIndexes[i]=0;}
 717   for (i=0; i<3; i++)   { fV0Indexes[i]=-1;}
 718   for (i=0;i<kTRDnPlanes;i++) {
 719     fTRDTimBin[i]=0;
 720   }
 721   for (i=0;i<4;i++) {fITSdEdxSamples[i]=0.;}
 722   for (i=0;i<4;i++) {fTPCPoints[i]=0;}
 723   for (i=0;i<10;i++) {fTOFInfo[i]=0;}
 724   for (i=0;i<12;i++) {fITSModule[i]=-1;}
 725
 726   // Calculate the AliExternalTrackParam content
 727
 728   Double_t xref;
 729   Double_t alpha;
 730   Double_t param[5];
 731   Double_t covar[15];
 732
 733   // Calculate alpha: the rotation angle of the corresponding local system (TPC sector)
 734   alpha = part->Phi()*180./TMath::Pi();
 735   if (alpha<0) alpha+= 360.;
 736   if (alpha>360) alpha -= 360.;
 737
 738   Int_t sector = (Int_t)(alpha/20.);
 739   alpha = 10. + 20.*sector;
 740   alpha /= 180;
 741   alpha *= TMath::Pi();
 742
 743   // Covariance matrix: no errors, the parameters are exact
 744   for (i=0; i<15; i++) covar[i]=0.;
 745
 746   // Get the vertex of origin and the momentum
 747   TVector3 ver(part->Vx(),part->Vy(),part->Vz());
 748   TVector3 mom(part->Px(),part->Py(),part->Pz());
 749
 750   // Rotate to the local coordinate system (TPC sector)
 751   ver.RotateZ(-alpha);
 752   mom.RotateZ(-alpha);
 753
 754   // X of the referense plane
 755   xref = ver.X();
 756
 757   Int_t pdgCode = part->GetPdgCode();
 758
 759   Double_t charge =
 760     TDatabasePDG::Instance()->GetParticle(pdgCode)->Charge();
 761
 762   param[0] = ver.Y();
 763   param[1] = ver.Z();
 764   param[2] = TMath::Sin(mom.Phi());
 765   param[3] = mom.Pz()/mom.Pt();
 766   param[4] = TMath::Sign(1/mom.Pt(),charge);
 767
 768   // Set AliExternalTrackParam
 769   Set(xref, alpha, param, covar);
 770
 771   // Set the PID
 772   Int_t indexPID = 99;
 773   if (pdgCode<0) pdgCode = -pdgCode;
 774   for (i=0;i<AliPID::kSPECIES;i++) if (pdgCode==AliPID::ParticleCode(i)) {indexPID = i; break;}
 775
 776   // If the particle is not valid charged the PID probabilities are set to 0
 777   if (indexPID < AliPID::kSPECIES) {
 778     fR[indexPID]=1.;
 779     fITSr[indexPID]=1.;
 780     fTPCr[indexPID]=1.;
 781     fTRDr[indexPID]=1.;
 782     fTOFr[indexPID]=1.;
 783     fHMPIDr[indexPID]=1.;
 784
 785   }
 786   // AliESD track label
 787   SetLabel(part->GetUniqueID());
 788
 789 }
 790
 791 //_______________________________________________________________________
 792 AliESDtrack::~AliESDtrack(){
 793   //
 794   // This is destructor according Coding Conventrions
 795   //
 796   //printf("Delete track\n");
 797   delete fIp;
 798   delete fTPCInner;
 799   delete fOp;
 800   delete fHMPIDp;
 801   delete fCp;
 802   delete fFriendTrack;
 803   delete fTPCdEdxInfo;
 804   if(fTRDnSlices)
 805     delete[] fTRDslices;
 806
 807   //Reset cached values - needed for TClonesArray in AliESDInputHandler
 808   fCacheNCrossedRows = -10.;
 809   fCacheChi2TPCConstrainedVsGlobal = -10.;
 810   if(fCacheChi2TPCConstrainedVsGlobalVertex) fCacheChi2TPCConstrainedVsGlobalVertex = 0;
 811
 812   if(fTOFcluster)
 813     delete[] fTOFcluster;
 814   fTOFcluster = NULL;
 815   fNtofClusters=0;
 816
 817   delete fDetectorPID;
 818
 819   if(fTrackTime) delete[] fTrackTime;
 820   if(fTOFLabel) delete[] fTOFLabel;
 821 }
 822
 823 AliESDtrack &AliESDtrack::operator=(const AliESDtrack &source){
 824
 825
 826   if(&source == this) return *this;
 827   AliExternalTrackParam::operator=(source);
 828
 829
 830   if(source.fCp){
 831     // we have the trackparam: assign or copy construct
 832     if(fCp)*fCp = *source.fCp;
 833     else fCp = new AliExternalTrackParam(*source.fCp);
 834   }
 835   else{
 836     // no track param delete the old one
 837     delete fCp;
 838     fCp = 0;
 839   }
 840
 841   if(source.fIp){
 842     // we have the trackparam: assign or copy construct
 843     if(fIp)*fIp = *source.fIp;
 844     else fIp = new AliExternalTrackParam(*source.fIp);
 845   }
 846   else{
 847     // no track param delete the old one
 848     delete fIp;
 849     fIp = 0;
 850   }
 851
 852
 853   if(source.fTPCInner){
 854     // we have the trackparam: assign or copy construct
 855     if(fTPCInner) *fTPCInner = *source.fTPCInner;
 856     else fTPCInner = new AliExternalTrackParam(*source.fTPCInner);
 857   }
 858   else{
 859     // no track param delete the old one
 860     delete fTPCInner;
 861     fTPCInner = 0;
 862   }
 863
 864   if(source.fTPCdEdxInfo) {
 865     if(fTPCdEdxInfo) *fTPCdEdxInfo = *source.fTPCdEdxInfo;
 866     fTPCdEdxInfo = new AliTPCdEdxInfo(*source.fTPCdEdxInfo);
 867   }
 868
 869   if(source.fOp){
 870     // we have the trackparam: assign or copy construct
 871     if(fOp) *fOp = *source.fOp;
 872     else fOp = new AliExternalTrackParam(*source.fOp);
 873   }
 874   else{
 875     // no track param delete the old one
 876     delete fOp;
 877     fOp = 0;
 878   }
 879
 880
 881   if(source.fHMPIDp){
 882     // we have the trackparam: assign or copy construct
 883     if(fHMPIDp) *fHMPIDp = *source.fHMPIDp;
 884     else fHMPIDp = new AliExternalTrackParam(*source.fHMPIDp);
 885   }
 886   else{
 887     // no track param delete the old one
 888     delete fHMPIDp;
 889     fHMPIDp = 0;
 890   }
 891
 892   // copy also the friend track
 893   // use copy constructor
 894   if(source.fFriendTrack){
 895     // we have the trackparam: assign or copy construct
 896     delete fFriendTrack; fFriendTrack=new AliESDfriendTrack(*source.fFriendTrack);
 897   }
 898   else{
 899     // no track param delete the old one
 900     delete fFriendTrack; fFriendTrack= 0;
 901   }
 902
 903   fTPCFitMap = source.fTPCFitMap;
 904   fTPCClusterMap = source.fTPCClusterMap;
 905   fTPCSharedMap  = source.fTPCSharedMap;
 906   // the simple stuff
 907   fFlags    = source.fFlags;
 908   fID       = source.fID;
 909   fLabel    = source.fLabel;
 910   fITSLabel = source.fITSLabel;
 911   for(int i = 0; i< 12;++i){
 912     fITSModule[i] = source.fITSModule[i];
 913   }
 914   fTPCLabel = source.fTPCLabel;
 915   fTRDLabel = source.fTRDLabel;
 916   if(source.fTOFLabel){
 917     if(!fTOFLabel) fTOFLabel = new Int_t[3];
 918     for(int i = 0; i< 3;++i){
 919       fTOFLabel[i] = source.fTOFLabel[i];
 920     }
 921   }
 922   fTOFCalChannel = source.fTOFCalChannel;
 923   fTOFindex      = source.fTOFindex;
 924   fHMPIDqn       = source.fHMPIDqn;
 925   fHMPIDcluIdx   = source.fHMPIDcluIdx;
 926   fCaloIndex    = source.fCaloIndex;
 927   for (int i=kNITSchi2Std;i--;) fITSchi2Std[i] = source.fITSchi2Std[i];
 928   for(int i = 0; i< 3;++i){
 929     fKinkIndexes[i] = source.fKinkIndexes[i];
 930     fV0Indexes[i]   = source.fV0Indexes[i];
 931   }
 932
 933   for(int i = 0; i< AliPID::kSPECIES;++i){
 934     fR[i]     = source.fR[i];
 935     fITSr[i]  = source.fITSr[i];
 936     fTPCr[i]  = source.fTPCr[i];
 937     fTRDr[i]  = source.fTRDr[i];
 938     fTOFr[i]  = source.fTOFr[i];
 939     fHMPIDr[i] = source.fHMPIDr[i];
 940   }
 941
 942   fMassForTracking = source.fMassForTracking;
 943
 944   fHMPIDtrkTheta = source.fHMPIDtrkTheta;
 945   fHMPIDtrkPhi   = source.fHMPIDtrkPhi;
 946   fHMPIDsignal   = source.fHMPIDsignal;
 947
 948
 949   if(fTrackTime){
 950     delete[] fTrackTime;
 951   }
 952   if(source.fTrackTime){
 953     fTrackTime = new Double32_t[AliPID::kSPECIESC];
 954     for(Int_t i=0;i < AliPID::kSPECIESC;i++)
 955       fTrackTime[i] = source.fTrackTime[i];
 956   }
 957
 958   fTrackLength   = source. fTrackLength;
 959   fdTPC  = source.fdTPC;
 960   fzTPC  = source.fzTPC;
 961   fCddTPC = source.fCddTPC;
 962   fCdzTPC = source.fCdzTPC;
 963   fCzzTPC = source.fCzzTPC;
 964   fCchi2TPC = source.fCchi2TPC;
 965
 966   fD  = source.fD;
 967   fZ  = source.fZ;
 968   fCdd = source.fCdd;
 969   fCdz = source.fCdz;
 970   fCzz = source.fCzz;
 971   fCchi2     = source.fCchi2;
 972
 973   fITSchi2   = source.fITSchi2;
 974   fTPCchi2   = source.fTPCchi2;
 975   fTPCchi2Iter1   = source.fTPCchi2Iter1;
 976   fTRDchi2   = source.fTRDchi2;
 977   fTOFchi2   = source.fTOFchi2;
 978   fHMPIDchi2 = source.fHMPIDchi2;
 979
 980   fGlobalChi2 = source.fGlobalChi2;
 981
 982   fITSsignal  = source.fITSsignal;
 983   for (Int_t i=0;i<4;i++) {fITSdEdxSamples[i]=source.fITSdEdxSamples[i];}
 984   fTPCsignal  = source.fTPCsignal;
 985   fTPCsignalTuned  = source.fTPCsignalTuned;
 986   fTPCsignalS = source.fTPCsignalS;
 987   for(int i = 0; i< 4;++i){
 988     fTPCPoints[i] = source.fTPCPoints[i];
 989   }
 990   fTRDsignal = source.fTRDsignal;
 991   fTRDNchamberdEdx = source.fTRDNchamberdEdx;
 992   fTRDNclusterdEdx = source.fTRDNclusterdEdx;
 993
 994   for(int i = 0;i < kTRDnPlanes;++i){
 995     fTRDTimBin[i] = source.fTRDTimBin[i];
 996   }
 997
 998   if(fTRDnSlices)
 999     delete[] fTRDslices;
1000   fTRDslices=0;
1001   fTRDnSlices=source.fTRDnSlices;
1002   if (fTRDnSlices) {
1003     fTRDslices=new Double32_t[fTRDnSlices];
1004     for(int j = 0;j < fTRDnSlices;++j) fTRDslices[j] = source.fTRDslices[j];
1005   }
1006
1007   fTRDQuality =   source.fTRDQuality;
1008   fTRDBudget  =   source.fTRDBudget;
1009   fTOFsignal  =   source.fTOFsignal;
1010   fTOFsignalTuned  = source.fTOFsignalTuned;
1011   fTOFsignalToT = source.fTOFsignalToT;
1012   fTOFsignalRaw = source.fTOFsignalRaw;
1013   fTOFsignalDz  = source.fTOFsignalDz;
1014   fTOFsignalDx  = source.fTOFsignalDx;
1015   fTOFdeltaBC   = source.fTOFdeltaBC;
1016   fTOFl0l1      = source.fTOFl0l1;
1017
1018   for(int i = 0;i<10;++i){
1019     fTOFInfo[i] = source.fTOFInfo[i];
1020   }
1021
1022   fHMPIDtrkX = source.fHMPIDtrkX;
1023   fHMPIDtrkY = source.fHMPIDtrkY;
1024   fHMPIDmipX = source.fHMPIDmipX;
1025   fHMPIDmipY = source.fHMPIDmipY;
1026
1027   fTPCncls    = source.fTPCncls;
1028   fTPCnclsF   = source.fTPCnclsF;
1029   fTPCsignalN = source.fTPCsignalN;
1030   fTPCnclsIter1    = source.fTPCnclsIter1;
1031   fTPCnclsFIter1   = source.fTPCnclsFIter1;
1032
1033   fITSncls = source.fITSncls;
1034   fITSClusterMap = source.fITSClusterMap;
1035   fITSSharedMap = source.fITSSharedMap;
1036   fTRDncls   = source.fTRDncls;
1037   fTRDncls0  = source.fTRDncls0;
1038   fTRDntracklets  = source.fTRDntracklets;
1039   fVertexID = source.fVertexID;
1040
1041   fCacheNCrossedRows = source.fCacheNCrossedRows;
1042   fCacheChi2TPCConstrainedVsGlobal = source.fCacheChi2TPCConstrainedVsGlobal;
1043   fCacheChi2TPCConstrainedVsGlobalVertex = source.fCacheChi2TPCConstrainedVsGlobalVertex;
1044
1045   delete fDetectorPID;
1046   fDetectorPID=0x0;
1047   if (source.fDetectorPID) fDetectorPID = new AliDetectorPID(*source.fDetectorPID);
1048
1049   fTrackPhiOnEMCal= source.fTrackPhiOnEMCal;
1050   fTrackEtaOnEMCal= source.fTrackEtaOnEMCal;
1051   fTrackPtOnEMCal= source.fTrackPtOnEMCal;
1052
1053   if(fTOFcluster){
1054     delete[] fTOFcluster;
1055   }
1056   fNtofClusters = source.fNtofClusters;
1057   if(fNtofClusters > 0){
1058     fTOFcluster = new Int_t[fNtofClusters];
1059         for(Int_t i=0;i < fNtofClusters;i++) fTOFcluster[i] = source.fTOFcluster[i];
1060   }
1061
1062   return *this;
1063 }
1064
1065
1066
1067 void AliESDtrack::Copy(TObject &obj) const {
1068
1069   // this overwrites the virtual TOBject::Copy()
1070   // to allow run time copying without casting
1071   // in AliESDEvent
1072
1073   if(this==&obj)return;
1074   AliESDtrack *robj = dynamic_cast<AliESDtrack*>(&obj);
1075   if(!robj)return; // not an AliESDtrack
1076   *robj = *this;
1077
1078 }
1079
1080
1081
1082 void AliESDtrack::AddCalibObject(TObject * object){
1083   //
1084   // add calib object to the list
1085   //
1086   if (!fFriendTrack) fFriendTrack  = new AliESDfriendTrack;
1087   if (!fFriendTrack) return;
1088   fFriendTrack->AddCalibObject(object);
1089 }
1090
1091 TObject *  AliESDtrack::GetCalibObject(Int_t index){
1092   //
1093   // return calib objct at given position
1094   //
1095   if (!fFriendTrack) return 0;
1096   return fFriendTrack->GetCalibObject(index);
1097 }
1098
1099
1100 Bool_t AliESDtrack::FillTPCOnlyTrack(AliESDtrack &track){
1101
1102   // Fills the information of the TPC-only first reconstruction pass
1103   // into the passed ESDtrack object. For consistency fTPCInner is also filled
1104   // again
1105
1106
1107
1108   // For data produced before r26675
1109   // RelateToVertexTPC was not properly called during reco
1110   // so you'll have to call it again, before FillTPCOnlyTrack
1111   //  Float_t p[2],cov[3];
1112   // track->GetImpactParametersTPC(p,cov);
1113   // if(p[0]==0&&p[1]==0) // <- Default values
1114   //  track->RelateToVertexTPC(esd->GetPrimaryVertexTPC(),esd->GetMagneticField(),kVeryBig);
1115
1116
1117   if(!fTPCInner)return kFALSE;
1118
1119   // fill the TPC track params to the global track parameters
1120   track.Set(fTPCInner->GetX(),fTPCInner->GetAlpha(),fTPCInner->GetParameter(),fTPCInner->GetCovariance());
1121   track.fD = fdTPC;
1122   track.fZ = fzTPC;
1123   track.fCdd = fCddTPC;
1124   track.fCdz = fCdzTPC;
1125   track.fCzz = fCzzTPC;
1126
1127   // copy the inner params
1128   if(track.fIp) *track.fIp = *fIp;
1129   else track.fIp = new AliExternalTrackParam(*fIp);
1130
1131   // copy the TPCinner parameters
1132   if(track.fTPCInner) *track.fTPCInner = *fTPCInner;
1133   else track.fTPCInner = new AliExternalTrackParam(*fTPCInner);
1134   track.fdTPC   = fdTPC;
1135   track.fzTPC   = fzTPC;
1136   track.fCddTPC = fCddTPC;
1137   track.fCdzTPC = fCdzTPC;
1138   track.fCzzTPC = fCzzTPC;
1139   track.fCchi2TPC = fCchi2TPC;
1140
1141   // copy all other TPC specific parameters
1142
1143   // replace label by TPC label
1144   track.fLabel    = fTPCLabel;
1145   track.fTPCLabel = fTPCLabel;
1146
1147   track.fTPCchi2 = fTPCchi2;
1148   track.fTPCchi2Iter1 = fTPCchi2Iter1;
1149   track.fTPCsignal = fTPCsignal;
1150   track.fTPCsignalTuned = fTPCsignalTuned;
1151   track.fTPCsignalS = fTPCsignalS;
1152   for(int i = 0;i<4;++i)track.fTPCPoints[i] = fTPCPoints[i];
1153
1154   track.fTPCncls    = fTPCncls;
1155   track.fTPCnclsF   = fTPCnclsF;
1156   track.fTPCsignalN =  fTPCsignalN;
1157   track.fTPCnclsIter1    = fTPCnclsIter1;
1158   track.fTPCnclsFIter1   = fTPCnclsFIter1;
1159
1160   // PID
1161   for(int i=0;i<AliPID::kSPECIES;++i){
1162     track.fTPCr[i] = fTPCr[i];
1163     // combined PID is TPC only!
1164     track.fR[i] = fTPCr[i];
1165   }
1166   track.fTPCFitMap = fTPCFitMap;
1167   track.fTPCClusterMap = fTPCClusterMap;
1168   track.fTPCSharedMap = fTPCSharedMap;
1169
1170
1171   // reset the flags
1172   track.fFlags = kTPCin;
1173   track.fID    = fID;
1174
1175   track.fFlags |= fFlags & kTPCpid; //copy the TPCpid status flag
1176
1177   for (Int_t i=0;i<3;i++) track.fKinkIndexes[i] = fKinkIndexes[i];
1178
1179   return kTRUE;
1180
1181 }
1182
1183 //_______________________________________________________________________
1184 void AliESDtrack::MakeMiniESDtrack(){
1185   // Resets everything except
1186   // fFlags: Reconstruction status flags
1187   // fLabel: Track label
1188   // fID:  Unique ID of the track
1189   // Impact parameter information
1190   // fR[AliPID::kSPECIES]: combined "detector response probability"
1191   // Running track parameters in the base class (AliExternalTrackParam)
1192
1193   fTrackLength = 0;
1194
1195   if(fTrackTime)
1196     for (Int_t i=0;i<AliPID::kSPECIESC;i++) fTrackTime[i] = 0;
1197
1198   // Reset track parameters constrained to the primary vertex
1199   delete fCp;fCp = 0;
1200
1201   // Reset track parameters at the inner wall of TPC
1202   delete fIp;fIp = 0;
1203   delete fTPCInner;fTPCInner=0;
1204   // Reset track parameters at the inner wall of the TRD
1205   delete fOp;fOp = 0;
1206   // Reset track parameters at the HMPID
1207   delete fHMPIDp;fHMPIDp = 0;
1208
1209
1210   // Reset ITS track related information
1211   fITSchi2 = 0;
1212   fITSncls = 0;
1213   fITSClusterMap=0;
1214   fITSSharedMap=0;
1215   fITSsignal = 0;
1216   for (Int_t i=0;i<4;i++) fITSdEdxSamples[i] = 0.;
1217   for (Int_t i=0;i<AliPID::kSPECIES;i++) fITSr[i]=0;
1218   fITSLabel = 0;
1219
1220   // Reset TPC related track information
1221   fTPCchi2 = 0;
1222   fTPCchi2Iter1 = 0;
1223   fTPCncls = 0;
1224   fTPCnclsF = 0;
1225   fTPCnclsIter1 = 0;
1226   fTPCnclsFIter1 = 0;
1227   fTPCFitMap = 0;
1228   fTPCClusterMap = 0;
1229   fTPCSharedMap = 0;
1230   fTPCsignal= 0;
1231   fTPCsignalTuned= 0;
1232   fTPCsignalS= 0;
1233   fTPCsignalN= 0;
1234   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTPCr[i]=0;
1235   fTPCLabel=0;
1236   for (Int_t i=0;i<4;i++) fTPCPoints[i] = 0;
1237   for (Int_t i=0; i<3;i++)   fKinkIndexes[i] = 0;
1238   for (Int_t i=0; i<3;i++)   fV0Indexes[i] = 0;
1239
1240   // Reset TRD related track information
1241   fTRDchi2 = 0;
1242   fTRDncls = 0;
1243   fTRDncls0 = 0;
1244   fTRDsignal = 0;
1245   fTRDNchamberdEdx = 0;
1246   fTRDNclusterdEdx = 0;
1247
1248   for (Int_t i=0;i<kTRDnPlanes;i++) {
1249     fTRDTimBin[i]  = 0;
1250   }
1251   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTRDr[i] = 0;
1252   fTRDLabel = 0;
1253   fTRDQuality  = 0;
1254   fTRDntracklets = 0;
1255   if(fTRDnSlices)
1256     delete[] fTRDslices;
1257   fTRDslices=0x0;
1258   fTRDnSlices=0;
1259   fTRDBudget  = 0;
1260
1261   // Reset TOF related track information
1262   fTOFchi2 = 0;
1263   fTOFindex = -1;
1264   fTOFsignal = 99999;
1265   fTOFCalChannel = -1;
1266   fTOFsignalToT = 99999;
1267   fTOFsignalRaw = 99999;
1268   fTOFsignalDz = 999;
1269   fTOFsignalDx = 999;
1270   fTOFdeltaBC = 999;
1271   fTOFl0l1 = 999;
1272   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTOFr[i] = 0;
1273   for (Int_t i=0;i<10;i++) fTOFInfo[i] = 0;
1274
1275   // Reset HMPID related track information
1276   fHMPIDchi2 = 0;
1277   fHMPIDqn = 0;
1278   fHMPIDcluIdx = -1;
1279   fHMPIDsignal = 0;
1280   for (Int_t i=0;i<AliPID::kSPECIES;i++) fHMPIDr[i] = 0;
1281   fMassForTracking = 0.13957;
1282   fHMPIDtrkTheta = 0;
1283   fHMPIDtrkPhi = 0;
1284   fHMPIDtrkX = 0;
1285   fHMPIDtrkY = 0;
1286   fHMPIDmipX = 0;
1287   fHMPIDmipY = 0;
1288   fCaloIndex = kEMCALNoMatch;
1289
1290   // reset global track chi2
1291   fGlobalChi2 = 0;
1292
1293   fVertexID = -2; // an orphan track
1294
1295   delete fFriendTrack; fFriendTrack = 0;
1296 }
1297
1298 //_______________________________________________________________________
1299 Int_t AliESDtrack::GetPID(Bool_t tpcOnly) const
1300 {
1301   // Returns the particle most probable id
1302   Int_t i;
1303   const Double32_t *prob = 0;
1304   if (tpcOnly) { // check if TPCpid is valid
1305     prob = fTPCr;
1306     for (i=0; i<AliPID::kSPECIES-1; i++) if (prob[i] != prob[i+1]) break;
1307     if (i == AliPID::kSPECIES-1) prob = 0; // not valid, try with combined pid
1308   }
1309   if (!prob) { // either requested TPCpid is not valid or comb.pid is requested
1310     prob = fR;
1311     for (i=0; i<AliPID::kSPECIES-1; i++) if (prob[i] != prob[i+1]) break;
1312     if (i == AliPID::kSPECIES-1) return AliPID::kPion;  // If all the probabilities are equal, return the pion mass
1313   }
1314   //
1315   Float_t max=0.;
1316   Int_t k=-1;
1317   for (i=0; i<AliPID::kSPECIES; i++) if (prob[i]>max) {k=i; max=prob[i];}
1318   //
1319   if (k==0) { // dE/dx "crossing points" in the TPC
1320     Double_t p=GetP();
1321     if ((p>0.38)&&(p<0.48))
1322       if (prob[0]<prob[3]*10.) return AliPID::kKaon;
1323     if ((p>0.75)&&(p<0.85))
1324       if (prob[0]<prob[4]*10.) return AliPID::kProton;
1325     return AliPID::kElectron;
1326   }
1327   if (k==1) return AliPID::kMuon;
1328   if (k==2||k==-1) return AliPID::kPion;
1329   if (k==3) return AliPID::kKaon;
1330   if (k==4) return AliPID::kProton;
1331   AliWarning("Undefined PID !");
1332   return AliPID::kPion;
1333 }
1334
1335 //_______________________________________________________________________
1336 Int_t AliESDtrack::GetTOFBunchCrossing(Double_t b, Bool_t pidTPConly) const
1337 {
1338   // Returns the number of bunch crossings after trigger (assuming 25ns spacing)
1339   const double kSpacing = 25e3; // min interbanch spacing
1340   const double kShift = 0;
1341   Int_t bcid = kTOFBCNA; // defualt one
1342   if (!IsOn(kTOFout) || !IsOn(kESDpid)) return bcid; // no info
1343   //
1344   double tdif = fTOFsignal;
1345   Double_t times[AliPID::kSPECIESC];
1346   GetIntegratedTimes(times);
1347   if (IsOn(kTIME)) { // integrated time info is there
1348     int pid = GetPID(pidTPConly);
1349
1350     tdif -= times[pid];
1351   }
1352   else { // assume integrated time info from TOF radius and momentum
1353     const double kRTOF = 385.;
1354     const double kCSpeed = 3.e-2; // cm/ps
1355     double p = GetP();
1356     if (p<0.01) return bcid;
1357     double m = GetMass(pidTPConly);
1358     double curv = GetC(b);
1359     double path = TMath::Abs(curv)>kAlmost0 ? // account for curvature
1360       2./curv*TMath::ASin(kRTOF*curv/2.)*TMath::Sqrt(1.+GetTgl()*GetTgl()) : kRTOF;
1361     tdif -= path/kCSpeed*TMath::Sqrt(1.+m*m/(p*p));
1362   }
1363   bcid = TMath::Nint((tdif - kShift)/kSpacing);
1364   return bcid;
1365 }
1366
1367 //______________________________________________________________________________
1368 Double_t AliESDtrack::M() const
1369 {
1370   // Returns the assumed mass
1371   // (the pion mass, if the particle can't be identified properly).
1372   static Bool_t printerr=kTRUE;
1373   if (printerr) {
1374      AliWarning("WARNING !!! ... THIS WILL BE PRINTED JUST ONCE !!!");
1375      printerr = kFALSE;
1376      AliWarning("This is the ESD mass. Use it with care !");
1377   }
1378   return GetMass();
1379 }
1380
1381 //______________________________________________________________________________
1382 Double_t AliESDtrack::E() const
1383 {
1384   // Returns the energy of the particle given its assumed mass.
1385   // Assumes the pion mass if the particle can't be identified properly.
1386
1387   Double_t m = M();
1388   Double_t p = P();
1389   return TMath::Sqrt(p*p + m*m);
1390 }
1391
1392 //______________________________________________________________________________
1393 Double_t AliESDtrack::Y() const
1394 {
1395   // Returns the rapidity of a particle given its assumed mass.
1396   // Assumes the pion mass if the particle can't be identified properly.
1397
1398   Double_t e = E();
1399   Double_t pz = Pz();
1400   if (e != TMath::Abs(pz)) { // energy was not equal to pz
1401     return 0.5*TMath::Log((e+pz)/(e-pz));
1402   } else { // energy was equal to pz
1403     return -999.;
1404   }
1405 }
1406
1407 //_______________________________________________________________________
1408 Bool_t AliESDtrack::UpdateTrackParams(const AliKalmanTrack *t, ULong_t flags){
1409   //
1410   // This function updates track's running parameters
1411   //
1412   Bool_t rc=kTRUE;
1413
1414   SetStatus(flags);
1415   fLabel=t->GetLabel();
1416
1417   if (t->IsStartedTimeIntegral()) {
1418     SetStatus(kTIME);
1419     Double_t times[AliPID::kSPECIESC];t->GetIntegratedTimes(times); SetIntegratedTimes(times);
1420     SetIntegratedLength(t->GetIntegratedLength());
1421   }
1422
1423   Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
1424   if (fFriendTrack) {
1425   if (flags==kITSout) fFriendTrack->SetITSOut(*t);
1426   if (flags==kTPCout) fFriendTrack->SetTPCOut(*t);
1427   if (flags==kTRDrefit) fFriendTrack->SetTRDIn(*t);
1428   }
1429
1430   switch (flags) {
1431
1432   case kITSin:
1433     fITSchi2Std[0] = t->GetChi2();
1434     //
1435   case kITSout:
1436     fITSchi2Std[1] = t->GetChi2();
1437   case kITSrefit:
1438     {
1439     fITSchi2Std[2] = t->GetChi2();
1440     fITSClusterMap=0;
1441     fITSncls=t->GetNumberOfClusters();
1442     if (fFriendTrack) {
1443     Int_t* indexITS = new Int_t[AliESDfriendTrack::kMaxITScluster];
1444     for (Int_t i=0;i<AliESDfriendTrack::kMaxITScluster;i++) {
1445         indexITS[i]=t->GetClusterIndex(i);
1446
1447         if (i<fITSncls) {
1448           Int_t l=(indexITS[i] & 0xf0000000) >> 28;
1449            SETBIT(fITSClusterMap,l);
1450         }
1451     }
1452     fFriendTrack->SetITSIndices(indexITS,AliESDfriendTrack::kMaxITScluster);
1453     delete [] indexITS;
1454     }
1455
1456     fITSchi2=t->GetChi2();
1457     fITSsignal=t->GetPIDsignal();
1458     fITSLabel = t->GetLabel();
1459     // keep in fOp the parameters outside ITS for ITS stand-alone tracks
1460     if (flags==kITSout) {
1461       if (!fOp) fOp=new AliExternalTrackParam(*t);
1462       else
1463         fOp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
1464     }
1465     }
1466     break;
1467
1468   case kTPCin: case kTPCrefit:
1469     {
1470     fTPCLabel = t->GetLabel();
1471     if (flags==kTPCin)  {
1472       fTPCInner=new AliExternalTrackParam(*t);
1473       fTPCnclsIter1=t->GetNumberOfClusters();
1474       fTPCchi2Iter1=t->GetChi2();
1475     }
1476     if (!fIp) fIp=new AliExternalTrackParam(*t);
1477     else
1478       fIp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
1479     }
1480     // Intentionally no break statement; need to set general TPC variables as well
1481   case kTPCout:
1482     {
1483     if (flags & kTPCout){
1484       if (!fOp) fOp=new AliExternalTrackParam(*t);
1485       else
1486         fOp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
1487     }
1488     fTPCncls=t->GetNumberOfClusters();
1489     fTPCchi2=t->GetChi2();
1490
1491     if (fFriendTrack) {  // Copy cluster indices
1492       Int_t* indexTPC = new Int_t[AliESDfriendTrack::kMaxTPCcluster];
1493       for (Int_t i=0;i<AliESDfriendTrack::kMaxTPCcluster;i++)
1494         indexTPC[i]=t->GetClusterIndex(i);
1495       fFriendTrack->SetTPCIndices(indexTPC,AliESDfriendTrack::kMaxTPCcluster);
1496       delete [] indexTPC;
1497     }
1498     fTPCsignal=t->GetPIDsignal();
1499     }
1500     break;
1501
1502   case kTRDin: case kTRDrefit:
1503     break;
1504   case kTRDout:
1505     {
1506     fTRDLabel = t->GetLabel();
1507     fTRDchi2  = t->GetChi2();
1508     fTRDncls  = t->GetNumberOfClusters();
1509     if (fFriendTrack) {
1510       Int_t* indexTRD = new Int_t[AliESDfriendTrack::kMaxTRDcluster];
1511       for (Int_t i=0;i<AliESDfriendTrack::kMaxTRDcluster;i++) indexTRD[i]=-2;
1512       for (Int_t i=0;i<6;i++) indexTRD[i]=t->GetTrackletIndex(i);
1513       fFriendTrack->SetTRDIndices(indexTRD,AliESDfriendTrack::kMaxTRDcluster);
1514       delete [] indexTRD;
1515     }
1516
1517     //commented out by Xianguo
1518     //fTRDsignal=t->GetPIDsignal();
1519     }
1520     break;
1521   case kTRDbackup:
1522     if (!fOp) fOp=new AliExternalTrackParam(*t);
1523     else
1524       fOp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
1525     fTRDncls0 = t->GetNumberOfClusters();
1526     break;
1527   case kTOFin:
1528     break;
1529   case kTOFout:
1530     break;
1531   case kTRDStop:
1532     break;
1533   case kHMPIDout:
1534   if (!fHMPIDp) fHMPIDp=new AliExternalTrackParam(*t);
1535     else
1536       fHMPIDp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
1537     break;
1538   default:
1539     AliError("Wrong flag !");
1540     return kFALSE;
1541   }
1542
1543   return rc;
1544 }
1545
1546 //_______________________________________________________________________
1547 void AliESDtrack::GetExternalParameters(Double_t &x, Double_t p[5]) const {
1548   //---------------------------------------------------------------------
1549   // This function returns external representation of the track parameters
1550   //---------------------------------------------------------------------
1551   x=GetX();
1552   for (Int_t i=0; i<5; i++) p[i]=GetParameter()[i];
1553 }
1554
1555 //_______________________________________________________________________
1556 void AliESDtrack::GetExternalCovariance(Double_t cov[15]) const {
1557   //---------------------------------------------------------------------
1558   // This function returns external representation of the cov. matrix
1559   //---------------------------------------------------------------------
1560   for (Int_t i=0; i<15; i++) cov[i]=AliExternalTrackParam::GetCovariance()[i];
1561 }
1562
1563 //_______________________________________________________________________
1564 Bool_t AliESDtrack::GetConstrainedExternalParameters
1565                  (Double_t &alpha, Double_t &x, Double_t p[5]) const {
1566   //---------------------------------------------------------------------
1567   // This function returns the constrained external track parameters
1568   //---------------------------------------------------------------------
1569   if (!fCp) return kFALSE;
1570   alpha=fCp->GetAlpha();
1571   x=fCp->GetX();
1572   for (Int_t i=0; i<5; i++) p[i]=fCp->GetParameter()[i];
1573   return kTRUE;
1574 }
1575
1576 //_______________________________________________________________________
1577 Bool_t
1578 AliESDtrack::GetConstrainedExternalCovariance(Double_t c[15]) const {
1579   //---------------------------------------------------------------------
1580   // This function returns the constrained external cov. matrix
1581   //---------------------------------------------------------------------
1582   if (!fCp) return kFALSE;
1583   for (Int_t i=0; i<15; i++) c[i]=fCp->GetCovariance()[i];
1584   return kTRUE;
1585 }
1586
1587 Bool_t
1588 AliESDtrack::GetInnerExternalParameters
1589                  (Double_t &alpha, Double_t &x, Double_t p[5]) const {
1590   //---------------------------------------------------------------------
1591   // This function returns external representation of the track parameters
1592   // at the inner layer of TPC
1593   //---------------------------------------------------------------------
1594   if (!fIp) return kFALSE;
1595   alpha=fIp->GetAlpha();
1596   x=fIp->GetX();
1597   for (Int_t i=0; i<5; i++) p[i]=fIp->GetParameter()[i];
1598   return kTRUE;
1599 }
1600
1601 Bool_t
1602 AliESDtrack::GetInnerExternalCovariance(Double_t cov[15]) const {
1603  //---------------------------------------------------------------------
1604  // This function returns external representation of the cov. matrix
1605  // at the inner layer of TPC
1606  //---------------------------------------------------------------------
1607   if (!fIp) return kFALSE;
1608   for (Int_t i=0; i<15; i++) cov[i]=fIp->GetCovariance()[i];
1609   return kTRUE;
1610 }
1611
1612 void
1613 AliESDtrack::SetOuterParam(const AliExternalTrackParam *p, ULong_t flags) {
1614   //
1615   // This is a direct setter for the outer track parameters
1616   //
1617   SetStatus(flags);
1618   if (fOp) delete fOp;
1619   fOp=new AliExternalTrackParam(*p);
1620 }
1621
1622 void
1623 AliESDtrack::SetOuterHmpParam(const AliExternalTrackParam *p, ULong_t flags) {
1624   //
1625   // This is a direct setter for the outer track parameters
1626   //
1627   SetStatus(flags);
1628   if (fHMPIDp) delete fHMPIDp;
1629   fHMPIDp=new AliExternalTrackParam(*p);
1630 }
1631
1632 Bool_t
1633 AliESDtrack::GetOuterExternalParameters
1634                  (Double_t &alpha, Double_t &x, Double_t p[5]) const {
1635   //---------------------------------------------------------------------
1636   // This function returns external representation of the track parameters
1637   // at the inner layer of TRD
1638   //---------------------------------------------------------------------
1639   if (!fOp) return kFALSE;
1640   alpha=fOp->GetAlpha();
1641   x=fOp->GetX();
1642   for (Int_t i=0; i<5; i++) p[i]=fOp->GetParameter()[i];
1643   return kTRUE;
1644 }
1645
1646 Bool_t
1647 AliESDtrack::GetOuterHmpExternalParameters
1648                  (Double_t &alpha, Double_t &x, Double_t p[5]) const {
1649   //---------------------------------------------------------------------
1650   // This function returns external representation of the track parameters
1651   // at the inner layer of TRD
1652   //---------------------------------------------------------------------
1653   if (!fHMPIDp) return kFALSE;
1654   alpha=fHMPIDp->GetAlpha();
1655   x=fHMPIDp->GetX();
1656   for (Int_t i=0; i<5; i++) p[i]=fHMPIDp->GetParameter()[i];
1657   return kTRUE;
1658 }
1659
1660 Bool_t
1661 AliESDtrack::GetOuterExternalCovariance(Double_t cov[15]) const {
1662  //---------------------------------------------------------------------
1663  // This function returns external representation of the cov. matrix
1664  // at the inner layer of TRD
1665  //---------------------------------------------------------------------
1666   if (!fOp) return kFALSE;
1667   for (Int_t i=0; i<15; i++) cov[i]=fOp->GetCovariance()[i];
1668   return kTRUE;
1669 }
1670
1671 Bool_t
1672 AliESDtrack::GetOuterHmpExternalCovariance(Double_t cov[15]) const {
1673  //---------------------------------------------------------------------
1674  // This function returns external representation of the cov. matrix
1675  // at the inner layer of TRD
1676  //---------------------------------------------------------------------
1677   if (!fHMPIDp) return kFALSE;
1678   for (Int_t i=0; i<15; i++) cov[i]=fHMPIDp->GetCovariance()[i];
1679   return kTRUE;
1680 }
1681
1682 Int_t AliESDtrack::GetNcls(Int_t idet) const
1683 {
1684   // Get number of clusters by subdetector index
1685   //
1686   Int_t ncls = 0;
1687   switch(idet){
1688   case 0:
1689     ncls = fITSncls;
1690     break;
1691   case 1:
1692     ncls = fTPCncls;
1693     break;
1694   case 2:
1695     ncls = fTRDncls;
1696     break;
1697   case 3:
1698     if (fTOFindex != -1)
1699       ncls = 1;
1700     break;
1701   case 4: //PHOS
1702     break;
1703   case 5: //HMPID
1704     if ((fHMPIDcluIdx >= 0) && (fHMPIDcluIdx < 7000000)) {
1705       if ((fHMPIDcluIdx%1000000 != 9999) && (fHMPIDcluIdx%1000000 != 99999)) {
1706         ncls = 1;
1707       }
1708     }
1709     break;
1710   default:
1711     break;
1712   }
1713   return ncls;
1714 }
1715
1716 Int_t AliESDtrack::GetClusters(Int_t idet, Int_t *idx) const
1717 {
1718   // Get cluster index array by subdetector index
1719   //
1720   Int_t ncls = 0;
1721   switch(idet){
1722   case 0:
1723     ncls = GetITSclusters(idx);
1724     break;
1725   case 1:
1726     ncls = GetTPCclusters(idx);
1727     break;
1728   case 2:
1729     ncls = GetTRDclusters(idx);
1730     break;
1731   case 3:
1732     if (fTOFindex != -1) {
1733       idx[0] = fTOFindex;
1734       ncls = 1;
1735     }
1736     break;
1737   case 4: //PHOS
1738     break;
1739   case 5:
1740     if ((fHMPIDcluIdx >= 0) && (fHMPIDcluIdx < 7000000)) {
1741       if ((fHMPIDcluIdx%1000000 != 9999) && (fHMPIDcluIdx%1000000 != 99999)) {
1742         idx[0] = GetHMPIDcluIdx();
1743         ncls = 1;
1744       }
1745     }
1746     break;
1747   case 6: //EMCAL
1748     break;
1749   default:
1750     break;
1751   }
1752   return ncls;
1753 }
1754
1755 //_______________________________________________________________________
1756 void AliESDtrack::GetIntegratedTimes(Double_t *times) const {
1757   Int_t index = -1;
1758
1759   if(fNtofClusters>0 && GetESDEvent()){
1760     TObjArray *tofclArray = GetESDEvent()->GetTOFcluster();
1761     AliESDTOFcluster *tofcl = (AliESDTOFcluster *) tofclArray->At(fTOFcluster[0]);
1762
1763     for(Int_t i=0;i < tofcl->GetNMatchableTracks();i++){
1764       if(tofcl->GetTrackIndex(i) == GetID()) index = i;
1765     }
1766     if(fNtofClusters>0 && index > -1){
1767       for (Int_t i=0; i<AliPID::kSPECIESC; i++) times[i]=tofcl->GetIntegratedTime(i,index);
1768       return;
1769     }
1770   }
1771   else if(fNtofClusters>0)
1772     AliInfo("No AliESDEvent available here!\n");
1773
1774   // Returns the array with integrated times for each particle hypothesis
1775   if(fTrackTime)
1776     for (Int_t i=0; i<AliPID::kSPECIESC; i++) times[i]=fTrackTime[i];
1777   else
1778     for (Int_t i=0; i<AliPID::kSPECIESC; i++) times[i]=0.0;
1779 }
1780 //_______________________________________________________________________
1781 Double_t AliESDtrack::GetIntegratedLength() const{
1782   Int_t index = -1;
1783   if(fNtofClusters>0 && GetESDEvent()){
1784     TObjArray *tofclArray = GetESDEvent()->GetTOFcluster();
1785     AliESDTOFcluster *tofcl = (AliESDTOFcluster *) tofclArray->At(fTOFcluster[0]);
1786
1787     for(Int_t i=0;i < tofcl->GetNMatchableTracks();i++){
1788       if(tofcl->GetTrackIndex(i) == GetID()) index = i;
1789     }
1790
1791     if(fNtofClusters>0 && index > -1)
1792       return tofcl->GetLength(index);
1793   }
1794   else if(fNtofClusters>0) AliInfo("No AliESDEvent available here!\n");
1795
1796   return fTrackLength;
1797 }
1798
1799 //_______________________________________________________________________
1800 void AliESDtrack::SetIntegratedTimes(const Double_t *times) {
1801   // Sets the array with integrated times for each particle hypotesis
1802   if(!fTrackTime)
1803     fTrackTime = new Double32_t[AliPID::kSPECIESC];
1804   for (Int_t i=0; i<AliPID::kSPECIESC; i++) fTrackTime[i]=times[i];
1805 }
1806
1807 //_______________________________________________________________________
1808 void AliESDtrack::SetITSpid(const Double_t *p) {
1809   // Sets values for the probability of each particle type (in ITS)
1810   SetPIDValues(fITSr,p,AliPID::kSPECIES);
1811   SetStatus(AliESDtrack::kITSpid);
1812 }
1813
1814 //_______________________________________________________________________
1815 void AliESDtrack::GetITSpid(Double_t *p) const {
1816   // Gets the probability of each particle type (in ITS)
1817   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fITSr[i];
1818 }
1819
1820 //_______________________________________________________________________
1821 Char_t AliESDtrack::GetITSclusters(Int_t *idx) const {
1822   //---------------------------------------------------------------------
1823   // This function returns indices of the assgined ITS clusters
1824   //---------------------------------------------------------------------
1825   if (idx && fFriendTrack) {
1826     Int_t *index=fFriendTrack->GetITSindices();
1827     for (Int_t i=0; i<AliESDfriendTrack::kMaxITScluster; i++) {
1828       if ( (i>=fITSncls) && (i<6) ) idx[i]=-1;
1829       else {
1830         if (index) {
1831           idx[i]=index[i];
1832         }
1833         else idx[i]= -2;
1834       }
1835     }
1836   }
1837   return fITSncls;
1838 }
1839
1840 //_______________________________________________________________________
1841 Bool_t AliESDtrack::GetITSModuleIndexInfo(Int_t ilayer,Int_t &idet,Int_t &status,
1842                                          Float_t &xloc,Float_t &zloc) const {
1843   //----------------------------------------------------------------------
1844   // This function encodes in the module number also the status of cluster association
1845   // "status" can have the following values:
1846   // 1 "found" (cluster is associated),
1847   // 2 "dead" (module is dead from OCDB),
1848   // 3 "skipped" (module or layer forced to be skipped),
1849   // 4 "outinz" (track out of z acceptance),
1850   // 5 "nocls" (no clusters in the road),
1851   // 6 "norefit" (cluster rejected during refit),
1852   // 7 "deadzspd" (holes in z in SPD)
1853   // Also given are the coordinates of the crossing point of track and module
1854   // (in the local module ref. system)
1855   // WARNING: THIS METHOD HAS TO BE SYNCHRONIZED WITH AliITStrackV2::GetModuleIndexInfo()!
1856   //----------------------------------------------------------------------
1857
1858   if(fITSModule[ilayer]==-1) {
1859     idet = -1;
1860     status=0;
1861     xloc=-99.; zloc=-99.;
1862     return kFALSE;
1863   }
1864
1865   Int_t module = fITSModule[ilayer];
1866
1867   idet = Int_t(module/1000000);
1868
1869   module -= idet*1000000;
1870
1871   status = Int_t(module/100000);
1872
1873   module -= status*100000;
1874
1875   Int_t signs = Int_t(module/10000);
1876
1877   module-=signs*10000;
1878
1879   Int_t xInt = Int_t(module/100);
1880   module -= xInt*100;
1881
1882   Int_t zInt = module;
1883
1884   if(signs==1) { xInt*=1; zInt*=1; }
1885   if(signs==2) { xInt*=1; zInt*=-1; }
1886   if(signs==3) { xInt*=-1; zInt*=1; }
1887   if(signs==4) { xInt*=-1; zInt*=-1; }
1888
1889   xloc = 0.1*(Float_t)xInt;
1890   zloc = 0.1*(Float_t)zInt;
1891
1892   if(status==4) idet = -1;
1893
1894   return kTRUE;
1895 }
1896
1897 //_______________________________________________________________________
1898 UShort_t AliESDtrack::GetTPCclusters(Int_t *idx) const {
1899   //---------------------------------------------------------------------
1900   // This function returns indices of the assgined ITS clusters
1901   //---------------------------------------------------------------------
1902   if (idx && fFriendTrack) {
1903     Int_t *index=fFriendTrack->GetTPCindices();
1904
1905     if (index){
1906       for (Int_t i=0; i<AliESDfriendTrack::kMaxTPCcluster; i++) idx[i]=index[i];
1907     }
1908     else {
1909       for (Int_t i=0; i<AliESDfriendTrack::kMaxTPCcluster; i++) idx[i]=-2;
1910     }
1911   }
1912   return fTPCncls;
1913 }
1914
1915 //_______________________________________________________________________
1916 Float_t AliESDtrack::GetTPCCrossedRows() const
1917 {
1918   // This function calls GetTPCClusterInfo with some default parameters which are used in the track selection and caches the outcome
1919   // because GetTPCClusterInfo is quite time-consuming
1920
1921   if (fCacheNCrossedRows > -1)
1922     return fCacheNCrossedRows;
1923
1924   fCacheNCrossedRows = GetTPCClusterInfo(2, 1);
1925   return fCacheNCrossedRows;
1926 }
1927
1928 //_______________________________________________________________________
1929 Float_t AliESDtrack::GetTPCClusterInfo(Int_t nNeighbours/*=3*/, Int_t type/*=0*/, Int_t row0, Int_t row1, Int_t bitType ) const
1930 {
1931   //
1932   // TPC cluster information
1933   // type 0: get fraction of found/findable clusters with neighbourhood definition
1934   //      1: findable clusters with neighbourhood definition
1935   //      2: found clusters
1936   // bitType:
1937   //      0 - all cluster used
1938   //      1 - clusters  used for the kalman update
1939   // definition of findable clusters:
1940   //            a cluster is defined as findable if there is another cluster
1941   //           within +- nNeighbours pad rows. The idea is to overcome threshold
1942   //           effects with a very simple algorithm.
1943   //
1944
1945
1946   Int_t found=0;
1947   Int_t findable=0;
1948   Int_t last=-nNeighbours;
1949   const TBits & clusterMap = (bitType%2==0) ? fTPCClusterMap : fTPCFitMap;
1950
1951   Int_t upperBound=clusterMap.GetNbits();
1952   if (upperBound>row1) upperBound=row1;
1953   for (Int_t i=row0; i<upperBound; ++i){
1954     //look to current row
1955     if (clusterMap[i]) {
1956       last=i;
1957       ++found;
1958       ++findable;
1959       continue;
1960     }
1961     //look to nNeighbours before
1962     if ((i-last)<=nNeighbours) {
1963       ++findable;
1964       continue;
1965     }
1966     //look to nNeighbours after
1967     for (Int_t j=i+1; j<i+1+nNeighbours; ++j){
1968       if (clusterMap[j]){
1969         ++findable;
1970         break;
1971       }
1972     }
1973   }
1974   if (type==2) return found;
1975   if (type==1) return findable;
1976
1977   if (type==0){
1978     Float_t fraction=0;
1979     if (findable>0)
1980       fraction=(Float_t)found/(Float_t)findable;
1981     else
1982       fraction=0;
1983     return fraction;
1984   }
1985   return 0;  // undefined type - default value
1986 }
1987
1988 //_______________________________________________________________________
1989 Float_t AliESDtrack::GetTPCClusterDensity(Int_t nNeighbours/*=3*/, Int_t type/*=0*/, Int_t row0, Int_t row1, Int_t bitType ) const
1990 {
1991   //
1992   // TPC cluster density -  only rows where signal before and after given row are used
1993   //                     -  slower function
1994   // type 0: get fraction of found/findable clusters with neighbourhood definition
1995   //      1: findable clusters with neighbourhood definition
1996   //      2: found clusters
1997   // bitType:
1998   //      0 - all cluster used
1999   //      1 - clusters  used for the kalman update
2000   // definition of findable clusters:
2001   //            a cluster is defined as findable if there is another cluster
2002   //           within +- nNeighbours pad rows. The idea is to overcome threshold
2003   //           effects with a very simple algorithm.
2004   //
2005   Int_t found=0;
2006   Int_t findable=0;
2007   //  Int_t last=-nNeighbours;
2008   const TBits & clusterMap = (bitType%2==0) ? fTPCClusterMap : fTPCFitMap;
2009   Int_t upperBound=clusterMap.GetNbits();
2010   if (upperBound>row1) upperBound=row1;
2011   for (Int_t i=row0; i<upperBound; ++i){
2012     Bool_t isUp=kFALSE;
2013     Bool_t isDown=kFALSE;
2014     for (Int_t idelta=1; idelta<=nNeighbours; idelta++){
2015       if (i-idelta>=0 && clusterMap[i-idelta]) isDown=kTRUE;
2016       if (i+idelta<upperBound && clusterMap[i+idelta]) isUp=kTRUE;
2017     }
2018     if (isUp&&isDown){
2019       ++findable;
2020       if (clusterMap[i]) ++found;
2021     }
2022   }
2023   if (type==2) return found;
2024   if (type==1) return findable;
2025
2026   if (type==0){
2027     Float_t fraction=0;
2028     if (findable>0)
2029       fraction=(Float_t)found/(Float_t)findable;
2030     else
2031       fraction=0;
2032     return fraction;
2033   }
2034   return 0;  // undefined type - default value
2035 }
2036
2037
2038
2039
2040 //_______________________________________________________________________
2041 Double_t AliESDtrack::GetTPCdensity(Int_t row0, Int_t row1) const{
2042   //
2043   // GetDensity of the clusters on given region between row0 and row1
2044   // Dead zone effect takin into acoount
2045   //
2046   if (!fFriendTrack) return 0.0;
2047   Int_t good  = 0;
2048   Int_t found = 0;
2049   //
2050   Int_t *index=fFriendTrack->GetTPCindices();
2051   for (Int_t i=row0;i<=row1;i++){
2052     Int_t idx = index[i];
2053     if (idx!=-1)  good++;             // track outside of dead zone
2054     if (idx>0)    found++;
2055   }
2056   Float_t density=0.5;
2057   if (good>TMath::Max((row1-row0)*0.5,0.0)) density = Float_t(found)/Float_t(good);
2058   return density;
2059 }
2060
2061 //_______________________________________________________________________
2062 void AliESDtrack::SetTPCpid(const Double_t *p) {
2063   // Sets values for the probability of each particle type (in TPC)
2064   SetPIDValues(fTPCr,p,AliPID::kSPECIES);
2065   SetStatus(AliESDtrack::kTPCpid);
2066 }
2067
2068 //_______________________________________________________________________
2069 void AliESDtrack::GetTPCpid(Double_t *p) const {
2070   // Gets the probability of each particle type (in TPC)
2071   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fTPCr[i];
2072 }
2073
2074 //_______________________________________________________________________
2075 UChar_t AliESDtrack::GetTRDclusters(Int_t *idx) const {
2076   //---------------------------------------------------------------------
2077   // This function returns indices of the assgined TRD clusters
2078   //---------------------------------------------------------------------
2079   if (idx && fFriendTrack) {
2080     Int_t *index=fFriendTrack->GetTRDindices();
2081
2082     if (index) {
2083       for (Int_t i=0; i<AliESDfriendTrack::kMaxTRDcluster; i++) idx[i]=index[i];
2084     }
2085     else {
2086       for (Int_t i=0; i<AliESDfriendTrack::kMaxTRDcluster; i++) idx[i]=-2;
2087     }
2088   }
2089   return fTRDncls;
2090 }
2091
2092 //_______________________________________________________________________
2093 UChar_t AliESDtrack::GetTRDtracklets(Int_t *idx) const {
2094 //
2095 // This function returns the number of TRD tracklets used in tracking
2096 // and it fills the indices of these tracklets in the array "idx" as they
2097 // are registered in the TRD track list.
2098 //
2099 // Caution :
2100 //   1. The idx array has to be allocated with a size >= AliESDtrack::kTRDnPlanes
2101 //   2. The idx array store not only the index but also the layer of the tracklet.
2102 //      Therefore tracks with TRD gaps contain default values for indices [-1]
2103
2104   if (!fFriendTrack) return 0;
2105   if (!idx) return GetTRDntracklets();
2106   Int_t *index=fFriendTrack->GetTRDindices();
2107   Int_t n = 0;
2108   for (Int_t i=0; i<kTRDnPlanes; i++){
2109     if (index){
2110       if(index[i]>=0) n++;
2111       idx[i]=index[i];
2112     }
2113     else idx[i] = -2;
2114   }
2115   return n;
2116 }
2117
2118 //_______________________________________________________________________
2119 void AliESDtrack::SetTRDpid(const Double_t *p) {
2120   // Sets values for the probability of each particle type (in TRD)
2121   SetPIDValues(fTRDr,p,AliPID::kSPECIES);
2122   SetStatus(AliESDtrack::kTRDpid);
2123 }
2124
2125 //_______________________________________________________________________
2126 void AliESDtrack::GetTRDpid(Double_t *p) const {
2127   // Gets the probability of each particle type (in TRD)
2128   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fTRDr[i];
2129 }
2130
2131 //_______________________________________________________________________
2132 void    AliESDtrack::SetTRDpid(Int_t iSpecies, Float_t p)
2133 {
2134   // Sets the probability of particle type iSpecies to p (in TRD)
2135   fTRDr[iSpecies] = p;
2136 }
2137
2138 Double_t AliESDtrack::GetTRDpid(Int_t iSpecies) const
2139 {
2140   // Returns the probability of particle type iSpecies (in TRD)
2141   return fTRDr[iSpecies];
2142 }
2143
2144 //____________________________________________________
2145 Int_t AliESDtrack::GetNumberOfTRDslices() const
2146 {
2147   // built in backward compatibility
2148   Int_t idx = fTRDnSlices - (kTRDnPlanes<<1);
2149   return idx<18 ? fTRDnSlices/kTRDnPlanes : idx/kTRDnPlanes;
2150 }
2151
2152 //____________________________________________________
2153 Double_t AliESDtrack::GetTRDmomentum(Int_t plane, Double_t *sp) const
2154 {
2155 //Returns momentum estimation and optional its error (sp)
2156 // in TRD layer "plane".
2157
2158   if (!fTRDnSlices) {
2159     AliDebug(2, "No TRD info allocated for this track.");
2160     return -1.;
2161   }
2162   if ((plane<0) || (plane>=kTRDnPlanes)) {
2163     AliWarning(Form("Request for TRD plane[%d] outside range.", plane));
2164     return -1.;
2165   }
2166
2167   Int_t idx = fTRDnSlices-(kTRDnPlanes<<1)+plane;
2168   // Protection for backward compatibility
2169   if(idx<(GetNumberOfTRDslices()*kTRDnPlanes)) return -1.;
2170
2171   if(sp) (*sp) = fTRDslices[idx+kTRDnPlanes];
2172   return fTRDslices[idx];
2173 }
2174
2175 //____________________________________________________
2176 Double_t  AliESDtrack::GetTRDslice(Int_t plane, Int_t slice) const {
2177   //Gets the charge from the slice of the plane
2178
2179   if(!fTRDslices) {
2180     //AliError("No TRD slices allocated for this track !");
2181     return -1.;
2182   }
2183   if ((plane<0) || (plane>=kTRDnPlanes)) {
2184     AliError("Info for TRD plane not available !");
2185     return -1.;
2186   }
2187   Int_t ns=GetNumberOfTRDslices();
2188   if ((slice<-1) || (slice>=ns)) {
2189     //AliError("Wrong TRD slice !");
2190     return -1.;
2191   }
2192
2193   if(slice>=0) return fTRDslices[plane*ns + slice];
2194
2195   // return average of the dEdx measurements
2196   Double_t q=0.; Double32_t *s = &fTRDslices[plane*ns];
2197   for (Int_t i=0; i<ns; i++, s++) if((*s)>0.) q+=(*s);
2198   return q/ns;
2199 }
2200
2201 //____________________________________________________
2202 void  AliESDtrack::SetNumberOfTRDslices(Int_t n) {
2203   //Sets the number of slices used for PID
2204   if (fTRDnSlices) return;
2205
2206   fTRDnSlices=n;
2207   fTRDslices=new Double32_t[fTRDnSlices];
2208
2209   // set-up correctly the allocated memory
2210   memset(fTRDslices, 0, n*sizeof(Double32_t));
2211   for (Int_t i=GetNumberOfTRDslices(); i--;) fTRDslices[i]=-1.;
2212 }
2213
2214 //____________________________________________________
2215 void  AliESDtrack::SetTRDslice(Double_t q, Int_t plane, Int_t slice) {
2216   //Sets the charge q in the slice of the plane
2217   if(!fTRDslices) {
2218     AliError("No TRD slices allocated for this track !");
2219     return;
2220   }
2221   if ((plane<0) || (plane>=kTRDnPlanes)) {
2222     AliError("Info for TRD plane not allocated !");
2223     return;
2224   }
2225   Int_t ns=GetNumberOfTRDslices();
2226   if ((slice<0) || (slice>=ns)) {
2227     AliError(Form("Wrong TRD slice %d/%d, NSlices=%d",plane,slice,ns));
2228     return;
2229   }
2230   Int_t n=plane*ns + slice;
2231   fTRDslices[n]=q;
2232 }
2233
2234
2235 //____________________________________________________
2236 void AliESDtrack::SetTRDmomentum(Double_t p, Int_t plane, Double_t *sp)
2237 {
2238   if(!fTRDslices) {
2239     AliError("No TRD slices allocated for this track !");
2240     return;
2241   }
2242   if ((plane<0) || (plane>=kTRDnPlanes)) {
2243     AliError("Info for TRD plane not allocated !");
2244     return;
2245   }
2246
2247   Int_t idx = fTRDnSlices-(kTRDnPlanes<<1)+plane;
2248   // Protection for backward compatibility
2249   if(idx<GetNumberOfTRDslices()*kTRDnPlanes) return;
2250
2251   if(sp) fTRDslices[idx+kTRDnPlanes] = (*sp);
2252   fTRDslices[idx] = p;
2253 }
2254
2255
2256 //_______________________________________________________________________
2257 void AliESDtrack::SetTOFpid(const Double_t *p) {
2258   // Sets the probability of each particle type (in TOF)
2259   SetPIDValues(fTOFr,p,AliPID::kSPECIES);
2260   SetStatus(AliESDtrack::kTOFpid);
2261 }
2262
2263 //_______________________________________________________________________
2264 void AliESDtrack::SetTOFLabel(const Int_t *p) {
2265   // Sets  (in TOF)
2266   if(!fTOFLabel) fTOFLabel = new Int_t[3];
2267   for (Int_t i=0; i<3; i++) fTOFLabel[i]=p[i];
2268 }
2269
2270 //_______________________________________________________________________
2271 void AliESDtrack::GetTOFpid(Double_t *p) const {
2272   // Gets probabilities of each particle type (in TOF)
2273   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fTOFr[i];
2274 }
2275
2276 //_______________________________________________________________________
2277 void AliESDtrack::GetTOFLabel(Int_t *p) const {
2278   // Gets (in TOF)
2279   if(fNtofClusters>0){
2280     TObjArray *tofclArray = GetESDEvent()->GetTOFcluster();
2281     AliESDTOFcluster *tofcl = (AliESDTOFcluster *) tofclArray->At(fTOFcluster[0]);
2282
2283     for (Int_t i=0; i<3; i++) p[i]=tofcl->GetLabel(i);
2284   }
2285   else{
2286     if(fTOFLabel)
2287       for (Int_t i=0; i<3; i++) p[i]=fTOFLabel[i];
2288   }
2289 }
2290
2291 //_______________________________________________________________________
2292 void AliESDtrack::GetTOFInfo(Float_t *info) const {
2293   // Gets (in TOF)
2294   for (Int_t i=0; i<10; i++) info[i]=fTOFInfo[i];
2295 }
2296
2297 //_______________________________________________________________________
2298 void AliESDtrack::SetTOFInfo(Float_t*info) {
2299   // Gets (in TOF)
2300   for (Int_t i=0; i<10; i++) fTOFInfo[i]=info[i];
2301 }
2302
2303
2304
2305 //_______________________________________________________________________
2306 void AliESDtrack::SetHMPIDpid(const Double_t *p) {
2307   // Sets the probability of each particle type (in HMPID)
2308   SetPIDValues(fHMPIDr,p,AliPID::kSPECIES);
2309   SetStatus(AliESDtrack::kHMPIDpid);
2310 }
2311
2312 //_______________________________________________________________________
2313 void  AliESDtrack::SetTPCdEdxInfo(AliTPCdEdxInfo * dEdxInfo){
2314   if(fTPCdEdxInfo) delete fTPCdEdxInfo;
2315   fTPCdEdxInfo = dEdxInfo;
2316 }
2317
2318 //_______________________________________________________________________
2319 void AliESDtrack::GetHMPIDpid(Double_t *p) const {
2320   // Gets probabilities of each particle type (in HMPID)
2321   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fHMPIDr[i];
2322 }
2323
2324
2325
2326 //_______________________________________________________________________
2327 void AliESDtrack::SetESDpid(const Double_t *p) {
2328   // Sets the probability of each particle type for the ESD track
2329   SetPIDValues(fR,p,AliPID::kSPECIES);
2330   SetStatus(AliESDtrack::kESDpid);
2331 }
2332
2333 //_______________________________________________________________________
2334 void AliESDtrack::GetESDpid(Double_t *p) const {
2335   // Gets probability of each particle type for the ESD track
2336   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fR[i];
2337 }
2338
2339 //_______________________________________________________________________
2340 Bool_t AliESDtrack::RelateToVertexTPC(const AliESDVertex *vtx,
2341 Double_t b, Double_t maxd, AliExternalTrackParam *cParam) {
2342   //
2343   // Try to relate the TPC-only track parameters to the vertex "vtx",
2344   // if the (rough) transverse impact parameter is not bigger then "maxd".
2345   //            Magnetic field is "b" (kG).
2346   //
2347   // a) The TPC-only paramters are extapolated to the DCA to the vertex.
2348   // b) The impact parameters and their covariance matrix are calculated.
2349   // c) An attempt to constrain the TPC-only params to the vertex is done.
2350   //    The constrained params are returned via "cParam".
2351   //
2352   // In the case of success, the returned value is kTRUE
2353   // otherwise, it's kFALSE)
2354   //
2355
2356   if (!fTPCInner) return kFALSE;
2357   if (!vtx) return kFALSE;
2358
2359   Double_t dz[2],cov[3];
2360   if (!fTPCInner->PropagateToDCA(vtx, b, maxd, dz, cov)) return kFALSE;
2361
2362   fdTPC = dz[0];
2363   fzTPC = dz[1];
2364   fCddTPC = cov[0];
2365   fCdzTPC = cov[1];
2366   fCzzTPC = cov[2];
2367
2368   Double_t covar[6]; vtx->GetCovMatrix(covar);
2369   Double_t p[2]={GetParameter()[0]-dz[0],GetParameter()[1]-dz[1]};
2370   Double_t c[3]={covar[2],0.,covar[5]};
2371
2372   Double_t chi2=GetPredictedChi2(p,c);
2373   if (chi2>kVeryBig) return kFALSE;
2374
2375   fCchi2TPC=chi2;
2376
2377   if (!cParam) return kTRUE;
2378
2379   *cParam = *fTPCInner;
2380   if (!cParam->Update(p,c)) return kFALSE;
2381
2382   return kTRUE;
2383 }
2384
2385 //_______________________________________________________________________
2386 Bool_t AliESDtrack::RelateToVertexTPCBxByBz(const AliESDVertex *vtx,
2387 Double_t b[3], Double_t maxd, AliExternalTrackParam *cParam) {
2388   //
2389   // Try to relate the TPC-only track parameters to the vertex "vtx",
2390   // if the (rough) transverse impact parameter is not bigger then "maxd".
2391   //
2392   // All three components of the magnetic field ,"b[3]" (kG),
2393   // are taken into account.
2394   //
2395   // a) The TPC-only paramters are extapolated to the DCA to the vertex.
2396   // b) The impact parameters and their covariance matrix are calculated.
2397   // c) An attempt to constrain the TPC-only params to the vertex is done.
2398   //    The constrained params are returned via "cParam".
2399   //
2400   // In the case of success, the returned value is kTRUE
2401   // otherwise, it's kFALSE)
2402   //
2403
2404   if (!fTPCInner) return kFALSE;
2405   if (!vtx) return kFALSE;
2406
2407   Double_t dz[2],cov[3];
2408   if (!fTPCInner->PropagateToDCABxByBz(vtx, b, maxd, dz, cov)) return kFALSE;
2409
2410   fdTPC = dz[0];
2411   fzTPC = dz[1];
2412   fCddTPC = cov[0];
2413   fCdzTPC = cov[1];
2414   fCzzTPC = cov[2];
2415
2416   Double_t covar[6]; vtx->GetCovMatrix(covar);
2417   Double_t p[2]={GetParameter()[0]-dz[0],GetParameter()[1]-dz[1]};
2418   Double_t c[3]={covar[2],0.,covar[5]};
2419
2420   Double_t chi2=GetPredictedChi2(p,c);
2421   if (chi2>kVeryBig) return kFALSE;
2422
2423   fCchi2TPC=chi2;
2424
2425   if (!cParam) return kTRUE;
2426
2427   *cParam = *fTPCInner;
2428   if (!cParam->Update(p,c)) return kFALSE;
2429
2430   return kTRUE;
2431 }
2432
2433 //_______________________________________________________________________
2434 Bool_t AliESDtrack::RelateToVertex(const AliESDVertex *vtx,
2435 Double_t b, Double_t maxd, AliExternalTrackParam *cParam) {
2436   //
2437   // Try to relate this track to the vertex "vtx",
2438   // if the (rough) transverse impact parameter is not bigger then "maxd".
2439   //            Magnetic field is "b" (kG).
2440   //
2441   // a) The track gets extapolated to the DCA to the vertex.
2442   // b) The impact parameters and their covariance matrix are calculated.
2443   // c) An attempt to constrain this track to the vertex is done.
2444   //    The constrained params are returned via "cParam".
2445   //
2446   // In the case of success, the returned value is kTRUE
2447   // (otherwise, it's kFALSE)
2448   //
2449
2450   if (!vtx) return kFALSE;
2451
2452   Double_t dz[2],cov[3];
2453   if (!PropagateToDCA(vtx, b, maxd, dz, cov)) return kFALSE;
2454
2455   fD = dz[0];
2456   fZ = dz[1];
2457   fCdd = cov[0];
2458   fCdz = cov[1];
2459   fCzz = cov[2];
2460
2461   Double_t covar[6]; vtx->GetCovMatrix(covar);
2462   Double_t p[2]={GetParameter()[0]-dz[0],GetParameter()[1]-dz[1]};
2463   Double_t c[3]={covar[2],0.,covar[5]};
2464
2465   Double_t chi2=GetPredictedChi2(p,c);
2466   if (chi2>kVeryBig) return kFALSE;
2467
2468   fCchi2=chi2;
2469
2470
2471   //--- Could now these lines be removed ? ---
2472   delete fCp;
2473   fCp=new AliExternalTrackParam(*this);
2474
2475   if (!fCp->Update(p,c)) {delete fCp; fCp=0; return kFALSE;}
2476   //----------------------------------------
2477
2478   fVertexID = vtx->GetID();
2479
2480   if (!cParam) return kTRUE;
2481
2482   *cParam = *this;
2483   if (!cParam->Update(p,c)) return kFALSE;
2484
2485   return kTRUE;
2486 }
2487
2488 //_______________________________________________________________________
2489 Bool_t AliESDtrack::RelateToVertexBxByBz(const AliESDVertex *vtx,
2490 Double_t b[3], Double_t maxd, AliExternalTrackParam *cParam) {
2491   //
2492   // Try to relate this track to the vertex "vtx",
2493   // if the (rough) transverse impact parameter is not bigger then "maxd".
2494   //            Magnetic field is "b" (kG).
2495   //
2496   // a) The track gets extapolated to the DCA to the vertex.
2497   // b) The impact parameters and their covariance matrix are calculated.
2498   // c) An attempt to constrain this track to the vertex is done.
2499   //    The constrained params are returned via "cParam".
2500   //
2501   // In the case of success, the returned value is kTRUE
2502   // (otherwise, it's kFALSE)
2503   //
2504
2505   if (!vtx) return kFALSE;
2506
2507   Double_t dz[2],cov[3];
2508   if (!PropagateToDCABxByBz(vtx, b, maxd, dz, cov)) return kFALSE;
2509
2510   fD = dz[0];
2511   fZ = dz[1];
2512   fCdd = cov[0];
2513   fCdz = cov[1];
2514   fCzz = cov[2];
2515
2516   Double_t covar[6]; vtx->GetCovMatrix(covar);
2517   Double_t p[2]={GetParameter()[0]-dz[0],GetParameter()[1]-dz[1]};
2518   Double_t c[3]={covar[2],0.,covar[5]};
2519
2520   Double_t chi2=GetPredictedChi2(p,c);
2521   if (chi2>kVeryBig) return kFALSE;
2522
2523   fCchi2=chi2;
2524
2525
2526   //--- Could now these lines be removed ? ---
2527   delete fCp;
2528   fCp=new AliExternalTrackParam(*this);
2529
2530   if (!fCp->Update(p,c)) {delete fCp; fCp=0; return kFALSE;}
2531   //----------------------------------------
2532
2533   fVertexID = vtx->GetID();
2534
2535   if (!cParam) return kTRUE;
2536
2537   *cParam = *this;
2538   if (!cParam->Update(p,c)) return kFALSE;
2539
2540   return kTRUE;
2541 }
2542
2543 //_______________________________________________________________________
2544 void AliESDtrack::Print(Option_t *) const {
2545   // Prints info on the track
2546   AliExternalTrackParam::Print();
2547   printf("ESD track info\n") ;
2548   Double_t p[AliPID::kSPECIES] ;
2549   Int_t index = 0 ;
2550   if( IsOn(kITSpid) ){
2551     printf("From ITS: ") ;
2552     GetITSpid(p) ;
2553     for(index = 0 ; index < AliPID::kSPECIES; index++)
2554       printf("%f, ", p[index]) ;
2555     printf("\n           signal = %f\n", GetITSsignal()) ;
2556   }
2557   if( IsOn(kTPCpid) ){
2558     printf("From TPC: ") ;
2559     GetTPCpid(p) ;
2560     for(index = 0 ; index < AliPID::kSPECIES; index++)
2561       printf("%f, ", p[index]) ;
2562     printf("\n           signal = %f\n", GetTPCsignal()) ;
2563   }
2564   if( IsOn(kTRDpid) ){
2565     printf("From TRD: ") ;
2566     GetTRDpid(p) ;
2567     for(index = 0 ; index < AliPID::kSPECIES; index++)
2568       printf("%f, ", p[index]) ;
2569       printf("\n           signal = %f\n", GetTRDsignal()) ;
2570       printf("\n           NchamberdEdx = %d\n", GetTRDNchamberdEdx()) ;
2571       printf("\n           NclusterdEdx = %d\n", GetTRDNclusterdEdx()) ;
2572   }
2573   if( IsOn(kTOFpid) ){
2574     printf("From TOF: ") ;
2575     GetTOFpid(p) ;
2576     for(index = 0 ; index < AliPID::kSPECIES; index++)
2577       printf("%f, ", p[index]) ;
2578     printf("\n           signal = %f\n", GetTOFsignal()) ;
2579   }
2580   if( IsOn(kHMPIDpid) ){
2581     printf("From HMPID: ") ;
2582     GetHMPIDpid(p) ;
2583     for(index = 0 ; index < AliPID::kSPECIES; index++)
2584       printf("%f, ", p[index]) ;
2585     printf("\n           signal = %f\n", GetHMPIDsignal()) ;
2586   }
2587 }
2588
2589
2590 //
2591 // Draw functionality
2592 // Origin: Marian Ivanov, Marian.Ivanov@cern.ch
2593 //
2594 void AliESDtrack::FillPolymarker(TPolyMarker3D *pol, Float_t magF, Float_t minR, Float_t maxR, Float_t stepR){
2595   //
2596   // Fill points in the polymarker
2597   //
2598   TObjArray arrayRef;
2599   arrayRef.AddLast(new AliExternalTrackParam(*this));
2600   if (fIp) arrayRef.AddLast(new AliExternalTrackParam(*fIp));
2601   if (fOp) arrayRef.AddLast(new AliExternalTrackParam(*fOp));
2602   if (fHMPIDp) arrayRef.AddLast(new AliExternalTrackParam(*fHMPIDp));
2603   //
2604   Double_t mpos[3]={0,0,0};
2605   Int_t entries=arrayRef.GetEntries();
2606   for (Int_t i=0;i<entries;i++){
2607     Double_t pos[3];
2608     ((AliExternalTrackParam*)arrayRef.At(i))->GetXYZ(pos);
2609     mpos[0]+=pos[0]/entries;
2610     mpos[1]+=pos[1]/entries;
2611     mpos[2]+=pos[2]/entries;
2612   }
2613   // Rotate to the mean position
2614   //
2615   Float_t fi= TMath::ATan2(mpos[1],mpos[0]);
2616   for (Int_t i=0;i<entries;i++){
2617     Bool_t res = ((AliExternalTrackParam*)arrayRef.At(i))->Rotate(fi);
2618     if (!res) delete arrayRef.RemoveAt(i);
2619   }
2620   Int_t counter=0;
2621   for (Double_t r=minR; r<maxR; r+=stepR){
2622     Double_t sweight=0;
2623     Double_t mlpos[3]={0,0,0};
2624     for (Int_t i=0;i<entries;i++){
2625       Double_t point[3]={0,0,0};
2626       AliExternalTrackParam *param = ((AliExternalTrackParam*)arrayRef.At(i));
2627       if (!param) continue;
2628       if (param->GetXYZAt(r,magF,point)){
2629         Double_t weight = 1./(10.+(r-param->GetX())*(r-param->GetX()));
2630         sweight+=weight;
2631         mlpos[0]+=point[0]*weight;
2632         mlpos[1]+=point[1]*weight;
2633         mlpos[2]+=point[2]*weight;
2634       }
2635     }
2636     if (sweight>0){
2637       mlpos[0]/=sweight;
2638       mlpos[1]/=sweight;
2639       mlpos[2]/=sweight;
2640       pol->SetPoint(counter,mlpos[0],mlpos[1], mlpos[2]);
2641       //      printf("xyz\t%f\t%f\t%f\n",mlpos[0], mlpos[1],mlpos[2]);
2642       counter++;
2643     }
2644   }
2645 }
2646
2647 //_______________________________________________________________________
2648 void AliESDtrack::SetITSdEdxSamples(const Double_t s[4]) {
2649   //
2650   // Store the dE/dx samples measured by the two SSD and two SDD layers.
2651   // These samples are corrected for the track segment length.
2652   //
2653   for (Int_t i=0; i<4; i++) fITSdEdxSamples[i]=s[i];
2654 }
2655
2656 //_______________________________________________________________________
2657 void AliESDtrack::GetITSdEdxSamples(Double_t s[4]) const {
2658   //
2659   // Get the dE/dx samples measured by the two SSD and two SDD layers.
2660   // These samples are corrected for the track segment length.
2661   //
2662   for (Int_t i=0; i<4; i++) s[i]=fITSdEdxSamples[i];
2663 }
2664
2665
2666 UShort_t   AliESDtrack::GetTPCnclsS(Int_t i0,Int_t i1) const{
2667   //
2668   // get number of shared TPC clusters
2669   //
2670   return  fTPCSharedMap.CountBits(i0)-fTPCSharedMap.CountBits(i1);
2671 }
2672
2673 UShort_t   AliESDtrack::GetTPCncls(Int_t i0,Int_t i1) const{
2674   //
2675   // get number of TPC clusters
2676   //
2677   return  fTPCClusterMap.CountBits(i0)-fTPCClusterMap.CountBits(i1);
2678 }
2679
2680 //____________________________________________________________________
2681 Double_t AliESDtrack::GetChi2TPCConstrainedVsGlobal(const AliESDVertex* vtx) const
2682 {
2683   // Calculates the chi2 between the TPC track (TPCinner) constrained to the primary vertex and the global track
2684   //
2685   // Returns -1 in case the calculation failed
2686   //
2687   // Value is cached as a non-persistent member.
2688   //
2689   // Code adapted from original code by GSI group (Jacek, Marian, Michael)
2690
2691   // cache, ignoring that a different vertex might be passed
2692   if (fCacheChi2TPCConstrainedVsGlobalVertex == vtx)
2693     return fCacheChi2TPCConstrainedVsGlobal;
2694
2695   fCacheChi2TPCConstrainedVsGlobal = -1;
2696   fCacheChi2TPCConstrainedVsGlobalVertex = vtx;
2697
2698   Double_t x[3];
2699   GetXYZ(x);
2700   Double_t b[3];
2701   AliTrackerBase::GetBxByBz(x,b);
2702
2703   if (!fTPCInner)  {
2704     AliWarning("Could not get TPC Inner Param.");
2705     return fCacheChi2TPCConstrainedVsGlobal;
2706   }
2707
2708   // clone for constraining
2709   AliExternalTrackParam* tpcInnerC = new AliExternalTrackParam(*fTPCInner);
2710   if (!tpcInnerC) {
2711     AliWarning("Clone of TPCInnerParam failed.");
2712     return fCacheChi2TPCConstrainedVsGlobal;
2713   }
2714
2715   // transform to the track reference frame
2716   Bool_t isOK = tpcInnerC->Rotate(GetAlpha());
2717   isOK &= tpcInnerC->PropagateTo(GetX(), b[2]);
2718   if (!isOK) {
2719     delete tpcInnerC;
2720     tpcInnerC = 0;
2721     AliWarning("Rotation/Propagation of track failed.") ;
2722     return fCacheChi2TPCConstrainedVsGlobal;
2723   }
2724
2725   // constrain TPCinner
2726   isOK = tpcInnerC->ConstrainToVertex(vtx, b);
2727
2728   // transform to the track reference frame
2729   isOK &= tpcInnerC->Rotate(GetAlpha());
2730   isOK &= tpcInnerC->PropagateTo(GetX(), b[2]);
2731
2732   if (!isOK) {
2733     AliWarning("ConstrainTPCInner failed.") ;
2734     delete tpcInnerC;
2735     tpcInnerC = 0;
2736     return fCacheChi2TPCConstrainedVsGlobal;
2737   }
2738
2739   // calculate chi2 between vi and vj vectors
2740   // with covi and covj covariance matrices
2741   // chi2ij = (vi-vj)^(T)*(covi+covj)^(-1)*(vi-vj)
2742   TMatrixD deltaT(5,1);
2743   TMatrixD delta(1,5);
2744   TMatrixD covarM(5,5);
2745
2746   for (Int_t ipar=0; ipar<5; ipar++) {
2747     deltaT(ipar,0) = tpcInnerC->GetParameter()[ipar] - GetParameter()[ipar];
2748     delta(0,ipar) = tpcInnerC->GetParameter()[ipar] - GetParameter()[ipar];
2749
2750     for (Int_t jpar=0; jpar<5; jpar++) {
2751       Int_t index = GetIndex(ipar,jpar);
2752       covarM(ipar,jpar) = GetCovariance()[index]+tpcInnerC->GetCovariance()[index];
2753     }
2754   }
2755   // chi2 distance TPC constrained and TPC+ITS
2756   TMatrixD covarMInv = covarM.Invert();
2757   TMatrixD mat2 = covarMInv*deltaT;
2758   TMatrixD chi2 = delta*mat2;
2759
2760   delete tpcInnerC;
2761   tpcInnerC = 0;
2762
2763   fCacheChi2TPCConstrainedVsGlobal = chi2(0,0);
2764   return fCacheChi2TPCConstrainedVsGlobal;
2765 }
2766
2767 void AliESDtrack::SetDetectorPID(const AliDetectorPID *pid)
2768 {
2769   //
2770   // Set the detector PID
2771   //
2772   if (fDetectorPID) delete fDetectorPID;
2773   fDetectorPID=pid;
2774
2775 }
2776
2777 Double_t AliESDtrack::GetLengthInActiveZone( Int_t mode, Double_t deltaY, Double_t deltaZ, Double_t bz, Double_t exbPhi , TTreeSRedirector * pcstream) const {
2778   //
2779   // Input parameters:
2780   //   mode  - type of external track parameters
2781   //   deltaY - user defined "dead region" in cm
2782   //   deltaZ - user defined "active region" in cm (250 cm drift lenght - 14 cm L1 delay
2783   //   bz     - magnetic field
2784   //   exbPhi - optional rotation due to the ExB effect
2785   // return value:
2786   //   the length of the track in cm in "active volume" of the TPC
2787   //
2788   if (mode==0) return GetLengthInActiveZone(this, deltaY,deltaZ,bz, exbPhi,pcstream);
2789   if (mode==1) return GetLengthInActiveZone(fIp, deltaY,deltaZ,bz, exbPhi,pcstream);
2790   if (mode==2) return GetLengthInActiveZone(fOp, deltaY,deltaZ,bz, exbPhi,pcstream);
2791   return 0;
2792 }
2793
2794 Double_t AliESDtrack::GetLengthInActiveZone(const AliExternalTrackParam  *paramT, Double_t deltaY, Double_t deltaZ, Double_t bz, Double_t exbPhi , TTreeSRedirector * pcstream) const {
2795   //
2796   // Numerical code to calculate the length of the track in active region of the TPC
2797   // ( can be speed up if somebody wants to invest time - analysical version shoult be possible)
2798   //
2799   // Input parameters:
2800   //   paramT - external track parameters
2801   //   deltaY - user defined "dead region" in cm
2802   //   deltaZ - user defined "active region" in cm (250 cm drift lenght - 14 cm L1 delay
2803   //   bz     - magnetic field
2804   //   exbPhi - optional rotation due to the ExB effect
2805   // return value:
2806   //   the length of the track in cm in "active volume" of the TPC
2807   //
2808   const Double_t rIn=85;
2809   const Double_t rOut=245;
2810   Double_t xyz[3], pxyz[3];
2811   if (paramT->GetXYZAt(rIn,bz,xyz)){
2812     paramT->GetPxPyPzAt(rIn,bz,pxyz);
2813   }else{
2814     paramT->GetXYZ(xyz);
2815     paramT->GetPxPyPz(pxyz);
2816   }
2817   //
2818   Double_t dca   = -paramT->GetD(0,0,bz);  // get impact parameter distance to point (0,0)
2819   Double_t radius= TMath::Abs(1/paramT->GetC(bz));  //
2820   Double_t sign  = paramT->GetSign();
2821   Double_t R0    = TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1]);   // radius at current point
2822   Double_t phiR0 = TMath::ATan2(xyz[1],xyz[0]);                // angle of given point
2823   Double_t dPhiR0= -TMath::ASin((dca*dca-2*dca*radius*sign+R0*R0)/(2*R0*(dca-radius*sign)));
2824   Double_t phi0  = phiR0-(dPhiR0);  // global phi offset to be added
2825   //
2826   //
2827   AliExternalTrackParam paramR=(*paramT);
2828   Double_t length=0;
2829   for (Double_t R=rIn; R<=rOut; R++){
2830     Double_t sinPhi=(dca*dca-2*dca*radius*sign+R*R)/(2*R*(dca-radius*sign));
2831     if (TMath::Abs(sinPhi)>=1) continue;
2832     Double_t dphi     = -TMath::ASin(sinPhi);
2833     Double_t phi      = phi0+dphi;                           // global phi
2834     Int_t    sector   = TMath::Nint(9*phi/(TMath::Pi()));
2835     Double_t dPhiEdge = phi-(sector*TMath::Pi()/9)+exbPhi;   // distance to sector boundary in rphi
2836     Double_t dX   = R*TMath::Cos(phi)-xyz[0];
2837     Double_t dY   = R*TMath::Sin(phi)-xyz[1];
2838     Double_t deltaPhi = 2*TMath::ASin(0.5*TMath::Sqrt(dX*dX+dY*dY)/radius);
2839     Double_t z = xyz[2]+deltaPhi*radius*paramT->GetTgl();
2840     if (TMath::Abs(dPhiEdge*R)>deltaY && TMath::Abs(z)<deltaZ){
2841       length++;
2842     }
2843     //    Double_t deltaZ= dphi*radius;
2844     if (pcstream){
2845       //should we keep debug possibility ?
2846       AliExternalTrackParam paramTcopy=(*paramT);
2847       paramR.Rotate(phi);
2848       paramR.PropagateTo(R,bz);
2849       (*pcstream)<<"debugEdge"<<
2850         "R="<<R<<                   // radius
2851         "dphiEdge="<<dPhiEdge<<     // distance to edge
2852         "phi0="<<phi0<<             // phi0 -phi at the track initial position
2853         "phi="<<phi<<               //
2854         "z="<<z<<
2855         "pT.="<<&paramTcopy<<
2856         "pR.="<<&paramR<<
2857         "\n";
2858     }
2859   }
2860   return length;
2861 }
2862
2863 void    AliESDtrack::SetTOFclusterArray(Int_t ncluster,Int_t *TOFcluster){
2864   AliInfo("Method has to be implemented!");
2865 //   fNtofClusters=ncluster;
2866 //   if(TOFcluster == fTOFcluster) return;
2867 //   if(fTOFcluster){ // reset previous content
2868 //     delete[] fTOFcluster;
2869 //     fTOFcluster = NULL;
2870 //     fNtofClusters=0;
2871 //   }
2872
2873 //   if(ncluster){ // set new content
2874 //     fTOFcluster = new Int_t[fNtofClusters];
2875 //     for(Int_t i=0;i < fNtofClusters;i++) fTOFcluster[i] = TOFcluster[i];
2876 //   }
2877 //   else
2878 //     fTOFcluster = 0;
2879 }
2880
2881 void    AliESDtrack::AddTOFcluster(Int_t icl){
2882   fNtofClusters++;
2883
2884   Int_t *old = fTOFcluster;
2885   fTOFcluster = new Int_t[fNtofClusters];
2886
2887   for(Int_t i=0;i < fNtofClusters-1;i++) fTOFcluster[i] = old[i];
2888   fTOFcluster[fNtofClusters-1] = icl;
2889
2890   if(fNtofClusters-1){ // delete previous content
2891     delete old;
2892     old = NULL;
2893   }
2894
2895 }
2896
2897 Double_t AliESDtrack::GetTOFsignal() const {
2898   if(fNtofClusters>0 && GetESDEvent()){
2899     TObjArray *tofclArray = GetESDEvent()->GetTOFcluster();
2900     AliESDTOFcluster *tofcl = (AliESDTOFcluster *) tofclArray->At(fTOFcluster[0]);
2901
2902     return tofcl->GetTime();
2903   }
2904   else if(fNtofClusters>0) AliInfo("No AliESDEvent available here!\n");
2905
2906   return fTOFsignal;
2907 }
2908
2909 Double_t AliESDtrack::GetTOFsignalToT() const
2910 {
2911   if(fNtofClusters>0 && GetESDEvent()){
2912     TObjArray *tofclArray = GetESDEvent()->GetTOFcluster();
2913     AliESDTOFcluster *tofcl = (AliESDTOFcluster *) tofclArray->At(fTOFcluster[0]);
2914
2915     return tofcl->GetTOT();
2916   }
2917   else if(fNtofClusters>0) AliInfo("No AliESDEvent available here!\n");
2918
2919   return fTOFsignalToT;
2920 }
2921
2922 Double_t AliESDtrack::GetTOFsignalRaw() const
2923 {
2924   if(fNtofClusters>0 && GetESDEvent()){
2925     TObjArray *tofclArray = GetESDEvent()->GetTOFcluster();
2926     AliESDTOFcluster *tofcl = (AliESDTOFcluster *) tofclArray->At(fTOFcluster[0]);
2927
2928     return tofcl->GetTimeRaw();
2929   }
2930   else if(fNtofClusters>0) AliInfo("No AliESDEvent available here!\n");
2931
2932   return fTOFsignalRaw;
2933 }
2934
2935 Double_t AliESDtrack::GetTOFsignalDz() const
2936 {
2937
2938   AliESDTOFcluster *tofcl;
2939
2940   Int_t index = -1;
2941   if(fNtofClusters>0 && GetESDEvent()){
2942     TObjArray *tofclArray = GetESDEvent()->GetTOFcluster();
2943     tofcl = (AliESDTOFcluster *) tofclArray->At(fTOFcluster[0]);
2944
2945     for(Int_t i=0;i < tofcl->GetNMatchableTracks();i++){
2946       if(tofcl->GetTrackIndex(i) == GetID()) index = i;
2947     }
2948   }
2949   else if(fNtofClusters>0) AliInfo("No AliESDEvent available here!\n");
2950
2951   if(fNtofClusters>0 && index > -1){
2952     return tofcl->GetDz(index);
2953   }
2954   return fTOFsignalDz;
2955 }
2956
2957 Double_t AliESDtrack::GetTOFsignalDx() const
2958 {
2959   AliESDTOFcluster *tofcl;
2960
2961   Int_t index = -1;
2962   if(fNtofClusters>0 && GetESDEvent()){
2963     TObjArray *tofclArray = GetESDEvent()->GetTOFcluster();
2964     tofcl = (AliESDTOFcluster *) tofclArray->At(fTOFcluster[0]);
2965     for(Int_t i=0;i < tofcl->GetNMatchableTracks();i++){
2966       if(tofcl->GetTrackIndex(i) == GetID()) index = i;
2967     }
2968   }
2969   else if(fNtofClusters>0) AliInfo("No AliESDEvent available here!\n");
2970   if(fNtofClusters>0 && index > -1){
2971     return tofcl->GetDx(index);
2972   }
2973   return fTOFsignalDx;
2974 }
2975
2976 Short_t  AliESDtrack::GetTOFDeltaBC() const
2977 {
2978   if(fNtofClusters>0 && GetESDEvent()){
2979     TObjArray *tofclArray = GetESDEvent()->GetTOFcluster();
2980     AliESDTOFcluster *tofcl = (AliESDTOFcluster *) tofclArray->At(fTOFcluster[0]);
2981     return tofcl->GetDeltaBC();
2982   }
2983   else if(fNtofClusters>0) AliInfo("No AliESDEvent available here!\n");
2984
2985   return fTOFdeltaBC;
2986 }
2987
2988 Short_t  AliESDtrack::GetTOFL0L1() const
2989 {
2990   if(fNtofClusters>0 && GetESDEvent()){
2991     TObjArray *tofclArray = GetESDEvent()->GetTOFcluster();
2992     AliESDTOFcluster *tofcl = (AliESDTOFcluster *) tofclArray->At(fTOFcluster[0]);
2993
2994     return tofcl->GetL0L1Latency();
2995   }
2996   else if(fNtofClusters>0) AliInfo("No AliESDEvent available here!\n");
2997
2998   return fTOFl0l1;
2999 }
3000
3001 Int_t   AliESDtrack::GetTOFCalChannel() const
3002 {
3003   if(fNtofClusters>0 && GetESDEvent()){
3004     TObjArray *tofclArray = GetESDEvent()->GetTOFcluster();
3005     AliESDTOFcluster *tofcl = (AliESDTOFcluster *) tofclArray->At(fTOFcluster[0]);
3006
3007     tofcl->GetTOFchannel();
3008   }
3009   else if(fNtofClusters>0) AliInfo("No AliESDEvent available here!\n");
3010
3011   return fTOFCalChannel;
3012 }
3013
3014 Int_t   AliESDtrack::GetTOFcluster() const
3015 {
3016   if(fNtofClusters>0 && GetESDEvent()){
3017     TObjArray *tofclArray = GetESDEvent()->GetTOFcluster();
3018     AliESDTOFcluster *tofcl = (AliESDTOFcluster *) tofclArray->At(fTOFcluster[0]);
3019
3020     tofcl->GetClusterIndex();
3021   }
3022   else if(fNtofClusters>0) AliInfo("No AliESDEvent available here!\n");
3023
3024   return fTOFindex;
3025 }
3026
3027 Int_t   AliESDtrack::GetTOFclusterN() const
3028 {
3029   return fNtofClusters;
3030 }
3031
3032 Bool_t  AliESDtrack::IsTOFHitAlreadyMatched() const{
3033   if(fNtofClusters>0 && GetESDEvent()){
3034     TObjArray *tofclArray = GetESDEvent()->GetTOFcluster();
3035     AliESDTOFcluster *tofcl = (AliESDTOFcluster *) tofclArray->At(fTOFcluster[0]);
3036
3037     if (tofcl->GetNMatchableTracks() > 1)
3038       return kTRUE;
3039   }
3040   else if(fNtofClusters>0) AliInfo("No AliESDEvent available here!\n");
3041
3042   return kFALSE;
3043 }
3044
3045 void AliESDtrack::ReMapTOFcluster(Int_t ncl,Int_t *mapping){
3046   for(Int_t i=0;i<fNtofClusters;i++){
3047     if(fTOFcluster[i]<ncl && fTOFcluster[i]>-1)
3048       fTOFcluster[i] = mapping[fTOFcluster[i]];
3049     else
3050       AliInfo(Form("TOF cluster re-mapping in AliESDtrack: out of range (%i > %i)\n",fTOFcluster[i],ncl));
3051   }
3052 }
3053
3054 void AliESDtrack::SortTOFcluster(){
3055   TObjArray *tofclArray = GetESDEvent()->GetTOFcluster();
3056
3057   for(Int_t i=0;i<fNtofClusters-1;i++){
3058     for(Int_t j=i+1;j<fNtofClusters;j++){
3059       AliESDTOFcluster *tofcl = (AliESDTOFcluster *) tofclArray->At(fTOFcluster[i]);
3060       Int_t index1 = -1;
3061       for(Int_t it=0;it < tofcl->GetNMatchableTracks();it++){
3062          if(tofcl->GetTrackIndex(it) == GetID()) index1 = it;
3063       }
3064       Double_t timedist1 = 10000;
3065       for(Int_t isp=0; isp< AliPID::kSPECIESC;isp++){
3066         Double_t timec = TMath::Abs(tofcl->GetTime() - tofcl->GetIntegratedTime(isp));
3067         if(timec < timedist1) timedist1 = timec;
3068       }
3069       timedist1 *= 0.03; // in cm
3070       Double_t radius1 = tofcl->GetDx(index1)*tofcl->GetDx(index1) + tofcl->GetDz(index1)*tofcl->GetDz(index1) + timedist1*timedist1;
3071
3072       AliESDTOFcluster *tofcl2 = (AliESDTOFcluster *) tofclArray->At(fTOFcluster[j]);
3073       Int_t index2 = -1;
3074       for(Int_t it=0;it < tofcl2->GetNMatchableTracks();it++){
3075          if(tofcl2->GetTrackIndex(it) == GetID()) index2 = it;
3076       }
3077       if(index1 == -1 || index2 == -1){
3078       }
3079       Double_t timedist2 = 10000;
3080       for(Int_t isp=0; isp< AliPID::kSPECIESC;isp++){
3081         Double_t timec = TMath::Abs(tofcl2->GetTime() - tofcl2->GetIntegratedTime(isp));
3082         if(timec < timedist2) timedist2 = timec;
3083       }
3084       timedist2 *= 0.03; // in cm
3085       Double_t radius2 = tofcl2->GetDx(index2)*tofcl2->GetDx(index2) + tofcl2->GetDz(index2)*tofcl2->GetDz(index2) + timedist2*timedist2;
3086
3087       if(radius2 < radius1){
3088         Int_t change = fTOFcluster[i];
3089         fTOFcluster[i] = fTOFcluster[j];
3090         fTOFcluster[j] = change;
3091       }
3092     }
3093   }
3094 }