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
 119 #include "AliESDVertex.h"
 120 #include "AliESDtrack.h"
 121 #include "AliKalmanTrack.h"
 122 #include "AliLog.h"
 123 #include "AliTrackPointArray.h"
 124 #include "TPolyMarker3D.h"
 125
 126 ClassImp(AliESDtrack)
 127
 128 void SetPIDValues(Double_t * dest, const Double_t * src, Int_t n) {
 129   // This function copies "n" PID weights from "scr" to "dest"
 130   // and normalizes their sum to 1 thus producing conditional probabilities.
 131   // The negative weights are set to 0.
 132   // In case all the weights are non-positive they are replaced by
 133   // uniform probabilities
 134
 135   if (n<=0) return;
 136
 137   Float_t uniform = 1./(Float_t)n;
 138
 139   Float_t sum = 0;
 140   for (Int_t i=0; i<n; i++)
 141     if (src[i]>=0) {
 142       sum+=src[i];
 143       dest[i] = src[i];
 144     }
 145     else {
 146       dest[i] = 0;
 147     }
 148
 149   if(sum>0)
 150     for (Int_t i=0; i<n; i++) dest[i] /= sum;
 151   else
 152     for (Int_t i=0; i<n; i++) dest[i] = uniform;
 153 }
 154
 155 //_______________________________________________________________________
 156 AliESDtrack::AliESDtrack() :
 157   AliExternalTrackParam(),
 158   fCp(0),
 159   fIp(0),
 160   fTPCInner(0),
 161   fOp(0),
 162   fFriendTrack(new AliESDfriendTrack()),
 163   fTPCClusterMap(159),//number of padrows
 164   fTPCSharedMap(159),//number of padrows
 165   fFlags(0),
 166   fID(0),
 167   fLabel(0),
 168   fITSLabel(0),
 169   fTPCLabel(0),
 170   fTRDLabel(0),
 171   fTOFCalChannel(0),
 172   fTOFindex(-1),
 173   fHMPIDqn(0),
 174   fHMPIDcluIdx(0),
 175   fEMCALindex(kEMCALNoMatch),
 176   fHMPIDtrkTheta(0),
 177   fHMPIDtrkPhi(0),
 178   fHMPIDsignal(0),
 179   fTrackLength(0),
 180   fdTPC(0),fzTPC(0),
 181   fCddTPC(0),fCdzTPC(0),fCzzTPC(0),
 182   fCchi2TPC(0),
 183   fD(0),fZ(0),
 184   fCdd(0),fCdz(0),fCzz(0),
 185   fCchi2(0),
 186   fITSchi2(0),
 187   fTPCchi2(0),
 188   fTRDchi2(0),
 189   fTOFchi2(0),
 190   fHMPIDchi2(0),
 191   fITSsignal(0),
 192   fTPCsignal(0),
 193   fTPCsignalS(0),
 194   fTRDsignal(0),
 195   fTRDQuality(0),
 196   fTRDBudget(0),
 197   fTOFsignal(0),
 198   fTOFsignalToT(0),
 199   fTOFsignalRaw(0),
 200   fTOFsignalDz(0),
 201   fHMPIDtrkX(0),
 202   fHMPIDtrkY(0),
 203   fHMPIDmipX(0),
 204   fHMPIDmipY(0),
 205   fTPCncls(0),
 206   fTPCnclsF(0),
 207   fTPCsignalN(0),
 208   fITSncls(0),
 209   fITSClusterMap(0),
 210   fTRDncls(0),
 211   fTRDncls0(0),
 212   fTRDpidQuality(0),
 213   fTRDnSlices(0),
 214   fTRDslices(0x0)
 215
 216 {
 217   //
 218   // The default ESD constructor
 219   //
 220   Int_t i;
 221   for (i=0; i<AliPID::kSPECIES; i++) {
 222     fTrackTime[i]=0.;
 223     fR[i]=0.;
 224     fITSr[i]=0.;
 225     fTPCr[i]=0.;
 226     fTRDr[i]=0.;
 227     fTOFr[i]=0.;
 228     fHMPIDr[i]=0.;
 229   }
 230
 231   for (i=0; i<3; i++)   { fKinkIndexes[i]=0;}
 232   for (i=0; i<3; i++)   { fV0Indexes[i]=0;}
 233   for (i=0;i<kTRDnPlanes;i++) {
 234     fTRDTimBin[i]=0;
 235   }
 236   for (i=0;i<4;i++) {fTPCPoints[i]=0;}
 237   for (i=0;i<3;i++) {fTOFLabel[i]=0;}
 238   for (i=0;i<10;i++) {fTOFInfo[i]=0;}
 239   for (i=0;i<12;i++) {fITSModule[i]=-1;}
 240 }
 241
 242 //_______________________________________________________________________
 243 AliESDtrack::AliESDtrack(const AliESDtrack& track):
 244   AliExternalTrackParam(track),
 245   fCp(0),
 246   fIp(0),
 247   fTPCInner(0),
 248   fOp(0),
 249   fFriendTrack(0),
 250   fTPCClusterMap(track.fTPCClusterMap),
 251   fTPCSharedMap(track.fTPCSharedMap),
 252   fFlags(track.fFlags),
 253   fID(track.fID),
 254   fLabel(track.fLabel),
 255   fITSLabel(track.fITSLabel),
 256   fTPCLabel(track.fTPCLabel),
 257   fTRDLabel(track.fTRDLabel),
 258   fTOFCalChannel(track.fTOFCalChannel),
 259   fTOFindex(track.fTOFindex),
 260   fHMPIDqn(track.fHMPIDqn),
 261   fHMPIDcluIdx(track.fHMPIDcluIdx),
 262   fEMCALindex(track.fEMCALindex),
 263   fHMPIDtrkTheta(track.fHMPIDtrkTheta),
 264   fHMPIDtrkPhi(track.fHMPIDtrkPhi),
 265   fHMPIDsignal(track.fHMPIDsignal),
 266   fTrackLength(track.fTrackLength),
 267   fdTPC(track.fdTPC),fzTPC(track.fzTPC),
 268   fCddTPC(track.fCddTPC),fCdzTPC(track.fCdzTPC),fCzzTPC(track.fCzzTPC),
 269   fCchi2TPC(track.fCchi2TPC),
 270   fD(track.fD),fZ(track.fZ),
 271   fCdd(track.fCdd),fCdz(track.fCdz),fCzz(track.fCzz),
 272   fCchi2(track.fCchi2),
 273   fITSchi2(track.fITSchi2),
 274   fTPCchi2(track.fTPCchi2),
 275   fTRDchi2(track.fTRDchi2),
 276   fTOFchi2(track.fTOFchi2),
 277   fHMPIDchi2(track.fHMPIDchi2),
 278   fITSsignal(track.fITSsignal),
 279   fTPCsignal(track.fTPCsignal),
 280   fTPCsignalS(track.fTPCsignalS),
 281   fTRDsignal(track.fTRDsignal),
 282   fTRDQuality(track.fTRDQuality),
 283   fTRDBudget(track.fTRDBudget),
 284   fTOFsignal(track.fTOFsignal),
 285   fTOFsignalToT(track.fTOFsignalToT),
 286   fTOFsignalRaw(track.fTOFsignalRaw),
 287   fTOFsignalDz(track.fTOFsignalDz),
 288   fHMPIDtrkX(track.fHMPIDtrkX),
 289   fHMPIDtrkY(track.fHMPIDtrkY),
 290   fHMPIDmipX(track.fHMPIDmipX),
 291   fHMPIDmipY(track.fHMPIDmipY),
 292   fTPCncls(track.fTPCncls),
 293   fTPCnclsF(track.fTPCnclsF),
 294   fTPCsignalN(track.fTPCsignalN),
 295   fITSncls(track.fITSncls),
 296   fITSClusterMap(track.fITSClusterMap),
 297   fTRDncls(track.fTRDncls),
 298   fTRDncls0(track.fTRDncls0),
 299   fTRDpidQuality(track.fTRDpidQuality),
 300   fTRDnSlices(track.fTRDnSlices),
 301   fTRDslices(0x0)
 302 {
 303   //
 304   //copy constructor
 305   //
 306   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTrackTime[i]=track.fTrackTime[i];
 307   for (Int_t i=0;i<AliPID::kSPECIES;i++)  fR[i]=track.fR[i];
 308   //
 309   for (Int_t i=0;i<AliPID::kSPECIES;i++) fITSr[i]=track.fITSr[i];
 310   //
 311   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTPCr[i]=track.fTPCr[i];
 312   for (Int_t i=0;i<4;i++) {fTPCPoints[i]=track.fTPCPoints[i];}
 313   for (Int_t i=0; i<3;i++)   { fKinkIndexes[i]=track.fKinkIndexes[i];}
 314   for (Int_t i=0; i<3;i++)   { fV0Indexes[i]=track.fV0Indexes[i];}
 315   //
 316   for (Int_t i=0;i<kTRDnPlanes;i++) {
 317     fTRDTimBin[i]=track.fTRDTimBin[i];
 318   }
 319
 320   if (fTRDnSlices) {
 321     fTRDslices=new Double32_t[fTRDnSlices];
 322     for (Int_t i=0; i<fTRDnSlices; i++) fTRDslices[i]=track.fTRDslices[i];
 323   }
 324
 325   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTRDr[i]=track.fTRDr[i];
 326   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTOFr[i]=track.fTOFr[i];
 327   for (Int_t i=0;i<3;i++) fTOFLabel[i]=track.fTOFLabel[i];
 328   for (Int_t i=0;i<10;i++) fTOFInfo[i]=track.fTOFInfo[i];
 329   for (Int_t i=0;i<12;i++) fITSModule[i]=track.fITSModule[i];
 330   for (Int_t i=0;i<AliPID::kSPECIES;i++) fHMPIDr[i]=track.fHMPIDr[i];
 331
 332   if (track.fCp) fCp=new AliExternalTrackParam(*track.fCp);
 333   if (track.fIp) fIp=new AliExternalTrackParam(*track.fIp);
 334   if (track.fTPCInner) fTPCInner=new AliExternalTrackParam(*track.fTPCInner);
 335   if (track.fOp) fOp=new AliExternalTrackParam(*track.fOp);
 336
 337   if (track.fFriendTrack) fFriendTrack=new AliESDfriendTrack(*(track.fFriendTrack));
 338 }
 339
 340 //_______________________________________________________________________
 341 AliESDtrack::AliESDtrack(TParticle * part) :
 342   AliExternalTrackParam(),
 343   fCp(0),
 344   fIp(0),
 345   fTPCInner(0),
 346   fOp(0),
 347   fFriendTrack(0),
 348   fTPCClusterMap(159),//number of padrows
 349   fTPCSharedMap(159),//number of padrows
 350   fFlags(0),
 351   fID(0),
 352   fLabel(0),
 353   fITSLabel(0),
 354   fTPCLabel(0),
 355   fTRDLabel(0),
 356   fTOFCalChannel(0),
 357   fTOFindex(-1),
 358   fHMPIDqn(0),
 359   fHMPIDcluIdx(0),
 360   fEMCALindex(kEMCALNoMatch),
 361   fHMPIDtrkTheta(0),
 362   fHMPIDtrkPhi(0),
 363   fHMPIDsignal(0),
 364   fTrackLength(0),
 365   fdTPC(0),fzTPC(0),
 366   fCddTPC(0),fCdzTPC(0),fCzzTPC(0),
 367   fCchi2TPC(0),
 368   fD(0),fZ(0),
 369   fCdd(0),fCdz(0),fCzz(0),
 370   fCchi2(0),
 371   fITSchi2(0),
 372   fTPCchi2(0),
 373   fTRDchi2(0),
 374   fTOFchi2(0),
 375   fHMPIDchi2(0),
 376   fITSsignal(0),
 377   fTPCsignal(0),
 378   fTPCsignalS(0),
 379   fTRDsignal(0),
 380   fTRDQuality(0),
 381   fTRDBudget(0),
 382   fTOFsignal(0),
 383   fTOFsignalToT(0),
 384   fTOFsignalRaw(0),
 385   fTOFsignalDz(0),
 386   fHMPIDtrkX(0),
 387   fHMPIDtrkY(0),
 388   fHMPIDmipX(0),
 389   fHMPIDmipY(0),
 390   fTPCncls(0),
 391   fTPCnclsF(0),
 392   fTPCsignalN(0),
 393   fITSncls(0),
 394   fITSClusterMap(0),
 395   fTRDncls(0),
 396   fTRDncls0(0),
 397   fTRDpidQuality(0),
 398   fTRDnSlices(0),
 399   fTRDslices(0x0)
 400 {
 401   //
 402   // ESD track from TParticle
 403   //
 404
 405   // Reset all the arrays
 406   Int_t i;
 407   for (i=0; i<AliPID::kSPECIES; i++) {
 408     fTrackTime[i]=0.;
 409     fR[i]=0.;
 410     fITSr[i]=0.;
 411     fTPCr[i]=0.;
 412     fTRDr[i]=0.;
 413     fTOFr[i]=0.;
 414     fHMPIDr[i]=0.;
 415   }
 416
 417   for (i=0; i<3; i++)   { fKinkIndexes[i]=0;}
 418   for (i=0; i<3; i++)   { fV0Indexes[i]=-1;}
 419   for (i=0;i<kTRDnPlanes;i++) {
 420     fTRDTimBin[i]=0;
 421   }
 422   for (i=0;i<4;i++) {fTPCPoints[i]=0;}
 423   for (i=0;i<3;i++) {fTOFLabel[i]=0;}
 424   for (i=0;i<10;i++) {fTOFInfo[i]=0;}
 425   for (i=0;i<12;i++) {fITSModule[i]=-1;}
 426
 427   // Calculate the AliExternalTrackParam content
 428
 429   Double_t xref;
 430   Double_t alpha;
 431   Double_t param[5];
 432   Double_t covar[15];
 433
 434   // Calculate alpha: the rotation angle of the corresponding local system (TPC sector)
 435   alpha = part->Phi()*180./TMath::Pi();
 436   if (alpha<0) alpha+= 360.;
 437   if (alpha>360) alpha -= 360.;
 438
 439   Int_t sector = (Int_t)(alpha/20.);
 440   alpha = 10. + 20.*sector;
 441   alpha /= 180;
 442   alpha *= TMath::Pi();
 443
 444   // Covariance matrix: no errors, the parameters are exact
 445   for (i=0; i<15; i++) covar[i]=0.;
 446
 447   // Get the vertex of origin and the momentum
 448   TVector3 ver(part->Vx(),part->Vy(),part->Vz());
 449   TVector3 mom(part->Px(),part->Py(),part->Pz());
 450
 451   // Rotate to the local coordinate system (TPC sector)
 452   ver.RotateZ(-alpha);
 453   mom.RotateZ(-alpha);
 454
 455   // X of the referense plane
 456   xref = ver.X();
 457
 458   Int_t pdgCode = part->GetPdgCode();
 459
 460   Double_t charge =
 461     TDatabasePDG::Instance()->GetParticle(pdgCode)->Charge();
 462
 463   param[0] = ver.Y();
 464   param[1] = ver.Z();
 465   param[2] = TMath::Sin(mom.Phi());
 466   param[3] = mom.Pz()/mom.Pt();
 467   param[4] = TMath::Sign(1/mom.Pt(),charge);
 468
 469   // Set AliExternalTrackParam
 470   Set(xref, alpha, param, covar);
 471
 472   // Set the PID
 473   Int_t indexPID = 99;
 474
 475   switch (TMath::Abs(pdgCode)) {
 476
 477   case  11: // electron
 478     indexPID = 0;
 479     break;
 480
 481   case 13: // muon
 482     indexPID = 1;
 483     break;
 484
 485   case 211: // pion
 486     indexPID = 2;
 487     break;
 488
 489   case 321: // kaon
 490     indexPID = 3;
 491     break;
 492
 493   case 2212: // proton
 494     indexPID = 4;
 495     break;
 496
 497   default:
 498     break;
 499   }
 500
 501   // If the particle is not e,mu,pi,K or p the PID probabilities are set to 0
 502   if (indexPID < AliPID::kSPECIES) {
 503     fR[indexPID]=1.;
 504     fITSr[indexPID]=1.;
 505     fTPCr[indexPID]=1.;
 506     fTRDr[indexPID]=1.;
 507     fTOFr[indexPID]=1.;
 508     fHMPIDr[indexPID]=1.;
 509
 510   }
 511   // AliESD track label
 512   SetLabel(part->GetUniqueID());
 513
 514 }
 515
 516 //_______________________________________________________________________
 517 AliESDtrack::~AliESDtrack(){
 518   //
 519   // This is destructor according Coding Conventrions
 520   //
 521   //printf("Delete track\n");
 522   delete fIp;
 523   delete fTPCInner;
 524   delete fOp;
 525   delete fCp;
 526   delete fFriendTrack;
 527   if(fTRDnSlices)
 528     delete[] fTRDslices;
 529 }
 530
 531 AliESDtrack &AliESDtrack::operator=(const AliESDtrack &source){
 532
 533
 534   if(&source == this) return *this;
 535   AliExternalTrackParam::operator=(source);
 536
 537
 538   if(source.fCp){
 539     // we have the trackparam: assign or copy construct
 540     if(fCp)*fCp = *source.fCp;
 541     else fCp = new AliExternalTrackParam(*source.fCp);
 542   }
 543   else{
 544     // no track param delete the old one
 545     if(fCp)delete fCp;
 546     fCp = 0;
 547   }
 548
 549   if(source.fIp){
 550     // we have the trackparam: assign or copy construct
 551     if(fIp)*fIp = *source.fIp;
 552     else fIp = new AliExternalTrackParam(*source.fIp);
 553   }
 554   else{
 555     // no track param delete the old one
 556     if(fIp)delete fIp;
 557     fIp = 0;
 558   }
 559
 560
 561   if(source.fTPCInner){
 562     // we have the trackparam: assign or copy construct
 563     if(fTPCInner) *fTPCInner = *source.fTPCInner;
 564     else fTPCInner = new AliExternalTrackParam(*source.fTPCInner);
 565   }
 566   else{
 567     // no track param delete the old one
 568     if(fTPCInner)delete fTPCInner;
 569     fTPCInner = 0;
 570   }
 571
 572
 573   if(source.fOp){
 574     // we have the trackparam: assign or copy construct
 575     if(fOp) *fOp = *source.fOp;
 576     else fOp = new AliExternalTrackParam(*source.fOp);
 577   }
 578   else{
 579     // no track param delete the old one
 580     if(fOp)delete fOp;
 581     fOp = 0;
 582   }
 583
 584   // copy also the friend track
 585   // use copy constructor
 586   if(source.fFriendTrack){
 587     // we have the trackparam: assign or copy construct
 588     delete fFriendTrack; fFriendTrack=new AliESDfriendTrack(*source.fFriendTrack);
 589   }
 590   else{
 591     // no track param delete the old one
 592     delete fFriendTrack; fFriendTrack= 0;
 593   }
 594
 595   fTPCClusterMap = source.fTPCClusterMap;
 596   fTPCSharedMap  = source.fTPCSharedMap;
 597   // the simple stuff
 598   fFlags    = source.fFlags;
 599   fID       = source.fID;
 600   fLabel    = source.fLabel;
 601   fITSLabel = source.fITSLabel;
 602   for(int i = 0; i< 12;++i){
 603     fITSModule[i] = source.fITSModule[i];
 604   }
 605   fTPCLabel = source.fTPCLabel;
 606   fTRDLabel = source.fTRDLabel;
 607   for(int i = 0; i< 3;++i){
 608     fTOFLabel[i] = source.fTOFLabel[i];
 609   }
 610   fTOFCalChannel = source.fTOFCalChannel;
 611   fTOFindex      = source.fTOFindex;
 612   fHMPIDqn       = source.fHMPIDqn;
 613   fHMPIDcluIdx   = source.fHMPIDcluIdx;
 614   fEMCALindex    = source.fEMCALindex;
 615
 616   for(int i = 0; i< 3;++i){
 617     fKinkIndexes[i] = source.fKinkIndexes[i];
 618     fV0Indexes[i]   = source.fV0Indexes[i];
 619   }
 620
 621   for(int i = 0; i< AliPID::kSPECIES;++i){
 622     fR[i]     = source.fR[i];
 623     fITSr[i]  = source.fITSr[i];
 624     fTPCr[i]  = source.fTPCr[i];
 625     fTRDr[i]  = source.fTRDr[i];
 626     fTOFr[i]  = source.fTOFr[i];
 627     fHMPIDr[i] = source.fHMPIDr[i];
 628     fTrackTime[i] = source.fTrackTime[i];
 629   }
 630
 631   fHMPIDtrkTheta = source.fHMPIDtrkTheta;
 632   fHMPIDtrkPhi   = source.fHMPIDtrkPhi;
 633   fHMPIDsignal   = source.fHMPIDsignal;
 634
 635
 636   fTrackLength   = source. fTrackLength;
 637   fdTPC  = source.fdTPC;
 638   fzTPC  = source.fzTPC;
 639   fCddTPC = source.fCddTPC;
 640   fCdzTPC = source.fCdzTPC;
 641   fCzzTPC = source.fCzzTPC;
 642   fCchi2TPC = source.fCchi2TPC;
 643
 644   fD  = source.fD;
 645   fZ  = source.fZ;
 646   fCdd = source.fCdd;
 647   fCdz = source.fCdz;
 648   fCzz = source.fCzz;
 649   fCchi2     = source.fCchi2;
 650
 651   fITSchi2   = source.fITSchi2;
 652   fTPCchi2   = source.fTPCchi2;
 653   fTRDchi2   = source.fTRDchi2;
 654   fTOFchi2   = source.fTOFchi2;
 655   fHMPIDchi2 = source.fHMPIDchi2;
 656
 657
 658   fITSsignal  = source.fITSsignal;
 659   fTPCsignal  = source.fTPCsignal;
 660   fTPCsignalS = source.fTPCsignalS;
 661   for(int i = 0; i< 4;++i){
 662     fTPCPoints[i] = source.fTPCPoints[i];
 663   }
 664   fTRDsignal = source.fTRDsignal;
 665
 666   for(int i = 0;i < kTRDnPlanes;++i){
 667     fTRDTimBin[i] = source.fTRDTimBin[i];
 668   }
 669
 670   if(fTRDnSlices)
 671     delete[] fTRDslices;
 672   fTRDslices=0;
 673   fTRDnSlices=source.fTRDnSlices;
 674   if (fTRDnSlices) {
 675     fTRDslices=new Double32_t[fTRDnSlices];
 676     for(int j = 0;j < fTRDnSlices;++j) fTRDslices[j] = source.fTRDslices[j];
 677   }
 678
 679   fTRDQuality =   source.fTRDQuality;
 680   fTRDBudget  =   source.fTRDBudget;
 681   fTOFsignal  =   source.fTOFsignal;
 682   fTOFsignalToT = source.fTOFsignalToT;
 683   fTOFsignalRaw = source.fTOFsignalRaw;
 684   fTOFsignalDz  = source.fTOFsignalDz;
 685
 686   for(int i = 0;i<10;++i){
 687     fTOFInfo[i] = source.fTOFInfo[i];
 688   }
 689
 690   fHMPIDtrkX = source.fHMPIDtrkX;
 691   fHMPIDtrkY = source.fHMPIDtrkY;
 692   fHMPIDmipX = source.fHMPIDmipX;
 693   fHMPIDmipY = source.fHMPIDmipY;
 694
 695   fTPCncls    = source.fTPCncls;
 696   fTPCnclsF   = source.fTPCnclsF;
 697   fTPCsignalN = source.fTPCsignalN;
 698
 699   fITSncls = source.fITSncls;
 700   fITSClusterMap = source.fITSClusterMap;
 701   fTRDncls   = source.fTRDncls;
 702   fTRDncls0  = source.fTRDncls0;
 703   fTRDpidQuality  = source.fTRDpidQuality;
 704   return *this;
 705 }
 706
 707
 708
 709 void AliESDtrack::Copy(TObject &obj) const {
 710
 711   // this overwrites the virtual TOBject::Copy()
 712   // to allow run time copying without casting
 713   // in AliESDEvent
 714
 715   if(this==&obj)return;
 716   AliESDtrack *robj = dynamic_cast<AliESDtrack*>(&obj);
 717   if(!robj)return; // not an AliESDtrack
 718   *robj = *this;
 719
 720 }
 721
 722
 723
 724 void AliESDtrack::AddCalibObject(TObject * object){
 725   //
 726   // add calib object to the list
 727   //
 728   if (!fFriendTrack) fFriendTrack  = new AliESDfriendTrack;
 729   fFriendTrack->AddCalibObject(object);
 730 }
 731
 732 TObject *  AliESDtrack::GetCalibObject(Int_t index){
 733   //
 734   // return calib objct at given position
 735   //
 736   if (!fFriendTrack) return 0;
 737   return fFriendTrack->GetCalibObject(index);
 738 }
 739
 740
 741 const Bool_t AliESDtrack::FillTPCOnlyTrack(AliESDtrack &track){
 742
 743   // Fills the information of the TPC-only first reconstruction pass
 744   // into the passed ESDtrack object. For consistency fTPCInner is also filled
 745   // again
 746
 747
 748
 749   // For data produced before r26675
 750   // RelateToVertexTPC was not properly called during reco
 751   // so you'll have to call it again, before FillTPCOnlyTrack
 752   //  Float_t p[2],cov[3];
 753   // track->GetImpactParametersTPC(p,cov);
 754   // if(p[0]==0&&p[1]==0) // <- Default values
 755   //  track->RelateToVertexTPC(esd->GetPrimaryVertexTPC(),esd->GetMagneticField(),kVeryBig);
 756
 757
 758   if(!fTPCInner)return kFALSE;
 759
 760   // fill the TPC track params to the global track parameters
 761   track.Set(fTPCInner->GetX(),fTPCInner->GetAlpha(),fTPCInner->GetParameter(),fTPCInner->GetCovariance());
 762   track.fD = fdTPC;
 763   track.fZ = fzTPC;
 764   track.fCdd = fCddTPC;
 765   track.fCdz = fCdzTPC;
 766   track.fCzz = fCzzTPC;
 767
 768   // copy the TPCinner parameters
 769   if(track.fTPCInner) *track.fTPCInner = *fTPCInner;
 770   else track.fTPCInner = new AliExternalTrackParam(*fTPCInner);
 771   track.fdTPC   = fdTPC;
 772   track.fzTPC   = fzTPC;
 773   track.fCddTPC = fCddTPC;
 774   track.fCdzTPC = fCdzTPC;
 775   track.fCzzTPC = fCzzTPC;
 776   track.fCchi2TPC = fCchi2TPC;
 777
 778
 779   // copy all other TPC specific parameters
 780
 781   // replace label by TPC label
 782   track.fLabel    = fTPCLabel;
 783   track.fTPCLabel = fTPCLabel;
 784
 785   track.fTPCchi2 = fTPCchi2;
 786   track.fTPCsignal = fTPCsignal;
 787   track.fTPCsignalS = fTPCsignalS;
 788   for(int i = 0;i<4;++i)track.fTPCPoints[i] = fTPCPoints[i];
 789
 790   track.fTPCncls    = fTPCncls;
 791   track.fTPCnclsF   = fTPCnclsF;
 792   track.fTPCsignalN =  fTPCsignalN;
 793
 794   // PID
 795   for(int i=0;i<AliPID::kSPECIES;++i){
 796     track.fTPCr[i] = fTPCr[i];
 797     // combined PID is TPC only!
 798     track.fR[i] = fTPCr[i];
 799   }
 800   track.fTPCClusterMap = fTPCClusterMap;
 801   track.fTPCSharedMap = fTPCSharedMap;
 802
 803
 804   // reset the flags
 805   track.fFlags = kTPCin;
 806   track.fID    = fID;
 807
 808
 809   for (Int_t i=0;i<3;i++) track.fKinkIndexes[i] = fKinkIndexes[i];
 810
 811   return kTRUE;
 812
 813 }
 814
 815 //_______________________________________________________________________
 816 void AliESDtrack::MakeMiniESDtrack(){
 817   // Resets everything except
 818   // fFlags: Reconstruction status flags
 819   // fLabel: Track label
 820   // fID:  Unique ID of the track
 821   // Impact parameter information
 822   // fR[AliPID::kSPECIES]: combined "detector response probability"
 823   // Running track parameters in the base class (AliExternalTrackParam)
 824
 825   fTrackLength = 0;
 826
 827   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTrackTime[i] = 0;
 828
 829   // Reset track parameters constrained to the primary vertex
 830   delete fCp;fCp = 0;
 831
 832   // Reset track parameters at the inner wall of TPC
 833   delete fIp;fIp = 0;
 834   delete fTPCInner;fTPCInner=0;
 835   // Reset track parameters at the inner wall of the TRD
 836   delete fOp;fOp = 0;
 837
 838
 839   // Reset ITS track related information
 840   fITSchi2 = 0;
 841   fITSncls = 0;
 842   fITSClusterMap=0;
 843   fITSsignal = 0;
 844   for (Int_t i=0;i<AliPID::kSPECIES;i++) fITSr[i]=0;
 845   fITSLabel = 0;
 846
 847   // Reset TPC related track information
 848   fTPCchi2 = 0;
 849   fTPCncls = 0;
 850   fTPCnclsF = 0;
 851   fTPCClusterMap = 0;
 852   fTPCSharedMap = 0;
 853   fTPCsignal= 0;
 854   fTPCsignalS= 0;
 855   fTPCsignalN= 0;
 856   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTPCr[i]=0;
 857   fTPCLabel=0;
 858   for (Int_t i=0;i<4;i++) fTPCPoints[i] = 0;
 859   for (Int_t i=0; i<3;i++)   fKinkIndexes[i] = 0;
 860   for (Int_t i=0; i<3;i++)   fV0Indexes[i] = 0;
 861
 862   // Reset TRD related track information
 863   fTRDchi2 = 0;
 864   fTRDncls = 0;
 865   fTRDncls0 = 0;
 866   fTRDsignal = 0;
 867   for (Int_t i=0;i<kTRDnPlanes;i++) {
 868     fTRDTimBin[i]  = 0;
 869   }
 870   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTRDr[i] = 0;
 871   fTRDLabel = 0;
 872   fTRDQuality  = 0;
 873   fTRDpidQuality = 0;
 874   if(fTRDnSlices)
 875     delete[] fTRDslices;
 876   fTRDslices=0x0;
 877   fTRDnSlices=0;
 878   fTRDBudget  = 0;
 879
 880   // Reset TOF related track information
 881   fTOFchi2 = 0;
 882   fTOFindex = -1;
 883   fTOFsignal = 0;
 884   fTOFCalChannel = 0;
 885   fTOFsignalToT = 0;
 886   fTOFsignalRaw = 0;
 887   fTOFsignalDz = 0;
 888   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTOFr[i] = 0;
 889   for (Int_t i=0;i<3;i++) fTOFLabel[i] = 0;
 890   for (Int_t i=0;i<10;i++) fTOFInfo[i] = 0;
 891
 892   // Reset HMPID related track information
 893   fHMPIDchi2 = 0;
 894   fHMPIDqn = 0;
 895   fHMPIDcluIdx = 0;
 896   fHMPIDsignal = 0;
 897   for (Int_t i=0;i<AliPID::kSPECIES;i++) fHMPIDr[i] = 0;
 898   fHMPIDtrkTheta = 0;
 899   fHMPIDtrkPhi = 0;
 900   fHMPIDtrkX = 0;
 901   fHMPIDtrkY = 0;
 902   fHMPIDmipX = 0;
 903   fHMPIDmipY = 0;
 904   fEMCALindex = kEMCALNoMatch;
 905
 906   delete fFriendTrack; fFriendTrack = 0;
 907 }
 908 //_______________________________________________________________________
 909 Double_t AliESDtrack::GetMass() const {
 910   // Returns the mass of the most probable particle type
 911   Float_t max=0.;
 912   Int_t k=-1;
 913   for (Int_t i=0; i<AliPID::kSPECIES; i++) {
 914     if (fR[i]>max) {k=i; max=fR[i];}
 915   }
 916   if (k==0) { // dE/dx "crossing points" in the TPC
 917      Double_t p=GetP();
 918      if ((p>0.38)&&(p<0.48))
 919         if (fR[0]<fR[3]*10.) return AliPID::ParticleMass(AliPID::kKaon);
 920      if ((p>0.75)&&(p<0.85))
 921         if (fR[0]<fR[4]*10.) return AliPID::ParticleMass(AliPID::kProton);
 922      return 0.00051;
 923   }
 924   if (k==1) return AliPID::ParticleMass(AliPID::kMuon);
 925   if (k==2||k==-1) return AliPID::ParticleMass(AliPID::kPion);
 926   if (k==3) return AliPID::ParticleMass(AliPID::kKaon);
 927   if (k==4) return AliPID::ParticleMass(AliPID::kProton);
 928   AliWarning("Undefined mass !");
 929   return AliPID::ParticleMass(AliPID::kPion);
 930 }
 931
 932 //______________________________________________________________________________
 933 Double_t AliESDtrack::E() const
 934 {
 935   // Returns the energy of the particle given its assumed mass.
 936   // Assumes the pion mass if the particle can't be identified properly.
 937
 938   Double_t m = M();
 939   Double_t p = P();
 940   return TMath::Sqrt(p*p + m*m);
 941 }
 942
 943 //______________________________________________________________________________
 944 Double_t AliESDtrack::Y() const
 945 {
 946   // Returns the rapidity of a particle given its assumed mass.
 947   // Assumes the pion mass if the particle can't be identified properly.
 948
 949   Double_t e = E();
 950   Double_t pz = Pz();
 951   if (e != TMath::Abs(pz)) { // energy was not equal to pz
 952     return 0.5*TMath::Log((e+pz)/(e-pz));
 953   } else { // energy was equal to pz
 954     return -999.;
 955   }
 956 }
 957
 958 //_______________________________________________________________________
 959 Bool_t AliESDtrack::UpdateTrackParams(const AliKalmanTrack *t, ULong_t flags){
 960   //
 961   // This function updates track's running parameters
 962   //
 963   Int_t *index=0;
 964   Bool_t rc=kTRUE;
 965
 966   SetStatus(flags);
 967   fLabel=t->GetLabel();
 968
 969   if (t->IsStartedTimeIntegral()) {
 970     SetStatus(kTIME);
 971     Double_t times[10];t->GetIntegratedTimes(times); SetIntegratedTimes(times);
 972     SetIntegratedLength(t->GetIntegratedLength());
 973   }
 974
 975   Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
 976
 977   switch (flags) {
 978
 979   case kITSin: case kITSout: case kITSrefit:
 980     fITSClusterMap=0;
 981     fITSncls=t->GetNumberOfClusters();
 982     index=fFriendTrack->GetITSindices();
 983     for (Int_t i=0;i<AliESDfriendTrack::kMaxITScluster;i++) {
 984         index[i]=t->GetClusterIndex(i);
 985         if (i<fITSncls) {
 986            Int_t l=(index[i] & 0xf0000000) >> 28;
 987            SETBIT(fITSClusterMap,l);
 988         }
 989     }
 990     fITSchi2=t->GetChi2();
 991     fITSsignal=t->GetPIDsignal();
 992     fITSLabel = t->GetLabel();
 993     // keep in fOp the parameters outside ITS for ITS stand-alone tracks
 994     if (flags==kITSout) {
 995       if (!fOp) fOp=new AliExternalTrackParam(*t);
 996       else
 997         fOp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
 998     }
 999     break;
1000
1001   case kTPCin: case kTPCrefit:
1002     fTPCLabel = t->GetLabel();
1003     if (flags==kTPCin)  fTPCInner=new AliExternalTrackParam(*t);
1004     if (!fIp) fIp=new AliExternalTrackParam(*t);
1005     else
1006       fIp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
1007   case kTPCout:
1008     index=fFriendTrack->GetTPCindices();
1009     if (flags & kTPCout){
1010       if (!fOp) fOp=new AliExternalTrackParam(*t);
1011       else
1012         fOp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
1013     }
1014     fTPCncls=t->GetNumberOfClusters();
1015     fTPCchi2=t->GetChi2();
1016
1017      {//prevrow must be declared in separate namespace, otherwise compiler cries:
1018       //"jump to case label crosses initialization of `Int_t prevrow'"
1019        Int_t prevrow = -1;
1020        //       for (Int_t i=0;i<fTPCncls;i++)
1021        for (Int_t i=0;i<AliESDfriendTrack::kMaxTPCcluster;i++)
1022         {
1023           index[i]=t->GetClusterIndex(i);
1024           Int_t idx = index[i];
1025
1026           if (idx<0) continue;
1027
1028           // Piotr's Cluster Map for HBT
1029           // ### please change accordingly if cluster array is changing
1030           // to "New TPC Tracking" style (with gaps in array)
1031           Int_t sect = (idx&0xff000000)>>24;
1032           Int_t row = (idx&0x00ff0000)>>16;
1033           if (sect > 18) row +=63; //if it is outer sector, add number of inner sectors
1034
1035           fTPCClusterMap.SetBitNumber(row,kTRUE);
1036
1037           //Fill the gap between previous row and this row with 0 bits
1038           //In case  ###  pleas change it as well - just set bit 0 in case there
1039           //is no associated clusters for current "i"
1040           if (prevrow < 0)
1041            {
1042              prevrow = row;//if previous bit was not assigned yet == this is the first one
1043            }
1044           else
1045            { //we don't know the order (inner to outer or reverse)
1046              //just to be save in case it is going to change
1047              Int_t n = 0, m = 0;
1048              if (prevrow < row)
1049               {
1050                 n = prevrow;
1051                 m = row;
1052               }
1053              else
1054               {
1055                 n = row;
1056                 m = prevrow;
1057               }
1058
1059              for (Int_t j = n+1; j < m; j++)
1060               {
1061                 fTPCClusterMap.SetBitNumber(j,kFALSE);
1062               }
1063              prevrow = row;
1064            }
1065           // End Of Piotr's Cluster Map for HBT
1066         }
1067      }
1068     fTPCsignal=t->GetPIDsignal();
1069     break;
1070
1071   case kTRDout: case kTRDin: case kTRDrefit:
1072     index     = fFriendTrack->GetTRDindices();
1073     fTRDLabel = t->GetLabel();
1074     fTRDchi2  = t->GetChi2();
1075     fTRDncls  = t->GetNumberOfClusters();
1076     for (Int_t i=0;i<6;i++) index[i]=t->GetTrackletIndex(i);
1077
1078     fTRDsignal=t->GetPIDsignal();
1079     break;
1080   case kTRDbackup:
1081     if (!fOp) fOp=new AliExternalTrackParam(*t);
1082     else
1083       fOp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
1084     fTRDncls0 = t->GetNumberOfClusters();
1085     break;
1086   case kTOFin:
1087     break;
1088   case kTOFout:
1089     break;
1090   case kTRDStop:
1091     break;
1092   default:
1093     AliError("Wrong flag !");
1094     return kFALSE;
1095   }
1096
1097   return rc;
1098 }
1099
1100 //_______________________________________________________________________
1101 void AliESDtrack::GetExternalParameters(Double_t &x, Double_t p[5]) const {
1102   //---------------------------------------------------------------------
1103   // This function returns external representation of the track parameters
1104   //---------------------------------------------------------------------
1105   x=GetX();
1106   for (Int_t i=0; i<5; i++) p[i]=GetParameter()[i];
1107 }
1108
1109 //_______________________________________________________________________
1110 void AliESDtrack::GetExternalCovariance(Double_t cov[15]) const {
1111   //---------------------------------------------------------------------
1112   // This function returns external representation of the cov. matrix
1113   //---------------------------------------------------------------------
1114   for (Int_t i=0; i<15; i++) cov[i]=AliExternalTrackParam::GetCovariance()[i];
1115 }
1116
1117 //_______________________________________________________________________
1118 Bool_t AliESDtrack::GetConstrainedExternalParameters
1119                  (Double_t &alpha, Double_t &x, Double_t p[5]) const {
1120   //---------------------------------------------------------------------
1121   // This function returns the constrained external track parameters
1122   //---------------------------------------------------------------------
1123   if (!fCp) return kFALSE;
1124   alpha=fCp->GetAlpha();
1125   x=fCp->GetX();
1126   for (Int_t i=0; i<5; i++) p[i]=fCp->GetParameter()[i];
1127   return kTRUE;
1128 }
1129
1130 //_______________________________________________________________________
1131 Bool_t
1132 AliESDtrack::GetConstrainedExternalCovariance(Double_t c[15]) const {
1133   //---------------------------------------------------------------------
1134   // This function returns the constrained external cov. matrix
1135   //---------------------------------------------------------------------
1136   if (!fCp) return kFALSE;
1137   for (Int_t i=0; i<15; i++) c[i]=fCp->GetCovariance()[i];
1138   return kTRUE;
1139 }
1140
1141 Bool_t
1142 AliESDtrack::GetInnerExternalParameters
1143                  (Double_t &alpha, Double_t &x, Double_t p[5]) const {
1144   //---------------------------------------------------------------------
1145   // This function returns external representation of the track parameters
1146   // at the inner layer of TPC
1147   //---------------------------------------------------------------------
1148   if (!fIp) return kFALSE;
1149   alpha=fIp->GetAlpha();
1150   x=fIp->GetX();
1151   for (Int_t i=0; i<5; i++) p[i]=fIp->GetParameter()[i];
1152   return kTRUE;
1153 }
1154
1155 Bool_t
1156 AliESDtrack::GetInnerExternalCovariance(Double_t cov[15]) const {
1157  //---------------------------------------------------------------------
1158  // This function returns external representation of the cov. matrix
1159  // at the inner layer of TPC
1160  //---------------------------------------------------------------------
1161   if (!fIp) return kFALSE;
1162   for (Int_t i=0; i<15; i++) cov[i]=fIp->GetCovariance()[i];
1163   return kTRUE;
1164 }
1165
1166 void
1167 AliESDtrack::SetOuterParam(const AliExternalTrackParam *p, ULong_t flags) {
1168   //
1169   // This is a direct setter for the outer track parameters
1170   //
1171   SetStatus(flags);
1172   if (fOp) delete fOp;
1173   fOp=new AliExternalTrackParam(*p);
1174 }
1175
1176 Bool_t
1177 AliESDtrack::GetOuterExternalParameters
1178                  (Double_t &alpha, Double_t &x, Double_t p[5]) const {
1179   //---------------------------------------------------------------------
1180   // This function returns external representation of the track parameters
1181   // at the inner layer of TRD
1182   //---------------------------------------------------------------------
1183   if (!fOp) return kFALSE;
1184   alpha=fOp->GetAlpha();
1185   x=fOp->GetX();
1186   for (Int_t i=0; i<5; i++) p[i]=fOp->GetParameter()[i];
1187   return kTRUE;
1188 }
1189
1190 Bool_t
1191 AliESDtrack::GetOuterExternalCovariance(Double_t cov[15]) const {
1192  //---------------------------------------------------------------------
1193  // This function returns external representation of the cov. matrix
1194  // at the inner layer of TRD
1195  //---------------------------------------------------------------------
1196   if (!fOp) return kFALSE;
1197   for (Int_t i=0; i<15; i++) cov[i]=fOp->GetCovariance()[i];
1198   return kTRUE;
1199 }
1200
1201 Int_t AliESDtrack::GetNcls(Int_t idet) const
1202 {
1203   // Get number of clusters by subdetector index
1204   //
1205   Int_t ncls = 0;
1206   switch(idet){
1207   case 0:
1208     ncls = fITSncls;
1209     break;
1210   case 1:
1211     ncls = fTPCncls;
1212     break;
1213   case 2:
1214     ncls = fTRDncls;
1215     break;
1216   case 3:
1217     if (fTOFindex != -1)
1218       ncls = 1;
1219     break;
1220   default:
1221     break;
1222   }
1223   return ncls;
1224 }
1225
1226 Int_t AliESDtrack::GetClusters(Int_t idet, Int_t *idx) const
1227 {
1228   // Get cluster index array by subdetector index
1229   //
1230   Int_t ncls = 0;
1231   switch(idet){
1232   case 0:
1233     ncls = GetITSclusters(idx);
1234     break;
1235   case 1:
1236     ncls = GetTPCclusters(idx);
1237     break;
1238   case 2:
1239     ncls = GetTRDclusters(idx);
1240     break;
1241   case 3:
1242     if (fTOFindex != -1) {
1243       idx[0] = fTOFindex;
1244       ncls = 1;
1245     }
1246     break;
1247   case 4: //PHOS
1248     break;
1249   case 5:
1250     if (fHMPIDcluIdx != 0) {
1251       idx[0] = GetHMPIDcluIdx();
1252       ncls = 1;
1253     }
1254     break;
1255   case 6: //EMCAL
1256     break;
1257   default:
1258     break;
1259   }
1260   return ncls;
1261 }
1262
1263 //_______________________________________________________________________
1264 void AliESDtrack::GetIntegratedTimes(Double_t *times) const {
1265   // Returns the array with integrated times for each particle hypothesis
1266   for (Int_t i=0; i<AliPID::kSPECIES; i++) times[i]=fTrackTime[i];
1267 }
1268
1269 //_______________________________________________________________________
1270 void AliESDtrack::SetIntegratedTimes(const Double_t *times) {
1271   // Sets the array with integrated times for each particle hypotesis
1272   for (Int_t i=0; i<AliPID::kSPECIES; i++) fTrackTime[i]=times[i];
1273 }
1274
1275 //_______________________________________________________________________
1276 void AliESDtrack::SetITSpid(const Double_t *p) {
1277   // Sets values for the probability of each particle type (in ITS)
1278   SetPIDValues(fITSr,p,AliPID::kSPECIES);
1279   SetStatus(AliESDtrack::kITSpid);
1280 }
1281
1282 //_______________________________________________________________________
1283 void AliESDtrack::GetITSpid(Double_t *p) const {
1284   // Gets the probability of each particle type (in ITS)
1285   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fITSr[i];
1286 }
1287
1288 //_______________________________________________________________________
1289 Char_t AliESDtrack::GetITSclusters(Int_t *idx) const {
1290   //---------------------------------------------------------------------
1291   // This function returns indices of the assgined ITS clusters
1292   //---------------------------------------------------------------------
1293   if (idx!=0) {
1294      Int_t *index=fFriendTrack->GetITSindices();
1295      for (Int_t i=0; i<AliESDfriendTrack::kMaxITScluster; i++) {
1296          if ( (i>=fITSncls) && (i<6) ) idx[i]=-1;
1297          else idx[i]=index[i];
1298      }
1299   }
1300   return fITSncls;
1301 }
1302
1303 //_______________________________________________________________________
1304 Bool_t AliESDtrack::GetITSModuleIndexInfo(Int_t ilayer,Int_t &idet,Int_t &status,
1305                                          Float_t &xloc,Float_t &zloc) const {
1306   //----------------------------------------------------------------------
1307   // This function encodes in the module number also the status of cluster association
1308   // "status" can have the following values:
1309   // 1 "found" (cluster is associated),
1310   // 2 "dead" (module is dead from OCDB),
1311   // 3 "skipped" (module or layer forced to be skipped),
1312   // 4 "outinz" (track out of z acceptance),
1313   // 5 "nocls" (no clusters in the road),
1314   // 6 "norefit" (cluster rejected during refit),
1315   // 7 "deadzspd" (holes in z in SPD)
1316   // Also given are the coordinates of the crossing point of track and module
1317   // (in the local module ref. system)
1318   // WARNING: THIS METHOD HAS TO BE SYNCHRONIZED WITH AliITStrackV2::GetModuleIndexInfo()!
1319   //----------------------------------------------------------------------
1320
1321   if(fITSModule[ilayer]==-1) {
1322     idet = -1;
1323     status=0;
1324     xloc=-99.; zloc=-99.;
1325     return kFALSE;
1326   }
1327
1328   Int_t module = fITSModule[ilayer];
1329
1330   idet = Int_t(module/1000000);
1331
1332   module -= idet*1000000;
1333
1334   status = Int_t(module/100000);
1335
1336   module -= status*100000;
1337
1338   Int_t signs = Int_t(module/10000);
1339
1340   module-=signs*10000;
1341
1342   Int_t xInt = Int_t(module/100);
1343   module -= xInt*100;
1344
1345   Int_t zInt = module;
1346
1347   if(signs==1) { xInt*=1; zInt*=1; }
1348   if(signs==2) { xInt*=1; zInt*=-1; }
1349   if(signs==3) { xInt*=-1; zInt*=1; }
1350   if(signs==4) { xInt*=-1; zInt*=-1; }
1351
1352   xloc = 0.1*(Float_t)xInt;
1353   zloc = 0.1*(Float_t)zInt;
1354
1355   if(status==4) idet = -1;
1356
1357   return kTRUE;
1358 }
1359
1360 //_______________________________________________________________________
1361 UShort_t AliESDtrack::GetTPCclusters(Int_t *idx) const {
1362   //---------------------------------------------------------------------
1363   // This function returns indices of the assgined ITS clusters
1364   //---------------------------------------------------------------------
1365   if (idx!=0) {
1366     Int_t *index=fFriendTrack->GetTPCindices();
1367     for (Int_t i=0; i<AliESDfriendTrack::kMaxTPCcluster; i++) idx[i]=index[i];
1368   }
1369   return fTPCncls;
1370 }
1371
1372 Double_t AliESDtrack::GetTPCdensity(Int_t row0, Int_t row1) const{
1373   //
1374   // GetDensity of the clusters on given region between row0 and row1
1375   // Dead zone effect takin into acoount
1376   //
1377   Int_t good  = 0;
1378   Int_t found = 0;
1379   //
1380   Int_t *index=fFriendTrack->GetTPCindices();
1381   for (Int_t i=row0;i<=row1;i++){
1382     Int_t idx = index[i];
1383     if (idx!=-1)  good++;             // track outside of dead zone
1384     if (idx>0)    found++;
1385   }
1386   Float_t density=0.5;
1387   if (good>(row1-row0)*0.5) density = Float_t(found)/Float_t(good);
1388   return density;
1389 }
1390
1391 //_______________________________________________________________________
1392 void AliESDtrack::SetTPCpid(const Double_t *p) {
1393   // Sets values for the probability of each particle type (in TPC)
1394   SetPIDValues(fTPCr,p,AliPID::kSPECIES);
1395   SetStatus(AliESDtrack::kTPCpid);
1396 }
1397
1398 //_______________________________________________________________________
1399 void AliESDtrack::GetTPCpid(Double_t *p) const {
1400   // Gets the probability of each particle type (in TPC)
1401   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fTPCr[i];
1402 }
1403
1404 //_______________________________________________________________________
1405 UChar_t AliESDtrack::GetTRDclusters(Int_t *idx) const {
1406   //---------------------------------------------------------------------
1407   // This function returns indices of the assgined TRD clusters
1408   //---------------------------------------------------------------------
1409   if (idx!=0) {
1410      Int_t *index=fFriendTrack->GetTRDindices();
1411      for (Int_t i=0; i<AliESDfriendTrack::kMaxTRDcluster; i++) idx[i]=index[i];
1412   }
1413   return fTRDncls;
1414 }
1415
1416 //_______________________________________________________________________
1417 UChar_t AliESDtrack::GetTRDtracklets(Int_t *idx) const {
1418   //---------------------------------------------------------------------
1419   // This function returns indices of the assigned TRD tracklets
1420   //---------------------------------------------------------------------
1421   if (idx!=0) {
1422      Int_t *index=fFriendTrack->GetTRDindices();
1423      for (Int_t i=0; i<6/*AliESDfriendTrack::kMaxTRDcluster*/; i++) idx[i]=index[i];
1424   }
1425   return fTRDncls;
1426 }
1427
1428 //_______________________________________________________________________
1429 void AliESDtrack::SetTRDpid(const Double_t *p) {
1430   // Sets values for the probability of each particle type (in TRD)
1431   SetPIDValues(fTRDr,p,AliPID::kSPECIES);
1432   SetStatus(AliESDtrack::kTRDpid);
1433 }
1434
1435 //_______________________________________________________________________
1436 void AliESDtrack::GetTRDpid(Double_t *p) const {
1437   // Gets the probability of each particle type (in TRD)
1438   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fTRDr[i];
1439 }
1440
1441 //_______________________________________________________________________
1442 void    AliESDtrack::SetTRDpid(Int_t iSpecies, Float_t p)
1443 {
1444   // Sets the probability of particle type iSpecies to p (in TRD)
1445   fTRDr[iSpecies] = p;
1446 }
1447
1448 Double_t AliESDtrack::GetTRDpid(Int_t iSpecies) const
1449 {
1450   // Returns the probability of particle type iSpecies (in TRD)
1451   return fTRDr[iSpecies];
1452 }
1453
1454 void  AliESDtrack::SetNumberOfTRDslices(Int_t n) {
1455   //Sets the number of slices used for PID
1456   if (fTRDnSlices != 0) return;
1457   fTRDnSlices=kTRDnPlanes*n;
1458   fTRDslices=new Double32_t[fTRDnSlices];
1459   for (Int_t i=0; i<fTRDnSlices; i++) fTRDslices[i]=-1.;
1460 }
1461
1462 void  AliESDtrack::SetTRDslice(Double_t q, Int_t plane, Int_t slice) {
1463   //Sets the charge q in the slice of the plane
1464   Int_t ns=GetNumberOfTRDslices();
1465   if (ns==0) {
1466     AliError("No TRD slices allocated for this track !");
1467     return;
1468   }
1469
1470   if ((plane<0) || (plane>=kTRDnPlanes)) {
1471     AliError("Wrong TRD plane !");
1472     return;
1473   }
1474   if ((slice<0) || (slice>=ns)) {
1475     AliError("Wrong TRD slice !");
1476     return;
1477   }
1478   Int_t n=plane*ns + slice;
1479   fTRDslices[n]=q;
1480 }
1481
1482 Double_t  AliESDtrack::GetTRDslice(Int_t plane, Int_t slice) const {
1483   //Gets the charge from the slice of the plane
1484   Int_t ns=GetNumberOfTRDslices();
1485   if (ns==0) {
1486     //AliError("No TRD slices allocated for this track !");
1487     return -1.;
1488   }
1489
1490   if ((plane<0) || (plane>=kTRDnPlanes)) {
1491     AliError("Wrong TRD plane !");
1492     return -1.;
1493   }
1494   if ((slice<-1) || (slice>=ns)) {
1495     //AliError("Wrong TRD slice !");
1496     return -1.;
1497   }
1498
1499   if (slice==-1) {
1500     Double_t q=0.;
1501     for (Int_t i=0; i<ns; i++) q+=fTRDslices[plane*ns + i];
1502     return q/ns;
1503   }
1504
1505   return fTRDslices[plane*ns + slice];
1506 }
1507
1508
1509 //_______________________________________________________________________
1510 void AliESDtrack::SetTOFpid(const Double_t *p) {
1511   // Sets the probability of each particle type (in TOF)
1512   SetPIDValues(fTOFr,p,AliPID::kSPECIES);
1513   SetStatus(AliESDtrack::kTOFpid);
1514 }
1515
1516 //_______________________________________________________________________
1517 void AliESDtrack::SetTOFLabel(const Int_t *p) {
1518   // Sets  (in TOF)
1519   for (Int_t i=0; i<3; i++) fTOFLabel[i]=p[i];
1520 }
1521
1522 //_______________________________________________________________________
1523 void AliESDtrack::GetTOFpid(Double_t *p) const {
1524   // Gets probabilities of each particle type (in TOF)
1525   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fTOFr[i];
1526 }
1527
1528 //_______________________________________________________________________
1529 void AliESDtrack::GetTOFLabel(Int_t *p) const {
1530   // Gets (in TOF)
1531   for (Int_t i=0; i<3; i++) p[i]=fTOFLabel[i];
1532 }
1533
1534 //_______________________________________________________________________
1535 void AliESDtrack::GetTOFInfo(Float_t *info) const {
1536   // Gets (in TOF)
1537   for (Int_t i=0; i<10; i++) info[i]=fTOFInfo[i];
1538 }
1539
1540 //_______________________________________________________________________
1541 void AliESDtrack::SetTOFInfo(Float_t*info) {
1542   // Gets (in TOF)
1543   for (Int_t i=0; i<10; i++) fTOFInfo[i]=info[i];
1544 }
1545
1546
1547
1548 //_______________________________________________________________________
1549 void AliESDtrack::SetHMPIDpid(const Double_t *p) {
1550   // Sets the probability of each particle type (in HMPID)
1551   SetPIDValues(fHMPIDr,p,AliPID::kSPECIES);
1552   SetStatus(AliESDtrack::kHMPIDpid);
1553 }
1554
1555 //_______________________________________________________________________
1556 void AliESDtrack::GetHMPIDpid(Double_t *p) const {
1557   // Gets probabilities of each particle type (in HMPID)
1558   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fHMPIDr[i];
1559 }
1560
1561
1562
1563 //_______________________________________________________________________
1564 void AliESDtrack::SetESDpid(const Double_t *p) {
1565   // Sets the probability of each particle type for the ESD track
1566   SetPIDValues(fR,p,AliPID::kSPECIES);
1567   SetStatus(AliESDtrack::kESDpid);
1568 }
1569
1570 //_______________________________________________________________________
1571 void AliESDtrack::GetESDpid(Double_t *p) const {
1572   // Gets probability of each particle type for the ESD track
1573   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fR[i];
1574 }
1575
1576 //_______________________________________________________________________
1577 Bool_t AliESDtrack::RelateToVertexTPC(const AliESDVertex *vtx,
1578 Double_t b, Double_t maxd, AliExternalTrackParam *cParam) {
1579   //
1580   // Try to relate the TPC-only track parameters to the vertex "vtx",
1581   // if the (rough) transverse impact parameter is not bigger then "maxd".
1582   //            Magnetic field is "b" (kG).
1583   //
1584   // a) The TPC-only paramters are extapolated to the DCA to the vertex.
1585   // b) The impact parameters and their covariance matrix are calculated.
1586   // c) An attempt to constrain the TPC-only params to the vertex is done.
1587   //    The constrained params are returned via "cParam".
1588   //
1589   // In the case of success, the returned value is kTRUE
1590   // otherwise, it's kFALSE)
1591   //
1592
1593   if (!fTPCInner) return kFALSE;
1594   if (!vtx) return kFALSE;
1595
1596   Double_t dz[2],cov[3];
1597   if (!fTPCInner->PropagateToDCA(vtx, b, maxd, dz, cov)) return kFALSE;
1598
1599   fdTPC = dz[0];
1600   fzTPC = dz[1];
1601   fCddTPC = cov[0];
1602   fCdzTPC = cov[1];
1603   fCzzTPC = cov[2];
1604
1605   Double_t covar[6]; vtx->GetCovMatrix(covar);
1606   Double_t p[2]={GetParameter()[0]-dz[0],GetParameter()[1]-dz[1]};
1607   Double_t c[3]={covar[2],0.,covar[5]};
1608
1609   Double_t chi2=GetPredictedChi2(p,c);
1610   if (chi2>kVeryBig) return kFALSE;
1611
1612   fCchi2TPC=chi2;
1613
1614   if (!cParam) return kTRUE;
1615
1616   *cParam = *fTPCInner;
1617   if (!cParam->Update(p,c)) return kFALSE;
1618
1619   return kTRUE;
1620 }
1621
1622 //_______________________________________________________________________
1623 Bool_t AliESDtrack::RelateToVertex(const AliESDVertex *vtx,
1624 Double_t b, Double_t maxd, AliExternalTrackParam *cParam) {
1625   //
1626   // Try to relate this track to the vertex "vtx",
1627   // if the (rough) transverse impact parameter is not bigger then "maxd".
1628   //            Magnetic field is "b" (kG).
1629   //
1630   // a) The track gets extapolated to the DCA to the vertex.
1631   // b) The impact parameters and their covariance matrix are calculated.
1632   // c) An attempt to constrain this track to the vertex is done.
1633   //    The constrained params are returned via "cParam".
1634   //
1635   // In the case of success, the returned value is kTRUE
1636   // (otherwise, it's kFALSE)
1637   //
1638
1639   if (!vtx) return kFALSE;
1640
1641   Double_t dz[2],cov[3];
1642   if (!PropagateToDCA(vtx, b, maxd, dz, cov)) return kFALSE;
1643
1644   fD = dz[0];
1645   fZ = dz[1];
1646   fCdd = cov[0];
1647   fCdz = cov[1];
1648   fCzz = cov[2];
1649
1650   Double_t covar[6]; vtx->GetCovMatrix(covar);
1651   Double_t p[2]={GetParameter()[0]-dz[0],GetParameter()[1]-dz[1]};
1652   Double_t c[3]={covar[2],0.,covar[5]};
1653
1654   Double_t chi2=GetPredictedChi2(p,c);
1655   if (chi2>kVeryBig) return kFALSE;
1656
1657   fCchi2=chi2;
1658
1659
1660   //--- Could now these lines be removed ? ---
1661   delete fCp;
1662   fCp=new AliExternalTrackParam(*this);
1663
1664   if (!fCp->Update(p,c)) {delete fCp; fCp=0; return kFALSE;}
1665   //----------------------------------------
1666
1667
1668   if (!cParam) return kTRUE;
1669
1670   *cParam = *this;
1671   if (!cParam->Update(p,c)) return kFALSE;
1672
1673   return kTRUE;
1674 }
1675
1676 //_______________________________________________________________________
1677 void AliESDtrack::Print(Option_t *) const {
1678   // Prints info on the track
1679   AliExternalTrackParam::Print();
1680   printf("ESD track info\n") ;
1681   Double_t p[AliPID::kSPECIESN] ;
1682   Int_t index = 0 ;
1683   if( IsOn(kITSpid) ){
1684     printf("From ITS: ") ;
1685     GetITSpid(p) ;
1686     for(index = 0 ; index < AliPID::kSPECIES; index++)
1687       printf("%f, ", p[index]) ;
1688     printf("\n           signal = %f\n", GetITSsignal()) ;
1689   }
1690   if( IsOn(kTPCpid) ){
1691     printf("From TPC: ") ;
1692     GetTPCpid(p) ;
1693     for(index = 0 ; index < AliPID::kSPECIES; index++)
1694       printf("%f, ", p[index]) ;
1695     printf("\n           signal = %f\n", GetTPCsignal()) ;
1696   }
1697   if( IsOn(kTRDpid) ){
1698     printf("From TRD: ") ;
1699     GetTRDpid(p) ;
1700     for(index = 0 ; index < AliPID::kSPECIES; index++)
1701       printf("%f, ", p[index]) ;
1702       printf("\n           signal = %f\n", GetTRDsignal()) ;
1703   }
1704   if( IsOn(kTOFpid) ){
1705     printf("From TOF: ") ;
1706     GetTOFpid(p) ;
1707     for(index = 0 ; index < AliPID::kSPECIES; index++)
1708       printf("%f, ", p[index]) ;
1709     printf("\n           signal = %f\n", GetTOFsignal()) ;
1710   }
1711   if( IsOn(kHMPIDpid) ){
1712     printf("From HMPID: ") ;
1713     GetHMPIDpid(p) ;
1714     for(index = 0 ; index < AliPID::kSPECIES; index++)
1715       printf("%f, ", p[index]) ;
1716     printf("\n           signal = %f\n", GetHMPIDsignal()) ;
1717   }
1718 }
1719
1720
1721 //
1722 // Draw functionality
1723 // Origin: Marian Ivanov, Marian.Ivanov@cern.ch
1724 //
1725 void AliESDtrack::FillPolymarker(TPolyMarker3D *pol, Float_t magF, Float_t minR, Float_t maxR, Float_t stepR){
1726   //
1727   // Fill points in the polymarker
1728   //
1729   TObjArray arrayRef;
1730   arrayRef.AddLast(new AliExternalTrackParam(*this));
1731   if (fIp) arrayRef.AddLast(new AliExternalTrackParam(*fIp));
1732   if (fOp) arrayRef.AddLast(new AliExternalTrackParam(*fOp));
1733   //
1734   Double_t mpos[3]={0,0,0};
1735   Int_t entries=arrayRef.GetEntries();
1736   for (Int_t i=0;i<entries;i++){
1737     Double_t pos[3];
1738     ((AliExternalTrackParam*)arrayRef.At(i))->GetXYZ(pos);
1739     mpos[0]+=pos[0]/entries;
1740     mpos[1]+=pos[1]/entries;
1741     mpos[2]+=pos[2]/entries;
1742   }
1743   // Rotate to the mean position
1744   //
1745   Float_t fi= TMath::ATan2(mpos[1],mpos[0]);
1746   for (Int_t i=0;i<entries;i++){
1747     Bool_t res = ((AliExternalTrackParam*)arrayRef.At(i))->Rotate(fi);
1748     if (!res) delete arrayRef.RemoveAt(i);
1749   }
1750   Int_t counter=0;
1751   for (Double_t r=minR; r<maxR; r+=stepR){
1752     Double_t sweight=0;
1753     Double_t mlpos[3]={0,0,0};
1754     for (Int_t i=0;i<entries;i++){
1755       Double_t point[3]={0,0,0};
1756       AliExternalTrackParam *param = ((AliExternalTrackParam*)arrayRef.At(i));
1757       if (!param) continue;
1758       if (param->GetXYZAt(r,magF,point)){
1759         Double_t weight = 1./(10.+(r-param->GetX())*(r-param->GetX()));
1760         sweight+=weight;
1761         mlpos[0]+=point[0]*weight;
1762         mlpos[1]+=point[1]*weight;
1763         mlpos[2]+=point[2]*weight;
1764       }
1765     }
1766     if (sweight>0){
1767       mlpos[0]/=sweight;
1768       mlpos[1]/=sweight;
1769       mlpos[2]/=sweight;
1770       pol->SetPoint(counter,mlpos[0],mlpos[1], mlpos[2]);
1771       printf("xyz\t%f\t%f\t%f\n",mlpos[0], mlpos[1],mlpos[2]);
1772       counter++;
1773     }
1774   }
1775 }