STEER/AOD/AliAODTrack.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-2007, 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 /* $Id$ */
  17
  18 //-------------------------------------------------------------------------
  19 //     AOD track implementation of AliVTrack
  20 //     Author: Markus Oldenburg, CERN
  21 //     Markus.Oldenburg@cern.ch
  22 //-------------------------------------------------------------------------
  23
  24 #include <TVector3.h>
  25 #include "AliLog.h"
  26 #include "AliExternalTrackParam.h"
  27 #include "AliVVertex.h"
  28 #include "AliAODTrack.h"
  29 #include "AliAODEvent.h"
  30
  31 ClassImp(AliAODTrack)
  32
  33 //______________________________________________________________________________
  34 AliAODTrack::AliAODTrack() :
  35   AliVTrack(),
  36   fRAtAbsorberEnd(0.),
  37   fChi2perNDF(-999.),
  38   fChi2MatchTrigger(0.),
  39   fFlags(0),
  40   fLabel(-999),
  41   fITSMuonClusterMap(0),
  42   fMUONtrigHitsMapTrg(0),
  43   fMUONtrigHitsMapTrk(0),
  44   fFilterMap(0),
  45   fTPCFitMap(),
  46   fTPCClusterMap(),
  47   fTPCSharedMap(),
  48   fTPCnclsF(0),
  49   fTPCNCrossedRows(0),
  50   fID(-999),
  51   fCharge(-99),
  52   fType(kUndef),
  53   fCaloIndex(kEMCALNoMatch),
  54   fCovMatrix(NULL),
  55   fDetPid(NULL),
  56   fProdVertex(NULL),
  57   fTrackPhiOnEMCal(-999),
  58   fTrackEtaOnEMCal(-999),
  59   fTPCsignalTuned(0),
  60   fAODEvent(NULL)
  61 {
  62   // default constructor
  63
  64   SetP();
  65   SetPosition((Float_t*)NULL);
  66   SetXYAtDCA(-999., -999.);
  67   SetPxPyPzAtDCA(-999., -999., -999.);
  68   SetPID((Float_t*)NULL);
  69 }
  70
  71 //______________________________________________________________________________
  72 AliAODTrack::AliAODTrack(Short_t id,
  73                          Int_t label,
  74                          Double_t p[3],
  75                          Bool_t cartesian,
  76                          Double_t x[3],
  77                          Bool_t isDCA,
  78                          Double_t covMatrix[21],
  79                          Short_t charge,
  80                          UChar_t itsClusMap,
  81                          Double_t pid[10],
  82                          AliAODVertex *prodVertex,
  83                          Bool_t usedForVtxFit,
  84                          Bool_t usedForPrimVtxFit,
  85                          AODTrk_t ttype,
  86                          UInt_t selectInfo,
  87                          Float_t chi2perNDF) :
  88   AliVTrack(),
  89   fRAtAbsorberEnd(0.),
  90   fChi2perNDF(chi2perNDF),
  91   fChi2MatchTrigger(0.),
  92   fFlags(0),
  93   fLabel(label),
  94   fITSMuonClusterMap(0),
  95   fMUONtrigHitsMapTrg(0),
  96   fMUONtrigHitsMapTrk(0),
  97   fFilterMap(selectInfo),
  98   fTPCFitMap(),
  99   fTPCClusterMap(),
 100   fTPCSharedMap(),
 101   fTPCnclsF(0),
 102   fTPCNCrossedRows(0),
 103   fID(id),
 104   fCharge(charge),
 105   fType(ttype),
 106   fCaloIndex(kEMCALNoMatch),
 107   fCovMatrix(NULL),
 108   fDetPid(NULL),
 109   fProdVertex(prodVertex),
 110   fTrackPhiOnEMCal(-999),
 111   fTrackEtaOnEMCal(-999),
 112   fTPCsignalTuned(0),
 113   fAODEvent(NULL)
 114 {
 115   // constructor
 116
 117   SetP(p, cartesian);
 118   SetPosition(x, isDCA);
 119   SetXYAtDCA(-999., -999.);
 120   SetPxPyPzAtDCA(-999., -999., -999.);
 121   SetUsedForVtxFit(usedForVtxFit);
 122   SetUsedForPrimVtxFit(usedForPrimVtxFit);
 123   if(covMatrix) SetCovMatrix(covMatrix);
 124   SetPID(pid);
 125   SetITSClusterMap(itsClusMap);
 126 }
 127
 128 //______________________________________________________________________________
 129 AliAODTrack::AliAODTrack(Short_t id,
 130                          Int_t label,
 131                          Float_t p[3],
 132                          Bool_t cartesian,
 133                          Float_t x[3],
 134                          Bool_t isDCA,
 135                          Float_t covMatrix[21],
 136                          Short_t charge,
 137                          UChar_t itsClusMap,
 138                          Float_t pid[10],
 139                          AliAODVertex *prodVertex,
 140                          Bool_t usedForVtxFit,
 141                          Bool_t usedForPrimVtxFit,
 142                          AODTrk_t ttype,
 143                          UInt_t selectInfo,
 144                          Float_t chi2perNDF) :
 145   AliVTrack(),
 146   fRAtAbsorberEnd(0.),
 147   fChi2perNDF(chi2perNDF),
 148   fChi2MatchTrigger(0.),
 149   fFlags(0),
 150   fLabel(label),
 151   fITSMuonClusterMap(0),
 152   fMUONtrigHitsMapTrg(0),
 153   fMUONtrigHitsMapTrk(0),
 154   fFilterMap(selectInfo),
 155   fTPCFitMap(),
 156   fTPCClusterMap(),
 157   fTPCSharedMap(),
 158   fTPCnclsF(0),
 159   fTPCNCrossedRows(0),
 160   fID(id),
 161   fCharge(charge),
 162   fType(ttype),
 163   fCaloIndex(kEMCALNoMatch),
 164   fCovMatrix(NULL),
 165   fDetPid(NULL),
 166   fProdVertex(prodVertex),
 167   fTrackPhiOnEMCal(-999),
 168   fTrackEtaOnEMCal(-999),
 169   fTPCsignalTuned(0),
 170   fAODEvent(NULL)
 171 {
 172   // constructor
 173
 174   SetP(p, cartesian);
 175   SetPosition(x, isDCA);
 176   SetXYAtDCA(-999., -999.);
 177   SetPxPyPzAtDCA(-999., -999., -999.);
 178   SetUsedForVtxFit(usedForVtxFit);
 179   SetUsedForPrimVtxFit(usedForPrimVtxFit);
 180   if(covMatrix) SetCovMatrix(covMatrix);
 181   SetPID(pid);
 182   SetITSClusterMap(itsClusMap);
 183 }
 184
 185 //______________________________________________________________________________
 186 AliAODTrack::~AliAODTrack()
 187 {
 188   // destructor
 189   delete fCovMatrix;
 190   delete fDetPid;
 191 }
 192
 193
 194 //______________________________________________________________________________
 195 AliAODTrack::AliAODTrack(const AliAODTrack& trk) :
 196   AliVTrack(trk),
 197   fRAtAbsorberEnd(trk.fRAtAbsorberEnd),
 198   fChi2perNDF(trk.fChi2perNDF),
 199   fChi2MatchTrigger(trk.fChi2MatchTrigger),
 200   fFlags(trk.fFlags),
 201   fLabel(trk.fLabel),
 202   fITSMuonClusterMap(trk.fITSMuonClusterMap),
 203   fMUONtrigHitsMapTrg(trk.fMUONtrigHitsMapTrg),
 204   fMUONtrigHitsMapTrk(trk.fMUONtrigHitsMapTrk),
 205   fFilterMap(trk.fFilterMap),
 206   fTPCFitMap(trk.fTPCFitMap),
 207   fTPCClusterMap(trk.fTPCClusterMap),
 208   fTPCSharedMap(trk.fTPCSharedMap),
 209   fTPCnclsF(trk.fTPCnclsF),
 210   fTPCNCrossedRows(trk.fTPCNCrossedRows),
 211   fID(trk.fID),
 212   fCharge(trk.fCharge),
 213   fType(trk.fType),
 214   fCaloIndex(trk.fCaloIndex),
 215   fCovMatrix(NULL),
 216   fDetPid(NULL),
 217   fProdVertex(trk.fProdVertex),
 218   fTrackPhiOnEMCal(trk.fTrackPhiOnEMCal),
 219   fTrackEtaOnEMCal(trk.fTrackEtaOnEMCal),
 220   fTPCsignalTuned(trk.fTPCsignalTuned),
 221   fAODEvent(trk.fAODEvent)
 222 {
 223   // Copy constructor
 224
 225   trk.GetP(fMomentum);
 226   trk.GetPosition(fPosition);
 227   SetXYAtDCA(trk.XAtDCA(), trk.YAtDCA());
 228   SetPxPyPzAtDCA(trk.PxAtDCA(), trk.PyAtDCA(), trk.PzAtDCA());
 229   SetUsedForVtxFit(trk.GetUsedForVtxFit());
 230   SetUsedForPrimVtxFit(trk.GetUsedForPrimVtxFit());
 231   if(trk.fCovMatrix) fCovMatrix=new AliAODRedCov<6>(*trk.fCovMatrix);
 232   if(trk.fDetPid) fDetPid=new AliAODPid(*trk.fDetPid);
 233   SetPID(trk.fPID);
 234 }
 235
 236 //______________________________________________________________________________
 237 AliAODTrack& AliAODTrack::operator=(const AliAODTrack& trk)
 238 {
 239   // Assignment operator
 240   if(this!=&trk) {
 241
 242     AliVTrack::operator=(trk);
 243
 244     trk.GetP(fMomentum);
 245     trk.GetPosition(fPosition);
 246     SetXYAtDCA(trk.XAtDCA(), trk.YAtDCA());
 247     SetPxPyPzAtDCA(trk.PxAtDCA(), trk.PyAtDCA(), trk.PzAtDCA());
 248     fRAtAbsorberEnd    = trk.fRAtAbsorberEnd;
 249     fChi2perNDF        = trk.fChi2perNDF;
 250     fChi2MatchTrigger  = trk.fChi2MatchTrigger;
 251     trk.GetPID(fPID);
 252     fFlags             = trk.fFlags;
 253     fLabel             = trk.fLabel;
 254     fITSMuonClusterMap = trk.fITSMuonClusterMap;
 255     fMUONtrigHitsMapTrg = trk.fMUONtrigHitsMapTrg;
 256     fMUONtrigHitsMapTrk = trk.fMUONtrigHitsMapTrk;
 257     fFilterMap         = trk.fFilterMap;
 258     fTPCFitMap         = trk.fTPCFitMap;
 259     fTPCClusterMap     = trk.fTPCClusterMap;
 260     fTPCSharedMap      = trk.fTPCSharedMap;
 261     fTPCnclsF          = trk.fTPCnclsF;
 262     fTPCNCrossedRows   = trk.fTPCNCrossedRows;
 263     fID                = trk.fID;
 264     fCharge            = trk.fCharge;
 265     fType              = trk.fType;
 266     fCaloIndex         = trk.fCaloIndex;
 267     fTrackPhiOnEMCal   = trk.fTrackPhiOnEMCal;
 268     fTrackEtaOnEMCal   = trk.fTrackEtaOnEMCal;
 269     fTPCsignalTuned    = trk.fTPCsignalTuned;
 270
 271     delete fCovMatrix;
 272     if(trk.fCovMatrix) fCovMatrix=new AliAODRedCov<6>(*trk.fCovMatrix);
 273     else fCovMatrix=NULL;
 274
 275
 276     fProdVertex        = trk.fProdVertex;
 277     SetUsedForVtxFit(trk.GetUsedForVtxFit());
 278     SetUsedForPrimVtxFit(trk.GetUsedForPrimVtxFit());
 279
 280     delete fDetPid;
 281     if(trk.fDetPid) fDetPid=new AliAODPid(*trk.fDetPid);
 282     else fDetPid=NULL;
 283
 284
 285   }
 286
 287   return *this;
 288 }
 289
 290 //______________________________________________________________________________
 291 Double_t AliAODTrack::M(AODTrkPID_t pid) const
 292 {
 293   // Returns the mass.
 294   // Masses for nuclei don't exist in the PDG tables, therefore they were put by hand.
 295
 296   switch (pid) {
 297
 298   case kElectron :
 299     return 0.000510999; //TDatabasePDG::Instance()->GetParticle(11/*::kElectron*/)->Mass();
 300     break;
 301
 302   case kMuon :
 303     return 0.1056584; //TDatabasePDG::Instance()->GetParticle(13/*::kMuonMinus*/)->Mass();
 304     break;
 305
 306   case kPion :
 307     return 0.13957; //TDatabasePDG::Instance()->GetParticle(211/*::kPiPlus*/)->Mass();
 308     break;
 309
 310   case kKaon :
 311     return 0.4937; //TDatabasePDG::Instance()->GetParticle(321/*::kKPlus*/)->Mass();
 312     break;
 313
 314   case kProton :
 315     return 0.9382720; //TDatabasePDG::Instance()->GetParticle(2212/*::kProton*/)->Mass();
 316     break;
 317
 318   case kDeuteron :
 319     return 1.8756; //TDatabasePDG::Instance()->GetParticle(1000010020)->Mass();
 320     break;
 321
 322   case kTriton :
 323     return 2.8089; //TDatabasePDG::Instance()->GetParticle(1000010030)->Mass();
 324     break;
 325
 326   case kHelium3 :
 327     return 2.8084; //TDatabasePDG::Instance()->GetParticle(1000020030)->Mass();
 328     break;
 329
 330   case kAlpha :
 331     return 3.7274; //TDatabasePDG::Instance()->GetParticle(1000020040)->Mass();
 332     break;
 333
 334   case kUnknown :
 335     return -999.;
 336     break;
 337
 338   default :
 339     return -999.;
 340   }
 341 }
 342
 343 //______________________________________________________________________________
 344 Double_t AliAODTrack::E(AODTrkPID_t pid) const
 345 {
 346   // Returns the energy of the particle of a given pid.
 347
 348   if (pid != kUnknown) { // particle was identified
 349     Double_t m = M(pid);
 350     return TMath::Sqrt(P()*P() + m*m);
 351   } else { // pid unknown
 352     return -999.;
 353   }
 354 }
 355
 356 //______________________________________________________________________________
 357 Double_t AliAODTrack::Y(AODTrkPID_t pid) const
 358 {
 359   // Returns the rapidity of a particle of a given pid.
 360
 361   if (pid != kUnknown) { // particle was identified
 362     Double_t e = E(pid);
 363     Double_t pz = Pz();
 364     if (e>=0 && e!=pz) { // energy was positive (e.g. not -999.) and not equal to pz
 365       return 0.5*TMath::Log((e+pz)/(e-pz));
 366     } else { // energy not known or equal to pz
 367       return -999.;
 368     }
 369   } else { // pid unknown
 370     return -999.;
 371   }
 372 }
 373
 374 //______________________________________________________________________________
 375 Double_t AliAODTrack::Y(Double_t m) const
 376 {
 377   // Returns the rapidity of a particle of a given mass.
 378
 379   if (m >= 0.) { // mass makes sense
 380     Double_t e = E(m);
 381     Double_t pz = Pz();
 382     if (e>=0 && e!=pz) { // energy was positive (e.g. not -999.) and not equal to pz
 383       return 0.5*TMath::Log((e+pz)/(e-pz));
 384     } else { // energy not known or equal to pz
 385       return -999.;
 386     }
 387   } else { // pid unknown
 388     return -999.;
 389   }
 390 }
 391
 392 //______________________________________________________________________________
 393 AliAODTrack::AODTrkPID_t AliAODTrack::GetMostProbablePID() const
 394 {
 395   // Returns the most probable PID array element.
 396
 397   Int_t nPID = 10;
 398   AODTrkPID_t loc = kUnknown;
 399   Double_t max = 0.;
 400   Bool_t allTheSame = kTRUE;
 401
 402   for (Int_t iPID = 0; iPID < nPID; iPID++) {
 403     if (fPID[iPID] >= max) {
 404       if (fPID[iPID] > max) {
 405         allTheSame = kFALSE;
 406         max = fPID[iPID];
 407         loc = (AODTrkPID_t)iPID;
 408       } else {
 409         allTheSame = kTRUE;
 410       }
 411     }
 412   }
 413   return allTheSame ? kUnknown : loc;
 414 }
 415
 416 //______________________________________________________________________________
 417 void AliAODTrack::ConvertAliPIDtoAODPID()
 418 {
 419   // Converts AliPID array.
 420   // The numbering scheme is the same for electrons, muons, pions, kaons, and protons.
 421   // Everything else has to be set to zero.
 422
 423   fPID[kDeuteron] = 0.;
 424   fPID[kTriton]   = 0.;
 425   fPID[kHelium3]  = 0.;
 426   fPID[kAlpha]    = 0.;
 427   fPID[kUnknown]  = 0.;
 428
 429   return;
 430 }
 431
 432
 433 //______________________________________________________________________________
 434 template <typename T> void AliAODTrack::SetP(const T *p, const Bool_t cartesian)
 435 {
 436   // Set the momentum
 437
 438   if (p) {
 439     if (cartesian) {
 440       Double_t pt2 = p[0]*p[0] + p[1]*p[1];
 441       Double_t pp  = TMath::Sqrt(pt2 + p[2]*p[2]);
 442
 443       fMomentum[0] = TMath::Sqrt(pt2); // pt
 444       fMomentum[1] = (pt2 != 0.) ? TMath::Pi()+TMath::ATan2(-p[1], -p[0]) : -999; // phi
 445       fMomentum[2] = (pp != 0.) ? TMath::ACos(p[2] / pp) : -999.; // theta
 446     } else {
 447       fMomentum[0] = p[0];  // pt
 448       fMomentum[1] = p[1];  // phi
 449       fMomentum[2] = p[2];  // theta
 450     }
 451   } else {
 452     fMomentum[0] = -999.;
 453     fMomentum[1] = -999.;
 454     fMomentum[2] = -999.;
 455   }
 456 }
 457
 458 /*
 459 //______________________________________________________________________________
 460 template <typename T> void AliAODTrack::SetPosition(const T *x, const Bool_t dca)
 461 {
 462   // set the position
 463
 464   if (x) {
 465     if (!dca) {
 466       ResetBit(kIsDCA);
 467
 468       fPosition[0] = x[0];
 469       fPosition[1] = x[1];
 470       fPosition[2] = x[2];
 471     } else {
 472       SetBit(kIsDCA);
 473       // don't know any better yet
 474       fPosition[0] = -999.;
 475       fPosition[1] = -999.;
 476       fPosition[2] = -999.;
 477     }
 478   } else {
 479     ResetBit(kIsDCA);
 480
 481     fPosition[0] = -999.;
 482     fPosition[1] = -999.;
 483     fPosition[2] = -999.;
 484   }
 485 }
 486 */
 487 //______________________________________________________________________________
 488 void AliAODTrack::SetDCA(Double_t d, Double_t z)
 489 {
 490   // set the dca
 491   fPosition[0] = d;
 492   fPosition[1] = z;
 493   fPosition[2] = 0.;
 494   SetBit(kIsDCA);
 495 }
 496
 497 //______________________________________________________________________________
 498 void AliAODTrack::Print(Option_t* /* option */) const
 499 {
 500   // prints information about AliAODTrack
 501
 502   printf("Object name: %s   Track type: %s\n", GetName(), GetTitle());
 503   printf("        px = %f\n", Px());
 504   printf("        py = %f\n", Py());
 505   printf("        pz = %f\n", Pz());
 506   printf("        pt = %f\n", Pt());
 507   printf("      1/pt = %f\n", OneOverPt());
 508   printf("     theta = %f\n", Theta());
 509   printf("       phi = %f\n", Phi());
 510   printf("  chi2/NDF = %f\n", Chi2perNDF());
 511   printf("    charge = %d\n", Charge());
 512 }
 513
 514 //______________________________________________________________________________
 515 void AliAODTrack::SetMatchTrigger(Int_t matchTrig)
 516 {
 517   // Set the MUON trigger information
 518   switch(matchTrig){
 519     case 0: // 0 track does not match trigger
 520       fITSMuonClusterMap=fITSMuonClusterMap&0x3fffffff;
 521       break;
 522     case 1: // 1 track match but does not pass pt cut
 523       fITSMuonClusterMap=(fITSMuonClusterMap&0x3fffffff)|0x40000000;
 524       break;
 525     case 2: // 2 track match Low pt cut
 526       fITSMuonClusterMap=(fITSMuonClusterMap&0x3fffffff)|0x80000000;
 527       break;
 528     case 3: // 3 track match High pt cut
 529       fITSMuonClusterMap=fITSMuonClusterMap|0xc0000000;
 530       break;
 531     default:
 532       fITSMuonClusterMap=fITSMuonClusterMap&0x3fffffff;
 533       AliWarning(Form("unknown case for matchTrig: %d\n",matchTrig));
 534   }
 535 }
 536
 537 //______________________________________________________________________________
 538 Bool_t AliAODTrack::HitsMuonChamber(Int_t MuonChamber, Int_t cathode) const
 539 {
 540   // return kTRUE if the track fires the given tracking or trigger chamber.
 541   // If the chamber is a trigger one:
 542   // - if cathode = 0 or 1, the track matches the corresponding cathode
 543   // - if cathode = -1, the track matches both cathodes
 544
 545   if (MuonChamber < 0) return kFALSE;
 546
 547   if (MuonChamber < 10) return TESTBIT(GetMUONClusterMap(), MuonChamber);
 548
 549   if (MuonChamber < 14) {
 550
 551     if (cathode < 0) return TESTBIT(GetHitsPatternInTrigCh(), 13-MuonChamber) &&
 552                             TESTBIT(GetHitsPatternInTrigCh(), 13-MuonChamber+4);
 553
 554     if (cathode < 2) return TESTBIT(GetHitsPatternInTrigCh(), 13-MuonChamber+(1-cathode)*4);
 555
 556   }
 557
 558   return kFALSE;
 559 }
 560
 561 //______________________________________________________________________________
 562 Bool_t AliAODTrack::MatchTriggerDigits() const
 563 {
 564   // return kTRUE if the track matches a digit on both planes of at least 2 trigger chambers
 565
 566   Int_t nMatchedChambers = 0;
 567   for (Int_t ich=10; ich<14; ich++) if (HitsMuonChamber(ich)) nMatchedChambers++;
 568
 569   return (nMatchedChambers >= 2);
 570 }
 571
 572 //______________________________________________________________________________
 573 Bool_t AliAODTrack::PropagateToDCA(const AliVVertex *vtx,
 574     Double_t b, Double_t maxd, Double_t dz[2], Double_t covar[3])
 575 {
 576   // compute impact parameters to the vertex vtx and their covariance matrix
 577   // b is the Bz, needed to propagate correctly the track to vertex
 578   // only the track parameters are update after the propagation (pos and mom),
 579   // not the covariance matrix. This is OK for propagation over short distance
 580   // inside the beam pipe.
 581   // return kFALSE is something went wrong
 582
 583   // convert to AliExternalTrackParam
 584   AliExternalTrackParam etp; etp.CopyFromVTrack(this);
 585
 586   Float_t xstart = etp.GetX();
 587   if(xstart>3.) {
 588     AliError("This method can be used only for propagation inside the beam pipe");
 589     return kFALSE;
 590   }
 591
 592   if(!etp.PropagateToDCA(vtx,b,maxd,dz,covar)) return kFALSE;
 593
 594   // update track position and momentum
 595   Double_t mom[3];
 596   etp.GetPxPyPz(mom);
 597   SetP(mom,kTRUE);
 598   etp.GetXYZ(mom);
 599   SetPosition(mom,kFALSE);
 600
 601
 602   return kTRUE;
 603 }
 604
 605 //______________________________________________________________________________
 606 Bool_t AliAODTrack::GetPxPyPz(Double_t p[3]) const
 607 {
 608     //---------------------------------------------------------------------
 609     // This function returns the global track momentum components
 610     //---------------------------------------------------------------------
 611   p[0]=Px(); p[1]=Py(); p[2]=Pz();
 612   return kTRUE;
 613 }
 614
 615
 616 //_______________________________________________________________________
 617 Float_t AliAODTrack::GetTPCClusterInfo(Int_t nNeighbours/*=3*/, Int_t type/*=0*/, Int_t row0, Int_t row1, Int_t bitType ) const
 618 {
 619   //
 620   // TPC cluster information
 621   // type 0: get fraction of found/findable clusters with neighbourhood definition
 622   //      1: findable clusters with neighbourhood definition
 623   //      2: found clusters
 624   // bitType:
 625   //      0 - all cluster used
 626   //      1 - clusters  used for the kalman update
 627   // definition of findable clusters:
 628   //            a cluster is defined as findable if there is another cluster
 629   //           within +- nNeighbours pad rows. The idea is to overcome threshold
 630   //           effects with a very simple algorithm.
 631   //
 632
 633
 634   Int_t found=0;
 635   Int_t findable=0;
 636   Int_t last=-nNeighbours;
 637   const TBits & clusterMap = (bitType%2==0) ? fTPCClusterMap : fTPCFitMap;
 638
 639   Int_t upperBound=clusterMap.GetNbits();
 640   if (upperBound>row1) upperBound=row1;
 641   for (Int_t i=row0; i<upperBound; ++i){
 642     //look to current row
 643     if (clusterMap[i]) {
 644       last=i;
 645       ++found;
 646       ++findable;
 647       continue;
 648     }
 649     //look to nNeighbours before
 650     if ((i-last)<=nNeighbours) {
 651       ++findable;
 652       continue;
 653     }
 654     //look to nNeighbours after
 655     for (Int_t j=i+1; j<i+1+nNeighbours; ++j){
 656       if (clusterMap[j]){
 657         ++findable;
 658         break;
 659       }
 660     }
 661   }
 662   if (type==2) return found;
 663   if (type==1) return findable;
 664
 665   if (type==0){
 666     Float_t fraction=0;
 667     if (findable>0)
 668       fraction=(Float_t)found/(Float_t)findable;
 669     else
 670       fraction=0;
 671     return fraction;
 672   }
 673   return 0;  // undefined type - default value
 674 }
 675
 676
 677 //______________________________________________________________________________
 678 Double_t  AliAODTrack::GetTRDslice(Int_t plane, Int_t slice) const {
 679   //
 680   // return TRD Pid information
 681   //
 682   if (!fDetPid) return -1;
 683   Double32_t *trdSlices=fDetPid->GetTRDsignal();
 684   if (!trdSlices) return -1;
 685   if ((plane<0) || (plane>=kTRDnPlanes)) {
 686     return -1.;
 687   }
 688
 689   Int_t ns=fDetPid->GetTRDnSlices()/kTRDnPlanes;
 690   if ((slice<-1) || (slice>=ns)) {
 691     return -1.;
 692   }
 693
 694   if(slice>=0) return trdSlices[plane*ns + slice];
 695
 696   // return average of the dEdx measurements
 697   Double_t q=0.; Double32_t *s = &trdSlices[plane*ns];
 698   for (Int_t i=0; i<ns; i++, s++) if((*s)>0.) q+=(*s);
 699   return q/ns;
 700 }
 701
 702 //______________________________________________________________________________
 703 UChar_t AliAODTrack::GetTRDntrackletsPID() const{
 704   //
 705   // return number of tracklets calculated from the slices
 706   //
 707   if(!fDetPid) return -1;
 708   return fDetPid->GetTRDntrackletsPID();
 709 }
 710
 711 //______________________________________________________________________________
 712 UChar_t AliAODTrack::GetTRDncls(Int_t layer) const {
 713   //
 714   // return number of TRD clusters
 715   //
 716   if(!fDetPid || layer > 5) return -1;
 717   if(layer < 0) return fDetPid->GetTRDncls();
 718   else return fDetPid->GetTRDncls(layer);
 719 }
 720
 721 //______________________________________________________________________________
 722 Double_t AliAODTrack::GetTRDmomentum(Int_t plane, Double_t */*sp*/) const
 723 {
 724   //Returns momentum estimation
 725   // in TRD layer "plane".
 726
 727   if (!fDetPid) return -1;
 728   const Double_t *trdMomentum=fDetPid->GetTRDmomentum();
 729
 730   if (!trdMomentum) {
 731     return -1.;
 732   }
 733   if ((plane<0) || (plane>=kTRDnPlanes)) {
 734     return -1.;
 735   }
 736
 737   return trdMomentum[plane];
 738 }
 739
 740 //_______________________________________________________________________
 741 Int_t AliAODTrack::GetTOFBunchCrossing(Double_t b, Bool_t) const
 742 {
 743   // Returns the number of bunch crossings after trigger (assuming 25ns spacing)
 744   const double kSpacing = 25e3; // min interbanch spacing
 745   const double kShift = 0;
 746   Int_t bcid = kTOFBCNA; // defualt one
 747   if (!IsOn(kTOFout) || !IsOn(kESDpid)) return bcid; // no info
 748   //
 749   double tdif = GetTOFsignal();
 750   if (IsOn(kTIME)) { // integrated time info is there
 751     int pid = (int)GetMostProbablePID();
 752     double ttimes[10];
 753     GetIntegratedTimes(ttimes);
 754     tdif -= ttimes[pid];
 755   }
 756   else { // assume integrated time info from TOF radius and momentum
 757     const double kRTOF = 385.;
 758     const double kCSpeed = 3.e-2; // cm/ps
 759     double p = P();
 760     if (p<0.001) p = 1.0;
 761     double m = M();
 762     double path =  kRTOF;     // mean TOF radius
 763     if (TMath::Abs(b)>kAlmost0) {  // account for curvature
 764       double curv = Pt()/(b*kB2C);
 765       if (curv>kAlmost0) {
 766         double tgl = Pz()/Pt();
 767         path = 2./curv*TMath::ASin(kRTOF*curv/2.)*TMath::Sqrt(1.+tgl*tgl);
 768       }
 769     }
 770     tdif -= path/kCSpeed*TMath::Sqrt(1.+m*m/(p*p));
 771   }
 772   bcid = TMath::Nint((tdif - kShift)/kSpacing);
 773   return bcid;
 774 }
 775
 776 //_____________________________________________________________________________
 777 Bool_t AliAODTrack::GetXYZAt(Double_t x, Double_t b, Double_t *r) const
 778 {
 779   //---------------------------------------------------------------------
 780   // This function returns the global track position extrapolated to
 781   // the radial position "x" (cm) in the magnetic field "b" (kG)
 782   //---------------------------------------------------------------------
 783
 784   //conversion of track parameter representation is
 785   //based on the implementation of AliExternalTrackParam::Set(...)
 786   //maybe some of this code can be moved to AliVTrack to avoid code duplication
 787   const double kSafe = 1e-5;
 788   Double_t alpha=0.0;
 789   Double_t radPos2 = fPosition[0]*fPosition[0]+fPosition[1]*fPosition[1];
 790   Double_t radMax  = 45.; // approximately ITS outer radius
 791   if (radPos2 < radMax*radMax) { // inside the ITS
 792      alpha = TMath::ATan2(fMomentum[1],fMomentum[0]);
 793   } else { // outside the ITS
 794      Float_t phiPos = TMath::Pi()+TMath::ATan2(-fPosition[1], -fPosition[0]);
 795      alpha =
 796      TMath::DegToRad()*(20*((((Int_t)(phiPos*TMath::RadToDeg()))/20))+10);
 797   }
 798   //
 799   Double_t cs=TMath::Cos(alpha), sn=TMath::Sin(alpha);
 800   // protection:  avoid alpha being too close to 0 or +-pi/2
 801   if (TMath::Abs(sn)<kSafe) {
 802     alpha = kSafe;
 803     cs=TMath::Cos(alpha);
 804     sn=TMath::Sin(alpha);
 805   }
 806   else if (cs<kSafe) {
 807     alpha -= TMath::Sign(kSafe, alpha);
 808     cs=TMath::Cos(alpha);
 809     sn=TMath::Sin(alpha);
 810   }
 811
 812   // Get the vertex of origin and the momentum
 813   TVector3 ver(fPosition[0],fPosition[1],fPosition[2]);
 814   TVector3 mom(fMomentum[0],fMomentum[1],fMomentum[2]);
 815   //
 816   // avoid momenta along axis
 817   if (TMath::Abs(mom[0])<kSafe) mom[0] = TMath::Sign(kSafe*TMath::Abs(mom[1]), mom[0]);
 818   if (TMath::Abs(mom[1])<kSafe) mom[1] = TMath::Sign(kSafe*TMath::Abs(mom[0]), mom[1]);
 819
 820   // Rotate to the local coordinate system
 821   ver.RotateZ(-alpha);
 822   mom.RotateZ(-alpha);
 823
 824   Double_t param0 = ver.Y();
 825   Double_t param1 = ver.Z();
 826   Double_t param2 = TMath::Sin(mom.Phi());
 827   Double_t param3 = mom.Pz()/mom.Pt();
 828   Double_t param4 = TMath::Sign(1/mom.Pt(),(Double_t)fCharge);
 829
 830   //calculate the propagated coordinates
 831   //this is based on AliExternalTrackParam::GetXYZAt(Double_t x, Double_t b, Double_t *r)
 832   Double_t dx=x-ver.X();
 833   if(TMath::Abs(dx)<=kAlmost0) return GetXYZ(r);
 834
 835   Double_t f1=param2;
 836   Double_t f2=f1 + dx*param4*b*kB2C;
 837
 838   if (TMath::Abs(f1) >= kAlmost1) return kFALSE;
 839   if (TMath::Abs(f2) >= kAlmost1) return kFALSE;
 840
 841   Double_t r1=TMath::Sqrt((1.-f1)*(1.+f1)), r2=TMath::Sqrt((1.-f2)*(1.+f2));
 842   r[0] = x;
 843   r[1] = param0 + dx*(f1+f2)/(r1+r2);
 844   r[2] = param1 + dx*(r2 + f2*(f1+f2)/(r1+r2))*param3;//Thanks to Andrea & Peter
 845
 846   return Local2GlobalPosition(r,alpha);
 847 }
 848
 849
 850