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