STEER/ESD/AliESDtrack.cxx

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