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