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