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 #include <TMath.h>
  24 #include <TParticle.h>
  25
  26 #include "AliESDVertex.h"
  27 #include "AliESDtrack.h"
  28 #include "AliKalmanTrack.h"
  29 #include "AliLog.h"
  30 #include "AliTrackPointArray.h"
  31
  32 ClassImp(AliESDtrack)
  33
  34 void SetPIDValues(Double_t * dest, const Double_t * src, Int_t n) {
  35   // This function copies "n" PID weights from "scr" to "dest"
  36   // and normalizes their sum to 1 thus producing conditional probabilities.
  37   // The negative weights are set to 0.
  38   // In case all the weights are non-positive they are replaced by
  39   // uniform probabilities
  40
  41   if (n<=0) return;
  42
  43   Float_t uniform = 1./(Float_t)n;
  44
  45   Float_t sum = 0;
  46   for (Int_t i=0; i<n; i++)
  47     if (src[i]>=0) {
  48       sum+=src[i];
  49       dest[i] = src[i];
  50     }
  51     else {
  52       dest[i] = 0;
  53     }
  54
  55   if(sum>0)
  56     for (Int_t i=0; i<n; i++) dest[i] /= sum;
  57   else
  58     for (Int_t i=0; i<n; i++) dest[i] = uniform;
  59 }
  60
  61 //_______________________________________________________________________
  62 AliESDtrack::AliESDtrack() :
  63   AliExternalTrackParam(),
  64   fCp(0),
  65   fIp(0),
  66   fTPCInner(0),
  67   fOp(0),
  68   fFriendTrack(new AliESDfriendTrack()),
  69   fTPCClusterMap(159),//number of padrows
  70   fTPCSharedMap(159),//number of padrows
  71   fFlags(0),
  72   fID(0),
  73   fLabel(0),
  74   fITSLabel(0),
  75   fTPCLabel(0),
  76   fTRDLabel(0),
  77   fTOFCalChannel(0),
  78   fTOFindex(0),
  79   fHMPIDqn(0),
  80   fHMPIDcluIdx(0),
  81   fEMCALindex(kEMCALNoMatch),
  82   fHMPIDtrkTheta(0),
  83   fHMPIDtrkPhi(0),
  84   fHMPIDsignal(0),
  85   fTrackLength(0),
  86   fD(0),fZ(0),
  87   fCdd(0),fCdz(0),fCzz(0),
  88   fCchi2(0),
  89   fITSchi2(0),
  90   fTPCchi2(0),
  91   fTRDchi2(0),
  92   fTOFchi2(0),
  93   fHMPIDchi2(0),
  94   fITSsignal(0),
  95   fTPCsignal(0),
  96   fTPCsignalS(0),
  97   fTRDsignal(0),
  98   fTRDQuality(0),
  99   fTRDBudget(0),
 100   fTOFsignal(0),
 101   fTOFsignalToT(0),
 102   fTOFsignalRaw(0),
 103   fTOFsignalDz(0),
 104   fHMPIDtrkX(0),
 105   fHMPIDtrkY(0),
 106   fHMPIDmipX(0),
 107   fHMPIDmipY(0),
 108   fTPCncls(0),
 109   fTPCnclsF(0),
 110   fTPCsignalN(0),
 111   fITSncls(0),
 112   fITSClusterMap(0),
 113   fTRDncls(0),
 114   fTRDncls0(0),
 115   fTRDpidQuality(0)
 116 {
 117   //
 118   // The default ESD constructor
 119   //
 120   Int_t i, j;
 121   for (i=0; i<AliPID::kSPECIES; i++) {
 122     fTrackTime[i]=0.;
 123     fR[i]=0.;
 124     fITSr[i]=0.;
 125     fTPCr[i]=0.;
 126     fTRDr[i]=0.;
 127     fTOFr[i]=0.;
 128     fHMPIDr[i]=0.;
 129   }
 130
 131   for (i=0; i<3; i++)   { fKinkIndexes[i]=0;}
 132   for (i=0; i<3; i++)   { fV0Indexes[i]=0;}
 133   for (i=0;i<kNPlane;i++) {
 134     for (j=0;j<kNSlice;j++) {
 135       fTRDsignals[i][j]=0.;
 136     }
 137     fTRDTimBin[i]=0;
 138   }
 139   for (i=0;i<4;i++) {fTPCPoints[i]=0;}
 140   for (i=0;i<3;i++) {fTOFLabel[i]=0;}
 141   for (i=0;i<10;i++) {fTOFInfo[i]=0;}
 142   for (i=0;i<12;i++) {fITSModule[i]=-1;}
 143 }
 144
 145 //_______________________________________________________________________
 146 AliESDtrack::AliESDtrack(const AliESDtrack& track):
 147   AliExternalTrackParam(track),
 148   fCp(0),
 149   fIp(0),
 150   fTPCInner(0),
 151   fOp(0),
 152   fFriendTrack(0),
 153   fTPCClusterMap(track.fTPCClusterMap),
 154   fTPCSharedMap(track.fTPCSharedMap),
 155   fFlags(track.fFlags),
 156   fID(track.fID),
 157   fLabel(track.fLabel),
 158   fITSLabel(track.fITSLabel),
 159   fTPCLabel(track.fTPCLabel),
 160   fTRDLabel(track.fTRDLabel),
 161   fTOFCalChannel(track.fTOFCalChannel),
 162   fTOFindex(track.fTOFindex),
 163   fHMPIDqn(track.fHMPIDqn),
 164   fHMPIDcluIdx(track.fHMPIDcluIdx),
 165   fEMCALindex(track.fEMCALindex),
 166   fHMPIDtrkTheta(track.fHMPIDtrkTheta),
 167   fHMPIDtrkPhi(track.fHMPIDtrkPhi),
 168   fHMPIDsignal(track.fHMPIDsignal),
 169   fTrackLength(track.fTrackLength),
 170   fD(track.fD),fZ(track.fZ),
 171   fCdd(track.fCdd),fCdz(track.fCdz),fCzz(track.fCzz),
 172   fCchi2(track.fCchi2),
 173   fITSchi2(track.fITSchi2),
 174   fTPCchi2(track.fTPCchi2),
 175   fTRDchi2(track.fTRDchi2),
 176   fTOFchi2(track.fTOFchi2),
 177   fHMPIDchi2(track.fHMPIDchi2),
 178   fITSsignal(track.fITSsignal),
 179   fTPCsignal(track.fTPCsignal),
 180   fTPCsignalS(track.fTPCsignalS),
 181   fTRDsignal(track.fTRDsignal),
 182   fTRDQuality(track.fTRDQuality),
 183   fTRDBudget(track.fTRDBudget),
 184   fTOFsignal(track.fTOFsignal),
 185   fTOFsignalToT(track.fTOFsignalToT),
 186   fTOFsignalRaw(track.fTOFsignalRaw),
 187   fTOFsignalDz(track.fTOFsignalDz),
 188   fHMPIDtrkX(track.fHMPIDtrkX),
 189   fHMPIDtrkY(track.fHMPIDtrkY),
 190   fHMPIDmipX(track.fHMPIDmipX),
 191   fHMPIDmipY(track.fHMPIDmipY),
 192   fTPCncls(track.fTPCncls),
 193   fTPCnclsF(track.fTPCnclsF),
 194   fTPCsignalN(track.fTPCsignalN),
 195   fITSncls(track.fITSncls),
 196   fITSClusterMap(track.fITSClusterMap),
 197   fTRDncls(track.fTRDncls),
 198   fTRDncls0(track.fTRDncls0),
 199   fTRDpidQuality(track.fTRDpidQuality)
 200 {
 201   //
 202   //copy constructor
 203   //
 204   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTrackTime[i]=track.fTrackTime[i];
 205   for (Int_t i=0;i<AliPID::kSPECIES;i++)  fR[i]=track.fR[i];
 206   //
 207   for (Int_t i=0;i<AliPID::kSPECIES;i++) fITSr[i]=track.fITSr[i];
 208   //
 209   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTPCr[i]=track.fTPCr[i];
 210   for (Int_t i=0;i<4;i++) {fTPCPoints[i]=track.fTPCPoints[i];}
 211   for (Int_t i=0; i<3;i++)   { fKinkIndexes[i]=track.fKinkIndexes[i];}
 212   for (Int_t i=0; i<3;i++)   { fV0Indexes[i]=track.fV0Indexes[i];}
 213   //
 214   for (Int_t i=0;i<kNPlane;i++) {
 215     for (Int_t j=0;j<kNSlice;j++) {
 216       fTRDsignals[i][j]=track.fTRDsignals[i][j];
 217     }
 218     fTRDTimBin[i]=track.fTRDTimBin[i];
 219   }
 220   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTRDr[i]=track.fTRDr[i];
 221   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTOFr[i]=track.fTOFr[i];
 222   for (Int_t i=0;i<3;i++) fTOFLabel[i]=track.fTOFLabel[i];
 223   for (Int_t i=0;i<10;i++) fTOFInfo[i]=track.fTOFInfo[i];
 224   for (Int_t i=0;i<12;i++) fITSModule[i]=track.fITSModule[i];
 225   for (Int_t i=0;i<AliPID::kSPECIES;i++) fHMPIDr[i]=track.fHMPIDr[i];
 226
 227   if (track.fCp) fCp=new AliExternalTrackParam(*track.fCp);
 228   if (track.fIp) fIp=new AliExternalTrackParam(*track.fIp);
 229   if (track.fTPCInner) fTPCInner=new AliExternalTrackParam(*track.fTPCInner);
 230   if (track.fOp) fOp=new AliExternalTrackParam(*track.fOp);
 231
 232   if (track.fFriendTrack) fFriendTrack=new AliESDfriendTrack(*(track.fFriendTrack));
 233 }
 234
 235 //_______________________________________________________________________
 236 AliESDtrack::AliESDtrack(TParticle * part) :
 237   AliExternalTrackParam(),
 238   fCp(0),
 239   fIp(0),
 240   fTPCInner(0),
 241   fOp(0),
 242   fFriendTrack(0),
 243   fTPCClusterMap(159),//number of padrows
 244   fTPCSharedMap(159),//number of padrows
 245   fFlags(0),
 246   fID(0),
 247   fLabel(0),
 248   fITSLabel(0),
 249   fTPCLabel(0),
 250   fTRDLabel(0),
 251   fTOFCalChannel(0),
 252   fTOFindex(0),
 253   fHMPIDqn(0),
 254   fHMPIDcluIdx(0),
 255   fEMCALindex(kEMCALNoMatch),
 256   fHMPIDtrkTheta(0),
 257   fHMPIDtrkPhi(0),
 258   fHMPIDsignal(0),
 259   fTrackLength(0),
 260   fD(0),fZ(0),
 261   fCdd(0),fCdz(0),fCzz(0),
 262   fCchi2(0),
 263   fITSchi2(0),
 264   fTPCchi2(0),
 265   fTRDchi2(0),
 266   fTOFchi2(0),
 267   fHMPIDchi2(0),
 268   fITSsignal(0),
 269   fTPCsignal(0),
 270   fTPCsignalS(0),
 271   fTRDsignal(0),
 272   fTRDQuality(0),
 273   fTRDBudget(0),
 274   fTOFsignal(0),
 275   fTOFsignalToT(0),
 276   fTOFsignalRaw(0),
 277   fTOFsignalDz(0),
 278   fHMPIDtrkX(0),
 279   fHMPIDtrkY(0),
 280   fHMPIDmipX(0),
 281   fHMPIDmipY(0),
 282   fTPCncls(0),
 283   fTPCnclsF(0),
 284   fTPCsignalN(0),
 285   fITSncls(0),
 286   fITSClusterMap(0),
 287   fTRDncls(0),
 288   fTRDncls0(0),
 289   fTRDpidQuality(0)
 290 {
 291   //
 292   // ESD track from TParticle
 293   //
 294
 295   // Reset all the arrays
 296   Int_t i, j;
 297   for (i=0; i<AliPID::kSPECIES; i++) {
 298     fTrackTime[i]=0.;
 299     fR[i]=0.;
 300     fITSr[i]=0.;
 301     fTPCr[i]=0.;
 302     fTRDr[i]=0.;
 303     fTOFr[i]=0.;
 304     fHMPIDr[i]=0.;
 305   }
 306
 307   for (i=0; i<3; i++)   { fKinkIndexes[i]=0;}
 308   for (i=0; i<3; i++)   { fV0Indexes[i]=-1;}
 309   for (i=0;i<kNPlane;i++) {
 310     for (j=0;j<kNSlice;j++) {
 311       fTRDsignals[i][j]=0.;
 312     }
 313     fTRDTimBin[i]=0;
 314   }
 315   for (i=0;i<4;i++) {fTPCPoints[i]=0;}
 316   for (i=0;i<3;i++) {fTOFLabel[i]=0;}
 317   for (i=0;i<10;i++) {fTOFInfo[i]=0;}
 318   for (i=0;i<12;i++) {fITSModule[i]=-1;}
 319
 320   // Calculate the AliExternalTrackParam content
 321
 322   Double_t xref;
 323   Double_t alpha;
 324   Double_t param[5];
 325   Double_t covar[15];
 326
 327   // Calculate alpha: the rotation angle of the corresponding local system (TPC sector)
 328   alpha = part->Phi()*180./TMath::Pi();
 329   if (alpha<0) alpha+= 360.;
 330   if (alpha>360) alpha -= 360.;
 331
 332   Int_t sector = (Int_t)(alpha/20.);
 333   alpha = 10. + 20.*sector;
 334   alpha /= 180;
 335   alpha *= TMath::Pi();
 336
 337   // Covariance matrix: no errors, the parameters are exact
 338   for (i=0; i<15; i++) covar[i]=0.;
 339
 340   // Get the vertex of origin and the momentum
 341   TVector3 ver(part->Vx(),part->Vy(),part->Vz());
 342   TVector3 mom(part->Px(),part->Py(),part->Pz());
 343
 344   // Rotate to the local coordinate system (TPC sector)
 345   ver.RotateZ(-alpha);
 346   mom.RotateZ(-alpha);
 347
 348   // X of the referense plane
 349   xref = ver.X();
 350
 351   Int_t pdgCode = part->GetPdgCode();
 352
 353   Double_t charge =
 354     TDatabasePDG::Instance()->GetParticle(pdgCode)->Charge();
 355
 356   param[0] = ver.Y();
 357   param[1] = ver.Z();
 358   param[2] = TMath::Sin(mom.Phi());
 359   param[3] = mom.Pz()/mom.Pt();
 360   param[4] = TMath::Sign(1/mom.Pt(),charge);
 361
 362   // Set AliExternalTrackParam
 363   Set(xref, alpha, param, covar);
 364
 365   // Set the PID
 366   Int_t indexPID = 99;
 367
 368   switch (TMath::Abs(pdgCode)) {
 369
 370   case  11: // electron
 371     indexPID = 0;
 372     break;
 373
 374   case 13: // muon
 375     indexPID = 1;
 376     break;
 377
 378   case 211: // pion
 379     indexPID = 2;
 380     break;
 381
 382   case 321: // kaon
 383     indexPID = 3;
 384     break;
 385
 386   case 2212: // proton
 387     indexPID = 4;
 388     break;
 389
 390   default:
 391     break;
 392   }
 393
 394   // If the particle is not e,mu,pi,K or p the PID probabilities are set to 0
 395   if (indexPID < AliPID::kSPECIES) {
 396     fR[indexPID]=1.;
 397     fITSr[indexPID]=1.;
 398     fTPCr[indexPID]=1.;
 399     fTRDr[indexPID]=1.;
 400     fTOFr[indexPID]=1.;
 401     fHMPIDr[indexPID]=1.;
 402
 403   }
 404   // AliESD track label
 405   SetLabel(part->GetUniqueID());
 406
 407 }
 408
 409 //_______________________________________________________________________
 410 AliESDtrack::~AliESDtrack(){
 411   //
 412   // This is destructor according Coding Conventrions
 413   //
 414   //printf("Delete track\n");
 415   delete fIp;
 416   delete fTPCInner;
 417   delete fOp;
 418   delete fCp;
 419   delete fFriendTrack;
 420 }
 421
 422 void AliESDtrack::AddCalibObject(TObject * object){
 423   //
 424   // add calib object to the list
 425   //
 426   if (!fFriendTrack) fFriendTrack  = new AliESDfriendTrack;
 427   fFriendTrack->AddCalibObject(object);
 428 }
 429
 430 TObject *  AliESDtrack::GetCalibObject(Int_t index){
 431   //
 432   // return calib objct at given position
 433   //
 434   if (!fFriendTrack) return 0;
 435   return fFriendTrack->GetCalibObject(index);
 436 }
 437
 438
 439 //_______________________________________________________________________
 440 void AliESDtrack::MakeMiniESDtrack(){
 441   // Resets everything except
 442   // fFlags: Reconstruction status flags
 443   // fLabel: Track label
 444   // fID:  Unique ID of the track
 445   // fD: Impact parameter in XY-plane
 446   // fZ: Impact parameter in Z
 447   // fR[AliPID::kSPECIES]: combined "detector response probability"
 448   // Running track parameters in the base class (AliExternalTrackParam)
 449
 450   fTrackLength = 0;
 451
 452   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTrackTime[i] = 0;
 453
 454   // Reset track parameters constrained to the primary vertex
 455   delete fCp;fCp = 0;
 456   fCchi2 = 0;
 457
 458   // Reset track parameters at the inner wall of TPC
 459   delete fIp;fIp = 0;
 460   delete fTPCInner;fTPCInner=0;
 461   // Reset track parameters at the inner wall of the TRD
 462   delete fOp;fOp = 0;
 463
 464
 465   // Reset ITS track related information
 466   fITSchi2 = 0;
 467   fITSncls = 0;
 468   fITSClusterMap=0;
 469   fITSsignal = 0;
 470   for (Int_t i=0;i<AliPID::kSPECIES;i++) fITSr[i]=0;
 471   fITSLabel = 0;
 472
 473   // Reset TPC related track information
 474   fTPCchi2 = 0;
 475   fTPCncls = 0;
 476   fTPCnclsF = 0;
 477   fTPCClusterMap = 0;
 478   fTPCSharedMap = 0;
 479   fTPCsignal= 0;
 480   fTPCsignalS= 0;
 481   fTPCsignalN= 0;
 482   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTPCr[i]=0;
 483   fTPCLabel=0;
 484   for (Int_t i=0;i<4;i++) fTPCPoints[i] = 0;
 485   for (Int_t i=0; i<3;i++)   fKinkIndexes[i] = 0;
 486   for (Int_t i=0; i<3;i++)   fV0Indexes[i] = 0;
 487
 488   // Reset TRD related track information
 489   fTRDchi2 = 0;
 490   fTRDncls = 0;
 491   fTRDncls0 = 0;
 492   fTRDsignal = 0;
 493   for (Int_t i=0;i<kNPlane;i++) {
 494     for (Int_t j=0;j<kNSlice;j++) {
 495       fTRDsignals[i][j] = 0;
 496     }
 497     fTRDTimBin[i]  = 0;
 498   }
 499   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTRDr[i] = 0;
 500   fTRDLabel = 0;
 501   fTRDQuality  = 0;
 502   fTRDpidQuality = 0;
 503   fTRDBudget  = 0;
 504
 505   // Reset TOF related track information
 506   fTOFchi2 = 0;
 507   fTOFindex = 0;
 508   fTOFsignal = 0;
 509   fTOFCalChannel = 0;
 510   fTOFsignalToT = 0;
 511   fTOFsignalRaw = 0;
 512   fTOFsignalDz = 0;
 513   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTOFr[i] = 0;
 514   for (Int_t i=0;i<3;i++) fTOFLabel[i] = 0;
 515   for (Int_t i=0;i<10;i++) fTOFInfo[i] = 0;
 516
 517   // Reset HMPID related track information
 518   fHMPIDchi2 = 0;
 519   fHMPIDqn = 0;
 520   fHMPIDcluIdx = 0;
 521   fHMPIDsignal = 0;
 522   for (Int_t i=0;i<AliPID::kSPECIES;i++) fHMPIDr[i] = 0;
 523   fHMPIDtrkTheta = 0;
 524   fHMPIDtrkPhi = 0;
 525   fHMPIDtrkX = 0;
 526   fHMPIDtrkY = 0;
 527   fHMPIDmipX = 0;
 528   fHMPIDmipY = 0;
 529   fEMCALindex = kEMCALNoMatch;
 530
 531   delete fFriendTrack; fFriendTrack = 0;
 532 }
 533 //_______________________________________________________________________
 534 Double_t AliESDtrack::GetMass() const {
 535   // Returns the mass of the most probable particle type
 536   Float_t max=0.;
 537   Int_t k=-1;
 538   for (Int_t i=0; i<AliPID::kSPECIES; i++) {
 539     if (fR[i]>max) {k=i; max=fR[i];}
 540   }
 541   if (k==0) { // dE/dx "crossing points" in the TPC
 542      Double_t p=GetP();
 543      if ((p>0.38)&&(p<0.48))
 544         if (fR[0]<fR[3]*10.) return AliPID::ParticleMass(AliPID::kKaon);
 545      if ((p>0.75)&&(p<0.85))
 546         if (fR[0]<fR[4]*10.) return AliPID::ParticleMass(AliPID::kProton);
 547      return 0.00051;
 548   }
 549   if (k==1) return AliPID::ParticleMass(AliPID::kMuon);
 550   if (k==2||k==-1) return AliPID::ParticleMass(AliPID::kPion);
 551   if (k==3) return AliPID::ParticleMass(AliPID::kKaon);
 552   if (k==4) return AliPID::ParticleMass(AliPID::kProton);
 553   AliWarning("Undefined mass !");
 554   return AliPID::ParticleMass(AliPID::kPion);
 555 }
 556
 557 //______________________________________________________________________________
 558 Double_t AliESDtrack::E() const
 559 {
 560   // Returns the energy of the particle given its assumed mass.
 561   // Assumes the pion mass if the particle can't be identified properly.
 562
 563   Double_t m = M();
 564   Double_t p = P();
 565   return TMath::Sqrt(p*p + m*m);
 566 }
 567
 568 //______________________________________________________________________________
 569 Double_t AliESDtrack::Y() const
 570 {
 571   // Returns the rapidity of a particle given its assumed mass.
 572   // Assumes the pion mass if the particle can't be identified properly.
 573
 574   Double_t e = E();
 575   Double_t pz = Pz();
 576   if (e != TMath::Abs(pz)) { // energy was not equal to pz
 577     return 0.5*TMath::Log((e+pz)/(e-pz));
 578   } else { // energy was equal to pz
 579     return -999.;
 580   }
 581 }
 582
 583 //_______________________________________________________________________
 584 Bool_t AliESDtrack::UpdateTrackParams(const AliKalmanTrack *t, ULong_t flags){
 585   //
 586   // This function updates track's running parameters
 587   //
 588   Int_t *index=0;
 589   Bool_t rc=kTRUE;
 590
 591   SetStatus(flags);
 592   fLabel=t->GetLabel();
 593
 594   if (t->IsStartedTimeIntegral()) {
 595     SetStatus(kTIME);
 596     Double_t times[10];t->GetIntegratedTimes(times); SetIntegratedTimes(times);
 597     SetIntegratedLength(t->GetIntegratedLength());
 598   }
 599
 600   Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
 601
 602   switch (flags) {
 603
 604   case kITSin: case kITSout: case kITSrefit:
 605     fITSClusterMap=0;
 606     fITSncls=t->GetNumberOfClusters();
 607     index=fFriendTrack->GetITSindices();
 608     for (Int_t i=0;i<AliESDfriendTrack::kMaxITScluster;i++) {
 609         index[i]=t->GetClusterIndex(i);
 610         if (i<fITSncls) {
 611            Int_t l=(index[i] & 0xf0000000) >> 28;
 612            SETBIT(fITSClusterMap,l);
 613         }
 614     }
 615     fITSchi2=t->GetChi2();
 616     fITSsignal=t->GetPIDsignal();
 617     fITSLabel = t->GetLabel();
 618     break;
 619
 620   case kTPCin: case kTPCrefit:
 621     fTPCLabel = t->GetLabel();
 622     if (flags==kTPCin)  fTPCInner=new AliExternalTrackParam(*t);
 623     if (!fIp) fIp=new AliExternalTrackParam(*t);
 624     else
 625       fIp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
 626   case kTPCout:
 627     index=fFriendTrack->GetTPCindices();
 628     if (flags & kTPCout){
 629       if (!fOp) fOp=new AliExternalTrackParam(*t);
 630       else
 631         fOp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
 632     }
 633     fTPCncls=t->GetNumberOfClusters();
 634     fTPCchi2=t->GetChi2();
 635
 636      {//prevrow must be declared in separate namespace, otherwise compiler cries:
 637       //"jump to case label crosses initialization of `Int_t prevrow'"
 638        Int_t prevrow = -1;
 639        //       for (Int_t i=0;i<fTPCncls;i++)
 640        for (Int_t i=0;i<AliESDfriendTrack::kMaxTPCcluster;i++)
 641         {
 642           index[i]=t->GetClusterIndex(i);
 643           Int_t idx = index[i];
 644
 645           if (idx<0) continue;
 646
 647           // Piotr's Cluster Map for HBT
 648           // ### please change accordingly if cluster array is changing
 649           // to "New TPC Tracking" style (with gaps in array)
 650           Int_t sect = (idx&0xff000000)>>24;
 651           Int_t row = (idx&0x00ff0000)>>16;
 652           if (sect > 18) row +=63; //if it is outer sector, add number of inner sectors
 653
 654           fTPCClusterMap.SetBitNumber(row,kTRUE);
 655
 656           //Fill the gap between previous row and this row with 0 bits
 657           //In case  ###  pleas change it as well - just set bit 0 in case there
 658           //is no associated clusters for current "i"
 659           if (prevrow < 0)
 660            {
 661              prevrow = row;//if previous bit was not assigned yet == this is the first one
 662            }
 663           else
 664            { //we don't know the order (inner to outer or reverse)
 665              //just to be save in case it is going to change
 666              Int_t n = 0, m = 0;
 667              if (prevrow < row)
 668               {
 669                 n = prevrow;
 670                 m = row;
 671               }
 672              else
 673               {
 674                 n = row;
 675                 m = prevrow;
 676               }
 677
 678              for (Int_t j = n+1; j < m; j++)
 679               {
 680                 fTPCClusterMap.SetBitNumber(j,kFALSE);
 681               }
 682              prevrow = row;
 683            }
 684           // End Of Piotr's Cluster Map for HBT
 685         }
 686      }
 687     fTPCsignal=t->GetPIDsignal();
 688     break;
 689
 690   case kTRDout: case kTRDin: case kTRDrefit:
 691     index=fFriendTrack->GetTRDindices();
 692     fTRDLabel = t->GetLabel();
 693     fTRDchi2=t->GetChi2();
 694     fTRDncls=6;//t->GetNumberOfTracklets(); //t->GetNumberOfClusters();
 695     //for (Int_t i=0;i<fTRDncls;i++) index[i]=t->GetClusterIndex(i);
 696     for (Int_t i=0;i<6;i++) index[i]=t->GetTrackletIndex(i);
 697
 698     fTRDsignal=t->GetPIDsignal();
 699     break;
 700   case kTRDbackup:
 701     if (!fOp) fOp=new AliExternalTrackParam(*t);
 702     else
 703       fOp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
 704     fTRDncls0 = t->GetNumberOfClusters();
 705     break;
 706   case kTOFin:
 707     break;
 708   case kTOFout:
 709     break;
 710   case kTRDStop:
 711     break;
 712   default:
 713     AliError("Wrong flag !");
 714     return kFALSE;
 715   }
 716
 717   return rc;
 718 }
 719
 720 //_______________________________________________________________________
 721 void AliESDtrack::GetExternalParameters(Double_t &x, Double_t p[5]) const {
 722   //---------------------------------------------------------------------
 723   // This function returns external representation of the track parameters
 724   //---------------------------------------------------------------------
 725   x=GetX();
 726   for (Int_t i=0; i<5; i++) p[i]=GetParameter()[i];
 727 }
 728
 729 //_______________________________________________________________________
 730 void AliESDtrack::GetExternalCovariance(Double_t cov[15]) const {
 731   //---------------------------------------------------------------------
 732   // This function returns external representation of the cov. matrix
 733   //---------------------------------------------------------------------
 734   for (Int_t i=0; i<15; i++) cov[i]=AliExternalTrackParam::GetCovariance()[i];
 735 }
 736
 737 //_______________________________________________________________________
 738 Bool_t AliESDtrack::GetConstrainedExternalParameters
 739                  (Double_t &alpha, Double_t &x, Double_t p[5]) const {
 740   //---------------------------------------------------------------------
 741   // This function returns the constrained external track parameters
 742   //---------------------------------------------------------------------
 743   if (!fCp) return kFALSE;
 744   alpha=fCp->GetAlpha();
 745   x=fCp->GetX();
 746   for (Int_t i=0; i<5; i++) p[i]=fCp->GetParameter()[i];
 747   return kTRUE;
 748 }
 749
 750 //_______________________________________________________________________
 751 Bool_t
 752 AliESDtrack::GetConstrainedExternalCovariance(Double_t c[15]) const {
 753   //---------------------------------------------------------------------
 754   // This function returns the constrained external cov. matrix
 755   //---------------------------------------------------------------------
 756   if (!fCp) return kFALSE;
 757   for (Int_t i=0; i<15; i++) c[i]=fCp->GetCovariance()[i];
 758   return kTRUE;
 759 }
 760
 761 Bool_t
 762 AliESDtrack::GetInnerExternalParameters
 763                  (Double_t &alpha, Double_t &x, Double_t p[5]) const {
 764   //---------------------------------------------------------------------
 765   // This function returns external representation of the track parameters
 766   // at the inner layer of TPC
 767   //---------------------------------------------------------------------
 768   if (!fIp) return kFALSE;
 769   alpha=fIp->GetAlpha();
 770   x=fIp->GetX();
 771   for (Int_t i=0; i<5; i++) p[i]=fIp->GetParameter()[i];
 772   return kTRUE;
 773 }
 774
 775 Bool_t
 776 AliESDtrack::GetInnerExternalCovariance(Double_t cov[15]) const {
 777  //---------------------------------------------------------------------
 778  // This function returns external representation of the cov. matrix
 779  // at the inner layer of TPC
 780  //---------------------------------------------------------------------
 781   if (!fIp) return kFALSE;
 782   for (Int_t i=0; i<15; i++) cov[i]=fIp->GetCovariance()[i];
 783   return kTRUE;
 784 }
 785
 786 Bool_t
 787 AliESDtrack::GetOuterExternalParameters
 788                  (Double_t &alpha, Double_t &x, Double_t p[5]) const {
 789   //---------------------------------------------------------------------
 790   // This function returns external representation of the track parameters
 791   // at the inner layer of TRD
 792   //---------------------------------------------------------------------
 793   if (!fOp) return kFALSE;
 794   alpha=fOp->GetAlpha();
 795   x=fOp->GetX();
 796   for (Int_t i=0; i<5; i++) p[i]=fOp->GetParameter()[i];
 797   return kTRUE;
 798 }
 799
 800 Bool_t
 801 AliESDtrack::GetOuterExternalCovariance(Double_t cov[15]) const {
 802  //---------------------------------------------------------------------
 803  // This function returns external representation of the cov. matrix
 804  // at the inner layer of TRD
 805  //---------------------------------------------------------------------
 806   if (!fOp) return kFALSE;
 807   for (Int_t i=0; i<15; i++) cov[i]=fOp->GetCovariance()[i];
 808   return kTRUE;
 809 }
 810
 811 Int_t AliESDtrack::GetNcls(Int_t idet) const
 812 {
 813   // Get number of clusters by subdetector index
 814   //
 815   Int_t ncls = 0;
 816   switch(idet){
 817   case 0:
 818     ncls = fITSncls;
 819     break;
 820   case 1:
 821     ncls = fTPCncls;
 822     break;
 823   case 2:
 824     ncls = fTRDncls;
 825     break;
 826   case 3:
 827     if (fTOFindex != 0)
 828       ncls = 1;
 829     break;
 830   default:
 831     break;
 832   }
 833   return ncls;
 834 }
 835
 836 Int_t AliESDtrack::GetClusters(Int_t idet, Int_t *idx) const
 837 {
 838   // Get cluster index array by subdetector index
 839   //
 840   Int_t ncls = 0;
 841   switch(idet){
 842   case 0:
 843     ncls = GetITSclusters(idx);
 844     break;
 845   case 1:
 846     ncls = GetTPCclusters(idx);
 847     break;
 848   case 2:
 849     ncls = GetTRDclusters(idx);
 850     break;
 851   case 3:
 852     if (fTOFindex != 0) {
 853       idx[0] = GetTOFcluster();
 854       ncls = 1;
 855     }
 856     break;
 857   default:
 858     break;
 859   }
 860   return ncls;
 861 }
 862
 863 //_______________________________________________________________________
 864 void AliESDtrack::GetIntegratedTimes(Double_t *times) const {
 865   // Returns the array with integrated times for each particle hypothesis
 866   for (Int_t i=0; i<AliPID::kSPECIES; i++) times[i]=fTrackTime[i];
 867 }
 868
 869 //_______________________________________________________________________
 870 void AliESDtrack::SetIntegratedTimes(const Double_t *times) {
 871   // Sets the array with integrated times for each particle hypotesis
 872   for (Int_t i=0; i<AliPID::kSPECIES; i++) fTrackTime[i]=times[i];
 873 }
 874
 875 //_______________________________________________________________________
 876 void AliESDtrack::SetITSpid(const Double_t *p) {
 877   // Sets values for the probability of each particle type (in ITS)
 878   SetPIDValues(fITSr,p,AliPID::kSPECIES);
 879   SetStatus(AliESDtrack::kITSpid);
 880 }
 881
 882 //_______________________________________________________________________
 883 void AliESDtrack::GetITSpid(Double_t *p) const {
 884   // Gets the probability of each particle type (in ITS)
 885   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fITSr[i];
 886 }
 887
 888 //_______________________________________________________________________
 889 Char_t AliESDtrack::GetITSclusters(Int_t *idx) const {
 890   //---------------------------------------------------------------------
 891   // This function returns indices of the assgined ITS clusters
 892   //---------------------------------------------------------------------
 893   if (idx!=0) {
 894      Int_t *index=fFriendTrack->GetITSindices();
 895      for (Int_t i=0; i<AliESDfriendTrack::kMaxITScluster; i++) idx[i]=index[i];
 896   }
 897   return fITSncls;
 898 }
 899
 900 //_______________________________________________________________________
 901 Bool_t AliESDtrack::GetITSModuleIndexInfo(Int_t ilayer,Int_t &idet,Int_t &status,
 902                                          Float_t &xloc,Float_t &zloc) const {
 903   //----------------------------------------------------------------------
 904   // This function encodes in the module number also the status of cluster association
 905   // "status" can have the following values:
 906   // 1 "found" (cluster is associated),
 907   // 2 "dead" (module is dead from OCDB),
 908   // 3 "skipped" (module or layer forced to be skipped),
 909   // 4 "outinz" (track out of z acceptance),
 910   // 5 "nocls" (no clusters in the road),
 911   // 6 "norefit" (cluster rejected during refit),
 912   // 7 "deadzspd" (holes in z in SPD)
 913   // Also given are the coordinates of the crossing point of track and module
 914   // (in the local module ref. system)
 915   // WARNING: THIS METHOD HAS TO BE SYNCHRONIZED WITH AliITStrackV2::GetModuleIndexInfo()!
 916   //----------------------------------------------------------------------
 917
 918   if(fITSModule[ilayer]==-1) {
 919     AliError("fModule was not set !");
 920     idet = -1;
 921     status=0;
 922     xloc=-99.; zloc=-99.;
 923     return kFALSE;
 924   }
 925
 926   Int_t module = fITSModule[ilayer];
 927
 928   idet = Int_t(module/1000000);
 929
 930   module -= idet*1000000;
 931
 932   status = Int_t(module/100000);
 933
 934   module -= status*100000;
 935
 936   Int_t signs = Int_t(module/10000);
 937
 938   module-=signs*10000;
 939
 940   Int_t xInt = Int_t(module/100);
 941   module -= xInt*100;
 942
 943   Int_t zInt = module;
 944
 945   if(signs==1) { xInt*=1; zInt*=1; }
 946   if(signs==2) { xInt*=1; zInt*=-1; }
 947   if(signs==3) { xInt*=-1; zInt*=1; }
 948   if(signs==4) { xInt*=-1; zInt*=-1; }
 949
 950   xloc = 0.1*(Float_t)xInt;
 951   zloc = 0.1*(Float_t)zInt;
 952
 953   if(status==4) idet = -1;
 954
 955   return kTRUE;
 956 }
 957
 958 //_______________________________________________________________________
 959 UShort_t AliESDtrack::GetTPCclusters(Int_t *idx) const {
 960   //---------------------------------------------------------------------
 961   // This function returns indices of the assgined ITS clusters
 962   //---------------------------------------------------------------------
 963   if (idx!=0) {
 964     Int_t *index=fFriendTrack->GetTPCindices();
 965     for (Int_t i=0; i<AliESDfriendTrack::kMaxTPCcluster; i++) idx[i]=index[i];
 966   }
 967   return fTPCncls;
 968 }
 969
 970 Double_t AliESDtrack::GetTPCdensity(Int_t row0, Int_t row1) const{
 971   //
 972   // GetDensity of the clusters on given region between row0 and row1
 973   // Dead zone effect takin into acoount
 974   //
 975   Int_t good  = 0;
 976   Int_t found = 0;
 977   //
 978   Int_t *index=fFriendTrack->GetTPCindices();
 979   for (Int_t i=row0;i<=row1;i++){
 980     Int_t idx = index[i];
 981     if (idx!=-1)  good++;             // track outside of dead zone
 982     if (idx>0)    found++;
 983   }
 984   Float_t density=0.5;
 985   if (good>(row1-row0)*0.5) density = Float_t(found)/Float_t(good);
 986   return density;
 987 }
 988
 989 //_______________________________________________________________________
 990 void AliESDtrack::SetTPCpid(const Double_t *p) {
 991   // Sets values for the probability of each particle type (in TPC)
 992   SetPIDValues(fTPCr,p,AliPID::kSPECIES);
 993   SetStatus(AliESDtrack::kTPCpid);
 994 }
 995
 996 //_______________________________________________________________________
 997 void AliESDtrack::GetTPCpid(Double_t *p) const {
 998   // Gets the probability of each particle type (in TPC)
 999   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fTPCr[i];
1000 }
1001
1002 //_______________________________________________________________________
1003 UChar_t AliESDtrack::GetTRDclusters(Int_t *idx) const {
1004   //---------------------------------------------------------------------
1005   // This function returns indices of the assgined TRD clusters
1006   //---------------------------------------------------------------------
1007   if (idx!=0) {
1008      Int_t *index=fFriendTrack->GetTRDindices();
1009      for (Int_t i=0; i<AliESDfriendTrack::kMaxTRDcluster; i++) idx[i]=index[i];
1010   }
1011   return fTRDncls;
1012 }
1013
1014 //_______________________________________________________________________
1015 UChar_t AliESDtrack::GetTRDtracklets(Int_t *idx) const {
1016   //---------------------------------------------------------------------
1017   // This function returns indices of the assigned TRD tracklets
1018   //---------------------------------------------------------------------
1019   if (idx!=0) {
1020      Int_t *index=fFriendTrack->GetTRDindices();
1021      for (Int_t i=0; i<6/*AliESDfriendTrack::kMaxTRDcluster*/; i++) idx[i]=index[i];
1022   }
1023   return fTRDncls;
1024 }
1025
1026 //_______________________________________________________________________
1027 void AliESDtrack::SetTRDpid(const Double_t *p) {
1028   // Sets values for the probability of each particle type (in TRD)
1029   SetPIDValues(fTRDr,p,AliPID::kSPECIES);
1030   SetStatus(AliESDtrack::kTRDpid);
1031 }
1032
1033 //_______________________________________________________________________
1034 void AliESDtrack::GetTRDpid(Double_t *p) const {
1035   // Gets the probability of each particle type (in TRD)
1036   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fTRDr[i];
1037 }
1038
1039 //_______________________________________________________________________
1040 void    AliESDtrack::SetTRDpid(Int_t iSpecies, Float_t p)
1041 {
1042   // Sets the probability of particle type iSpecies to p (in TRD)
1043   fTRDr[iSpecies] = p;
1044 }
1045
1046 Double_t AliESDtrack::GetTRDpid(Int_t iSpecies) const
1047 {
1048   // Returns the probability of particle type iSpecies (in TRD)
1049   return fTRDr[iSpecies];
1050 }
1051
1052 //_______________________________________________________________________
1053 void AliESDtrack::SetTOFpid(const Double_t *p) {
1054   // Sets the probability of each particle type (in TOF)
1055   SetPIDValues(fTOFr,p,AliPID::kSPECIES);
1056   SetStatus(AliESDtrack::kTOFpid);
1057 }
1058
1059 //_______________________________________________________________________
1060 void AliESDtrack::SetTOFLabel(const Int_t *p) {
1061   // Sets  (in TOF)
1062   for (Int_t i=0; i<3; i++) fTOFLabel[i]=p[i];
1063 }
1064
1065 //_______________________________________________________________________
1066 void AliESDtrack::GetTOFpid(Double_t *p) const {
1067   // Gets probabilities of each particle type (in TOF)
1068   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fTOFr[i];
1069 }
1070
1071 //_______________________________________________________________________
1072 void AliESDtrack::GetTOFLabel(Int_t *p) const {
1073   // Gets (in TOF)
1074   for (Int_t i=0; i<3; i++) p[i]=fTOFLabel[i];
1075 }
1076
1077 //_______________________________________________________________________
1078 void AliESDtrack::GetTOFInfo(Float_t *info) const {
1079   // Gets (in TOF)
1080   for (Int_t i=0; i<10; i++) info[i]=fTOFInfo[i];
1081 }
1082
1083 //_______________________________________________________________________
1084 void AliESDtrack::SetTOFInfo(Float_t*info) {
1085   // Gets (in TOF)
1086   for (Int_t i=0; i<10; i++) fTOFInfo[i]=info[i];
1087 }
1088
1089
1090
1091 //_______________________________________________________________________
1092 void AliESDtrack::SetHMPIDpid(const Double_t *p) {
1093   // Sets the probability of each particle type (in HMPID)
1094   SetPIDValues(fHMPIDr,p,AliPID::kSPECIES);
1095   SetStatus(AliESDtrack::kHMPIDpid);
1096 }
1097
1098 //_______________________________________________________________________
1099 void AliESDtrack::GetHMPIDpid(Double_t *p) const {
1100   // Gets probabilities of each particle type (in HMPID)
1101   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fHMPIDr[i];
1102 }
1103
1104
1105
1106 //_______________________________________________________________________
1107 void AliESDtrack::SetESDpid(const Double_t *p) {
1108   // Sets the probability of each particle type for the ESD track
1109   SetPIDValues(fR,p,AliPID::kSPECIES);
1110   SetStatus(AliESDtrack::kESDpid);
1111 }
1112
1113 //_______________________________________________________________________
1114 void AliESDtrack::GetESDpid(Double_t *p) const {
1115   // Gets probability of each particle type for the ESD track
1116   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fR[i];
1117 }
1118
1119 //_______________________________________________________________________
1120 Bool_t AliESDtrack::RelateToVertex
1121 (const AliESDVertex *vtx, Double_t b, Double_t maxd) {
1122   //
1123   // Try to relate this track to the vertex "vtx",
1124   // if the (rough) transverse impact parameter is not bigger then "maxd".
1125   //            Magnetic field is "b" (kG).
1126   //
1127   // a) The track gets extapolated to the DCA to the vertex.
1128   // b) The impact parameters and their covariance matrix are calculated.
1129   // c) An attempt to constrain this track to the vertex is done.
1130   //
1131   //    In the case of success, the returned value is kTRUE
1132   //    (otherwise, it's kFALSE)
1133   //
1134
1135   if (!vtx) return kFALSE;
1136
1137   Double_t dz[2],cov[3];
1138   if (!PropagateToDCA(vtx, b, maxd, dz, cov)) return kFALSE;
1139
1140   fD = dz[0];
1141   fZ = dz[1];
1142   fCdd = cov[0];
1143   fCdz = cov[1];
1144   fCzz = cov[2];
1145
1146   Double_t covar[6]; vtx->GetCovMatrix(covar);
1147   Double_t p[2]={GetParameter()[0]-dz[0],GetParameter()[1]-dz[1]};
1148   Double_t c[3]={covar[2],0.,covar[5]};
1149
1150   Double_t chi2=GetPredictedChi2(p,c);
1151   if (chi2>77.) return kFALSE;
1152
1153   delete fCp;
1154   fCp=new AliExternalTrackParam(*this);
1155
1156   if (!fCp->Update(p,c)) {delete fCp; fCp=0; return kFALSE;}
1157
1158   fCchi2=chi2;
1159   return kTRUE;
1160 }
1161
1162 //_______________________________________________________________________
1163 void AliESDtrack::Print(Option_t *) const {
1164   // Prints info on the track
1165
1166   printf("ESD track info\n") ;
1167   Double_t p[AliPID::kSPECIESN] ;
1168   Int_t index = 0 ;
1169   if( IsOn(kITSpid) ){
1170     printf("From ITS: ") ;
1171     GetITSpid(p) ;
1172     for(index = 0 ; index < AliPID::kSPECIES; index++)
1173       printf("%f, ", p[index]) ;
1174     printf("\n           signal = %f\n", GetITSsignal()) ;
1175   }
1176   if( IsOn(kTPCpid) ){
1177     printf("From TPC: ") ;
1178     GetTPCpid(p) ;
1179     for(index = 0 ; index < AliPID::kSPECIES; index++)
1180       printf("%f, ", p[index]) ;
1181     printf("\n           signal = %f\n", GetTPCsignal()) ;
1182   }
1183   if( IsOn(kTRDpid) ){
1184     printf("From TRD: ") ;
1185     GetTRDpid(p) ;
1186     for(index = 0 ; index < AliPID::kSPECIES; index++)
1187       printf("%f, ", p[index]) ;
1188     printf("\n           signal = %f\n", GetTRDsignal()) ;
1189   }
1190   if( IsOn(kTOFpid) ){
1191     printf("From TOF: ") ;
1192     GetTOFpid(p) ;
1193     for(index = 0 ; index < AliPID::kSPECIES; index++)
1194       printf("%f, ", p[index]) ;
1195     printf("\n           signal = %f\n", GetTOFsignal()) ;
1196   }
1197   if( IsOn(kHMPIDpid) ){
1198     printf("From HMPID: ") ;
1199     GetHMPIDpid(p) ;
1200     for(index = 0 ; index < AliPID::kSPECIES; index++)
1201       printf("%f, ", p[index]) ;
1202     printf("\n           signal = %f\n", GetHMPIDsignal()) ;
1203   }
1204 }
1205