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