1 /**************************************************************************
2 * Copyright(c) 1998-2007, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
18 //-------------------------------------------------------------------------
19 // AOD track implementation of AliVTrack
20 // Author: Markus Oldenburg, CERN
21 // Markus.Oldenburg@cern.ch
22 //-------------------------------------------------------------------------
25 #include "AliExternalTrackParam.h"
26 #include "AliVVertex.h"
27 #include "AliAODTrack.h"
31 //______________________________________________________________________________
32 AliAODTrack::AliAODTrack() :
36 fChi2MatchTrigger(0.),
39 fITSMuonClusterMap(0),
48 fCaloIndex(kEMCALNoMatch),
52 fTrackPhiOnEMCal(-999),
53 fTrackEtaOnEMCal(-999)
55 // default constructor
58 SetPosition((Float_t*)NULL);
59 SetXYAtDCA(-999., -999.);
60 SetPxPyPzAtDCA(-999., -999., -999.);
61 SetPID((Float_t*)NULL);
64 //______________________________________________________________________________
65 AliAODTrack::AliAODTrack(Short_t id,
71 Double_t covMatrix[21],
75 AliAODVertex *prodVertex,
77 Bool_t usedForPrimVtxFit,
83 fChi2perNDF(chi2perNDF),
84 fChi2MatchTrigger(0.),
87 fITSMuonClusterMap(0),
88 fFilterMap(selectInfo),
96 fCaloIndex(kEMCALNoMatch),
99 fProdVertex(prodVertex),
100 fTrackPhiOnEMCal(-999),
101 fTrackEtaOnEMCal(-999)
106 SetPosition(x, isDCA);
107 SetXYAtDCA(-999., -999.);
108 SetPxPyPzAtDCA(-999., -999., -999.);
109 SetUsedForVtxFit(usedForVtxFit);
110 SetUsedForPrimVtxFit(usedForPrimVtxFit);
111 if(covMatrix) SetCovMatrix(covMatrix);
113 SetITSClusterMap(itsClusMap);
116 //______________________________________________________________________________
117 AliAODTrack::AliAODTrack(Short_t id,
123 Float_t covMatrix[21],
127 AliAODVertex *prodVertex,
128 Bool_t usedForVtxFit,
129 Bool_t usedForPrimVtxFit,
132 Float_t chi2perNDF) :
135 fChi2perNDF(chi2perNDF),
136 fChi2MatchTrigger(0.),
139 fITSMuonClusterMap(0),
140 fFilterMap(selectInfo),
148 fCaloIndex(kEMCALNoMatch),
151 fProdVertex(prodVertex),
152 fTrackPhiOnEMCal(-999),
153 fTrackEtaOnEMCal(-999)
158 SetPosition(x, isDCA);
159 SetXYAtDCA(-999., -999.);
160 SetPxPyPzAtDCA(-999., -999., -999.);
161 SetUsedForVtxFit(usedForVtxFit);
162 SetUsedForPrimVtxFit(usedForPrimVtxFit);
163 if(covMatrix) SetCovMatrix(covMatrix);
165 SetITSClusterMap(itsClusMap);
168 //______________________________________________________________________________
169 AliAODTrack::~AliAODTrack()
177 //______________________________________________________________________________
178 AliAODTrack::AliAODTrack(const AliAODTrack& trk) :
180 fRAtAbsorberEnd(trk.fRAtAbsorberEnd),
181 fChi2perNDF(trk.fChi2perNDF),
182 fChi2MatchTrigger(trk.fChi2MatchTrigger),
185 fITSMuonClusterMap(trk.fITSMuonClusterMap),
186 fFilterMap(trk.fFilterMap),
187 fTPCFitMap(trk.fTPCFitMap),
188 fTPCClusterMap(trk.fTPCClusterMap),
189 fTPCSharedMap(trk.fTPCSharedMap),
190 fTPCnclsF(trk.fTPCnclsF),
192 fCharge(trk.fCharge),
194 fCaloIndex(trk.fCaloIndex),
197 fProdVertex(trk.fProdVertex),
198 fTrackPhiOnEMCal(trk.fTrackPhiOnEMCal),
199 fTrackEtaOnEMCal(trk.fTrackEtaOnEMCal)
204 trk.GetPosition(fPosition);
205 SetXYAtDCA(trk.XAtDCA(), trk.YAtDCA());
206 SetPxPyPzAtDCA(trk.PxAtDCA(), trk.PyAtDCA(), trk.PzAtDCA());
207 SetUsedForVtxFit(trk.GetUsedForVtxFit());
208 SetUsedForPrimVtxFit(trk.GetUsedForPrimVtxFit());
209 if(trk.fCovMatrix) fCovMatrix=new AliAODRedCov<6>(*trk.fCovMatrix);
210 if(trk.fDetPid) fDetPid=new AliAODPid(*trk.fDetPid);
214 //______________________________________________________________________________
215 AliAODTrack& AliAODTrack::operator=(const AliAODTrack& trk)
217 // Assignment operator
220 AliVTrack::operator=(trk);
223 trk.GetPosition(fPosition);
224 SetXYAtDCA(trk.XAtDCA(), trk.YAtDCA());
225 SetPxPyPzAtDCA(trk.PxAtDCA(), trk.PyAtDCA(), trk.PzAtDCA());
226 fRAtAbsorberEnd = trk.fRAtAbsorberEnd;
227 fChi2perNDF = trk.fChi2perNDF;
228 fChi2MatchTrigger = trk.fChi2MatchTrigger;
232 fITSMuonClusterMap = trk.fITSMuonClusterMap;
233 fFilterMap = trk.fFilterMap;
234 fTPCFitMap = trk.fTPCFitMap;
235 fTPCClusterMap = trk.fTPCClusterMap;
236 fTPCSharedMap = trk.fTPCSharedMap;
237 fTPCnclsF = trk.fTPCnclsF;
239 fCharge = trk.fCharge;
241 fCaloIndex = trk.fCaloIndex;
242 fTrackPhiOnEMCal = trk.fTrackPhiOnEMCal;
243 fTrackEtaOnEMCal = trk.fTrackEtaOnEMCal;
246 if(trk.fCovMatrix) fCovMatrix=new AliAODRedCov<6>(*trk.fCovMatrix);
247 else fCovMatrix=NULL;
250 fProdVertex = trk.fProdVertex;
251 SetUsedForVtxFit(trk.GetUsedForVtxFit());
252 SetUsedForPrimVtxFit(trk.GetUsedForPrimVtxFit());
255 if(trk.fDetPid) fDetPid=new AliAODPid(*trk.fDetPid);
264 //______________________________________________________________________________
265 Double_t AliAODTrack::M(AODTrkPID_t pid) const
268 // Masses for nuclei don't exist in the PDG tables, therefore they were put by hand.
273 return 0.000510999; //TDatabasePDG::Instance()->GetParticle(11/*::kElectron*/)->Mass();
277 return 0.1056584; //TDatabasePDG::Instance()->GetParticle(13/*::kMuonMinus*/)->Mass();
281 return 0.13957; //TDatabasePDG::Instance()->GetParticle(211/*::kPiPlus*/)->Mass();
285 return 0.4937; //TDatabasePDG::Instance()->GetParticle(321/*::kKPlus*/)->Mass();
289 return 0.9382720; //TDatabasePDG::Instance()->GetParticle(2212/*::kProton*/)->Mass();
293 return 1.8756; //TDatabasePDG::Instance()->GetParticle(1000010020)->Mass();
297 return 2.8089; //TDatabasePDG::Instance()->GetParticle(1000010030)->Mass();
301 return 2.8084; //TDatabasePDG::Instance()->GetParticle(1000020030)->Mass();
305 return 3.7274; //TDatabasePDG::Instance()->GetParticle(1000020040)->Mass();
317 //______________________________________________________________________________
318 Double_t AliAODTrack::E(AODTrkPID_t pid) const
320 // Returns the energy of the particle of a given pid.
322 if (pid != kUnknown) { // particle was identified
324 return TMath::Sqrt(P()*P() + m*m);
325 } else { // pid unknown
330 //______________________________________________________________________________
331 Double_t AliAODTrack::Y(AODTrkPID_t pid) const
333 // Returns the rapidity of a particle of a given pid.
335 if (pid != kUnknown) { // particle was identified
338 if (e>=0 && e!=pz) { // energy was positive (e.g. not -999.) and not equal to pz
339 return 0.5*TMath::Log((e+pz)/(e-pz));
340 } else { // energy not known or equal to pz
343 } else { // pid unknown
348 //______________________________________________________________________________
349 Double_t AliAODTrack::Y(Double_t m) const
351 // Returns the rapidity of a particle of a given mass.
353 if (m >= 0.) { // mass makes sense
356 if (e>=0 && e!=pz) { // energy was positive (e.g. not -999.) and not equal to pz
357 return 0.5*TMath::Log((e+pz)/(e-pz));
358 } else { // energy not known or equal to pz
361 } else { // pid unknown
366 //______________________________________________________________________________
367 AliAODTrack::AODTrkPID_t AliAODTrack::GetMostProbablePID() const
369 // Returns the most probable PID array element.
372 AODTrkPID_t loc = kUnknown;
374 Bool_t allTheSame = kTRUE;
376 for (Int_t iPID = 0; iPID < nPID; iPID++) {
377 if (fPID[iPID] >= max) {
378 if (fPID[iPID] > max) {
381 loc = (AODTrkPID_t)iPID;
387 return allTheSame ? kUnknown : loc;
390 //______________________________________________________________________________
391 void AliAODTrack::ConvertAliPIDtoAODPID()
393 // Converts AliPID array.
394 // The numbering scheme is the same for electrons, muons, pions, kaons, and protons.
395 // Everything else has to be set to zero.
397 fPID[kDeuteron] = 0.;
407 //______________________________________________________________________________
408 template <class T> void AliAODTrack::SetP(const T *p, const Bool_t cartesian)
414 Double_t pt2 = p[0]*p[0] + p[1]*p[1];
415 Double_t pp = TMath::Sqrt(pt2 + p[2]*p[2]);
417 fMomentum[0] = TMath::Sqrt(pt2); // pt
418 fMomentum[1] = (pt2 != 0.) ? TMath::Pi()+TMath::ATan2(-p[1], -p[0]) : -999; // phi
419 fMomentum[2] = (pp != 0.) ? TMath::ACos(p[2] / pp) : -999.; // theta
421 fMomentum[0] = p[0]; // pt
422 fMomentum[1] = p[1]; // phi
423 fMomentum[2] = p[2]; // theta
426 fMomentum[0] = -999.;
427 fMomentum[1] = -999.;
428 fMomentum[2] = -999.;
432 //______________________________________________________________________________
433 template <class T> void AliAODTrack::SetPosition(const T *x, const Bool_t dca)
446 // don't know any better yet
447 fPosition[0] = -999.;
448 fPosition[1] = -999.;
449 fPosition[2] = -999.;
454 fPosition[0] = -999.;
455 fPosition[1] = -999.;
456 fPosition[2] = -999.;
460 //______________________________________________________________________________
461 void AliAODTrack::SetDCA(Double_t d, Double_t z)
470 //______________________________________________________________________________
471 void AliAODTrack::Print(Option_t* /* option */) const
473 // prints information about AliAODTrack
475 printf("Object name: %s Track type: %s\n", GetName(), GetTitle());
476 printf(" px = %f\n", Px());
477 printf(" py = %f\n", Py());
478 printf(" pz = %f\n", Pz());
479 printf(" pt = %f\n", Pt());
480 printf(" 1/pt = %f\n", OneOverPt());
481 printf(" theta = %f\n", Theta());
482 printf(" phi = %f\n", Phi());
483 printf(" chi2/NDF = %f\n", Chi2perNDF());
484 printf(" charge = %d\n", Charge());
487 //______________________________________________________________________________
488 void AliAODTrack::SetMatchTrigger(Int_t matchTrig)
490 // Set the MUON trigger information
492 case 0: // 0 track does not match trigger
493 fITSMuonClusterMap=fITSMuonClusterMap&0x3fffffff;
495 case 1: // 1 track match but does not pass pt cut
496 fITSMuonClusterMap=(fITSMuonClusterMap&0x3fffffff)|0x40000000;
498 case 2: // 2 track match Low pt cut
499 fITSMuonClusterMap=(fITSMuonClusterMap&0x3fffffff)|0x80000000;
501 case 3: // 3 track match High pt cut
502 fITSMuonClusterMap=fITSMuonClusterMap|0xc0000000;
505 fITSMuonClusterMap=fITSMuonClusterMap&0x3fffffff;
506 AliWarning(Form("unknown case for matchTrig: %d\n",matchTrig));
510 //______________________________________________________________________________
511 Bool_t AliAODTrack::HitsMuonChamber(Int_t MuonChamber, Int_t cathode) const
513 // return kTRUE if the track fires the given tracking or trigger chamber.
514 // If the chamber is a trigger one:
515 // - if cathode = 0 or 1, the track matches the corresponding cathode
516 // - if cathode = -1, the track matches both cathodes
518 if (MuonChamber < 0) return kFALSE;
520 if (MuonChamber < 10) return TESTBIT(GetMUONClusterMap(), MuonChamber);
522 if (MuonChamber < 14) {
524 if (cathode < 0) return TESTBIT(GetHitsPatternInTrigCh(), 13-MuonChamber) &&
525 TESTBIT(GetHitsPatternInTrigCh(), 13-MuonChamber+4);
527 if (cathode < 2) return TESTBIT(GetHitsPatternInTrigCh(), 13-MuonChamber+(1-cathode)*4);
534 //______________________________________________________________________________
535 Bool_t AliAODTrack::MatchTriggerDigits() const
537 // return kTRUE if the track matches a digit on both planes of at least 2 trigger chambers
539 Int_t nMatchedChambers = 0;
540 for (Int_t ich=10; ich<14; ich++) if (HitsMuonChamber(ich)) nMatchedChambers++;
542 return (nMatchedChambers >= 2);
545 //______________________________________________________________________________
546 Bool_t AliAODTrack::PropagateToDCA(const AliVVertex *vtx,
547 Double_t b, Double_t maxd, Double_t dz[2], Double_t covar[3])
549 // compute impact parameters to the vertex vtx and their covariance matrix
550 // b is the Bz, needed to propagate correctly the track to vertex
551 // only the track parameters are update after the propagation (pos and mom),
552 // not the covariance matrix. This is OK for propagation over short distance
553 // inside the beam pipe.
554 // return kFALSE is something went wrong
556 // convert to AliExternalTrackParam
557 AliExternalTrackParam etp; etp.CopyFromVTrack(this);
559 Float_t xstart = etp.GetX();
561 AliError("This method can be used only for propagation inside the beam pipe");
565 if(!etp.PropagateToDCA(vtx,b,maxd,dz,covar)) return kFALSE;
567 // update track position and momentum
572 SetPosition(mom,kFALSE);
578 //______________________________________________________________________________
579 Bool_t AliAODTrack::GetPxPyPz(Double_t p[3]) const
581 //---------------------------------------------------------------------
582 // This function returns the global track momentum components
583 //---------------------------------------------------------------------
584 p[0]=Px(); p[1]=Py(); p[2]=Pz();
589 //_______________________________________________________________________
590 Float_t AliAODTrack::GetTPCClusterInfo(Int_t nNeighbours/*=3*/, Int_t type/*=0*/, Int_t row0, Int_t row1, Int_t bitType ) const
593 // TPC cluster information
594 // type 0: get fraction of found/findable clusters with neighbourhood definition
595 // 1: findable clusters with neighbourhood definition
598 // 0 - all cluster used
599 // 1 - clusters used for the kalman update
600 // definition of findable clusters:
601 // a cluster is defined as findable if there is another cluster
602 // within +- nNeighbours pad rows. The idea is to overcome threshold
603 // effects with a very simple algorithm.
609 Int_t last=-nNeighbours;
610 const TBits & clusterMap = (bitType%2==0) ? fTPCClusterMap : fTPCFitMap;
612 Int_t upperBound=clusterMap.GetNbits();
613 if (upperBound>row1) upperBound=row1;
614 for (Int_t i=row0; i<upperBound; ++i){
615 //look to current row
622 //look to nNeighbours before
623 if ((i-last)<=nNeighbours) {
627 //look to nNeighbours after
628 for (Int_t j=i+1; j<i+1+nNeighbours; ++j){
635 if (type==2) return found;
636 if (type==1) return findable;
641 fraction=(Float_t)found/(Float_t)findable;
646 return 0; // undefined type - default value
650 //______________________________________________________________________________
651 Double_t AliAODTrack::GetTRDslice(Int_t plane, Int_t slice) const {
653 // return TRD Pid information
655 if (!fDetPid) return -1;
656 Double32_t *trdSlices=fDetPid->GetTRDsignal();
657 if (!trdSlices) return -1;
658 if ((plane<0) || (plane>=kTRDnPlanes)) {
662 Int_t ns=fDetPid->GetTRDnSlices()/kTRDnPlanes;
663 if ((slice<-1) || (slice>=ns)) {
667 if(slice>=0) return trdSlices[plane*ns + slice];
669 // return average of the dEdx measurements
670 Double_t q=0.; Double32_t *s = &trdSlices[plane*ns];
671 for (Int_t i=0; i<ns; i++, s++) if((*s)>0.) q+=(*s);
675 //______________________________________________________________________________
676 UChar_t AliAODTrack::GetTRDntrackletsPID() const{
678 // return number of tracklets calculated from the slices
680 if(!fDetPid) return -1;
681 return fDetPid->GetTRDntrackletsPID();
684 //______________________________________________________________________________
685 UChar_t AliAODTrack::GetTRDncls(Int_t layer) const {
687 // return number of TRD clusters
689 if(!fDetPid || layer > 5) return -1;
690 if(layer < 0) return fDetPid->GetTRDncls();
691 else return fDetPid->GetTRDncls(layer);
694 //______________________________________________________________________________
695 Double_t AliAODTrack::GetTRDmomentum(Int_t plane, Double_t */*sp*/) const
697 //Returns momentum estimation
698 // in TRD layer "plane".
700 if (!fDetPid) return -1;
701 const Float_t *trdMomentum=fDetPid->GetTRDmomentum();
706 if ((plane<0) || (plane>=kTRDnPlanes)) {
710 return trdMomentum[plane];
713 //_______________________________________________________________________
714 Int_t AliAODTrack::GetTOFBunchCrossing(Double_t b, Bool_t) const
716 // Returns the number of bunch crossings after trigger (assuming 25ns spacing)
717 const double kSpacing = 25e3; // min interbanch spacing
718 const double kShift = 0;
719 Int_t bcid = kTOFBCNA; // defualt one
720 if (!IsOn(kTOFout) || !IsOn(kESDpid)) return bcid; // no info
722 double tdif = GetTOFsignal();
723 if (IsOn(kTIME)) { // integrated time info is there
724 int pid = (int)GetMostProbablePID();
726 GetIntegratedTimes(ttimes);
729 else { // assume integrated time info from TOF radius and momentum
730 const double kRTOF = 385.;
731 const double kCSpeed = 3.e-2; // cm/ps
733 if (p<0.001) p = 1.0;
735 double path = kRTOF; // mean TOF radius
736 if (TMath::Abs(b)>kAlmost0) { // account for curvature
737 double curv = Pt()/(b*kB2C);
739 double tgl = Pz()/Pt();
740 path = 2./curv*TMath::ASin(kRTOF*curv/2.)*TMath::Sqrt(1.+tgl*tgl);
743 tdif -= path/kCSpeed*TMath::Sqrt(1.+m*m/(p*p));
745 bcid = TMath::Nint((tdif - kShift)/kSpacing);