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),
51 // default constructor
54 SetPosition((Float_t*)NULL);
55 SetXYAtDCA(-999., -999.);
56 SetPxPyPzAtDCA(-999., -999., -999.);
57 SetPID((Float_t*)NULL);
60 //______________________________________________________________________________
61 AliAODTrack::AliAODTrack(Short_t id,
67 Double_t covMatrix[21],
71 AliAODVertex *prodVertex,
73 Bool_t usedForPrimVtxFit,
79 fChi2perNDF(chi2perNDF),
80 fChi2MatchTrigger(0.),
83 fITSMuonClusterMap(0),
84 fFilterMap(selectInfo),
93 fProdVertex(prodVertex)
98 SetPosition(x, isDCA);
99 SetXYAtDCA(-999., -999.);
100 SetPxPyPzAtDCA(-999., -999., -999.);
101 SetUsedForVtxFit(usedForVtxFit);
102 SetUsedForPrimVtxFit(usedForPrimVtxFit);
103 if(covMatrix) SetCovMatrix(covMatrix);
105 SetITSClusterMap(itsClusMap);
108 //______________________________________________________________________________
109 AliAODTrack::AliAODTrack(Short_t id,
115 Float_t covMatrix[21],
119 AliAODVertex *prodVertex,
120 Bool_t usedForVtxFit,
121 Bool_t usedForPrimVtxFit,
124 Float_t chi2perNDF) :
127 fChi2perNDF(chi2perNDF),
128 fChi2MatchTrigger(0.),
131 fITSMuonClusterMap(0),
132 fFilterMap(selectInfo),
141 fProdVertex(prodVertex)
146 SetPosition(x, isDCA);
147 SetXYAtDCA(-999., -999.);
148 SetPxPyPzAtDCA(-999., -999., -999.);
149 SetUsedForVtxFit(usedForVtxFit);
150 SetUsedForPrimVtxFit(usedForPrimVtxFit);
151 if(covMatrix) SetCovMatrix(covMatrix);
153 SetITSClusterMap(itsClusMap);
156 //______________________________________________________________________________
157 AliAODTrack::~AliAODTrack()
165 //______________________________________________________________________________
166 AliAODTrack::AliAODTrack(const AliAODTrack& trk) :
168 fRAtAbsorberEnd(trk.fRAtAbsorberEnd),
169 fChi2perNDF(trk.fChi2perNDF),
170 fChi2MatchTrigger(trk.fChi2MatchTrigger),
173 fITSMuonClusterMap(trk.fITSMuonClusterMap),
174 fFilterMap(trk.fFilterMap),
175 fTPCClusterMap(trk.fTPCClusterMap),
176 fTPCSharedMap(trk.fTPCSharedMap),
177 fTPCnclsF(trk.fTPCnclsF),
179 fCharge(trk.fCharge),
183 fProdVertex(trk.fProdVertex)
188 trk.GetPosition(fPosition);
189 SetXYAtDCA(trk.XAtDCA(), trk.YAtDCA());
190 SetPxPyPzAtDCA(trk.PxAtDCA(), trk.PyAtDCA(), trk.PzAtDCA());
191 SetUsedForVtxFit(trk.GetUsedForVtxFit());
192 SetUsedForPrimVtxFit(trk.GetUsedForPrimVtxFit());
193 if(trk.fCovMatrix) fCovMatrix=new AliAODRedCov<6>(*trk.fCovMatrix);
194 if(trk.fDetPid) fDetPid=new AliAODPid(*trk.fDetPid);
198 //______________________________________________________________________________
199 AliAODTrack& AliAODTrack::operator=(const AliAODTrack& trk)
201 // Assignment operator
204 AliVTrack::operator=(trk);
207 trk.GetPosition(fPosition);
210 SetXYAtDCA(trk.XAtDCA(), trk.YAtDCA());
211 SetPxPyPzAtDCA(trk.PxAtDCA(), trk.PyAtDCA(), trk.PzAtDCA());
213 fRAtAbsorberEnd = trk.fRAtAbsorberEnd;
215 fChi2perNDF = trk.fChi2perNDF;
216 fChi2MatchTrigger = trk.fChi2MatchTrigger;
221 fITSMuonClusterMap = trk.fITSMuonClusterMap;
222 fFilterMap = trk.fFilterMap;
226 fCharge = trk.fCharge;
230 if(trk.fCovMatrix) fCovMatrix=new AliAODRedCov<6>(*trk.fCovMatrix);
231 else fCovMatrix=NULL;
232 fProdVertex = trk.fProdVertex;
234 SetUsedForVtxFit(trk.GetUsedForVtxFit());
235 SetUsedForPrimVtxFit(trk.GetUsedForPrimVtxFit());
238 if(trk.fDetPid) fDetPid=new AliAODPid(*trk.fDetPid);
245 //______________________________________________________________________________
246 Double_t AliAODTrack::M(AODTrkPID_t pid) const
249 // Masses for nuclei don't exist in the PDG tables, therefore they were put by hand.
254 return 0.000510999; //TDatabasePDG::Instance()->GetParticle(11/*::kElectron*/)->Mass();
258 return 0.1056584; //TDatabasePDG::Instance()->GetParticle(13/*::kMuonMinus*/)->Mass();
262 return 0.13957; //TDatabasePDG::Instance()->GetParticle(211/*::kPiPlus*/)->Mass();
266 return 0.4937; //TDatabasePDG::Instance()->GetParticle(321/*::kKPlus*/)->Mass();
270 return 0.9382720; //TDatabasePDG::Instance()->GetParticle(2212/*::kProton*/)->Mass();
274 return 1.8756; //TDatabasePDG::Instance()->GetParticle(1000010020)->Mass();
278 return 2.8089; //TDatabasePDG::Instance()->GetParticle(1000010030)->Mass();
282 return 2.8084; //TDatabasePDG::Instance()->GetParticle(1000020030)->Mass();
286 return 3.7274; //TDatabasePDG::Instance()->GetParticle(1000020040)->Mass();
298 //______________________________________________________________________________
299 Double_t AliAODTrack::E(AODTrkPID_t pid) const
301 // Returns the energy of the particle of a given pid.
303 if (pid != kUnknown) { // particle was identified
305 return TMath::Sqrt(P()*P() + m*m);
306 } else { // pid unknown
311 //______________________________________________________________________________
312 Double_t AliAODTrack::Y(AODTrkPID_t pid) const
314 // Returns the rapidity of a particle of a given pid.
316 if (pid != kUnknown) { // particle was identified
319 if (e>=0 && e!=pz) { // energy was positive (e.g. not -999.) and not equal to pz
320 return 0.5*TMath::Log((e+pz)/(e-pz));
321 } else { // energy not known or equal to pz
324 } else { // pid unknown
329 //______________________________________________________________________________
330 Double_t AliAODTrack::Y(Double_t m) const
332 // Returns the rapidity of a particle of a given mass.
334 if (m >= 0.) { // mass makes sense
337 if (e>=0 && e!=pz) { // energy was positive (e.g. not -999.) and not equal to pz
338 return 0.5*TMath::Log((e+pz)/(e-pz));
339 } else { // energy not known or equal to pz
342 } else { // pid unknown
347 //______________________________________________________________________________
348 AliAODTrack::AODTrkPID_t AliAODTrack::GetMostProbablePID() const
350 // Returns the most probable PID array element.
353 AODTrkPID_t loc = kUnknown;
355 Bool_t allTheSame = kTRUE;
357 for (Int_t iPID = 0; iPID < nPID; iPID++) {
358 if (fPID[iPID] >= max) {
359 if (fPID[iPID] > max) {
362 loc = (AODTrkPID_t)iPID;
368 return allTheSame ? kUnknown : loc;
371 //______________________________________________________________________________
372 void AliAODTrack::ConvertAliPIDtoAODPID()
374 // Converts AliPID array.
375 // The numbering scheme is the same for electrons, muons, pions, kaons, and protons.
376 // Everything else has to be set to zero.
378 fPID[kDeuteron] = 0.;
388 //______________________________________________________________________________
389 template <class T> void AliAODTrack::SetP(const T *p, const Bool_t cartesian)
395 Double_t pt2 = p[0]*p[0] + p[1]*p[1];
396 Double_t pp = TMath::Sqrt(pt2 + p[2]*p[2]);
398 fMomentum[0] = TMath::Sqrt(pt2); // pt
399 fMomentum[1] = (pt2 != 0.) ? TMath::Pi()+TMath::ATan2(-p[1], -p[0]) : -999; // phi
400 fMomentum[2] = (pp != 0.) ? TMath::ACos(p[2] / pp) : -999.; // theta
402 fMomentum[0] = p[0]; // pt
403 fMomentum[1] = p[1]; // phi
404 fMomentum[2] = p[2]; // theta
407 fMomentum[0] = -999.;
408 fMomentum[1] = -999.;
409 fMomentum[2] = -999.;
413 //______________________________________________________________________________
414 template <class T> void AliAODTrack::SetPosition(const T *x, const Bool_t dca)
427 // don't know any better yet
428 fPosition[0] = -999.;
429 fPosition[1] = -999.;
430 fPosition[2] = -999.;
435 fPosition[0] = -999.;
436 fPosition[1] = -999.;
437 fPosition[2] = -999.;
441 //______________________________________________________________________________
442 void AliAODTrack::SetDCA(Double_t d, Double_t z)
451 //______________________________________________________________________________
452 void AliAODTrack::Print(Option_t* /* option */) const
454 // prints information about AliAODTrack
456 printf("Object name: %s Track type: %s\n", GetName(), GetTitle());
457 printf(" px = %f\n", Px());
458 printf(" py = %f\n", Py());
459 printf(" pz = %f\n", Pz());
460 printf(" pt = %f\n", Pt());
461 printf(" 1/pt = %f\n", OneOverPt());
462 printf(" theta = %f\n", Theta());
463 printf(" phi = %f\n", Phi());
464 printf(" chi2/NDF = %f\n", Chi2perNDF());
465 printf(" charge = %d\n", Charge());
468 //______________________________________________________________________________
469 void AliAODTrack::SetMatchTrigger(Int_t matchTrig)
471 // Set the MUON trigger information
473 case 0: // 0 track does not match trigger
474 fITSMuonClusterMap=fITSMuonClusterMap&0x3fffffff;
476 case 1: // 1 track match but does not pass pt cut
477 fITSMuonClusterMap=(fITSMuonClusterMap&0x3fffffff)|0x40000000;
479 case 2: // 2 track match Low pt cut
480 fITSMuonClusterMap=(fITSMuonClusterMap&0x3fffffff)|0x80000000;
482 case 3: // 3 track match High pt cut
483 fITSMuonClusterMap=fITSMuonClusterMap|0xc0000000;
486 fITSMuonClusterMap=fITSMuonClusterMap&0x3fffffff;
487 AliWarning(Form("unknown case for matchTrig: %d\n",matchTrig));
491 //______________________________________________________________________________
492 Bool_t AliAODTrack::HitsMuonChamber(Int_t MuonChamber, Int_t cathode) const
494 // return kTRUE if the track fires the given tracking or trigger chamber.
495 // If the chamber is a trigger one:
496 // - if cathode = 0 or 1, the track matches the corresponding cathode
497 // - if cathode = -1, the track matches both cathodes
499 if (MuonChamber < 0) return kFALSE;
501 if (MuonChamber < 10) return TESTBIT(GetMUONClusterMap(), MuonChamber);
503 if (MuonChamber < 14) {
505 if (cathode < 0) return TESTBIT(GetHitsPatternInTrigCh(), 13-MuonChamber) &&
506 TESTBIT(GetHitsPatternInTrigCh(), 13-MuonChamber+4);
508 if (cathode < 2) return TESTBIT(GetHitsPatternInTrigCh(), 13-MuonChamber+(1-cathode)*4);
515 //______________________________________________________________________________
516 Bool_t AliAODTrack::MatchTriggerDigits() const
518 // return kTRUE if the track matches a digit on both planes of at least 2 trigger chambers
520 Int_t nMatchedChambers = 0;
521 for (Int_t ich=10; ich<14; ich++) if (HitsMuonChamber(ich)) nMatchedChambers++;
523 return (nMatchedChambers >= 2);
526 //______________________________________________________________________________
527 Bool_t AliAODTrack::PropagateToDCA(const AliVVertex *vtx,
528 Double_t b, Double_t maxd, Double_t dz[2], Double_t covar[3])
530 // compute impact parameters to the vertex vtx and their covariance matrix
531 // b is the Bz, needed to propagate correctly the track to vertex
532 // only the track parameters are update after the propagation (pos and mom),
533 // not the covariance matrix. This is OK for propagation over short distance
534 // inside the beam pipe.
535 // return kFALSE is something went wrong
537 // convert to AliExternalTrackParam
538 AliExternalTrackParam etp(this);
540 Float_t xstart = etp.GetX();
542 AliError("This method can be used only for propagation inside the beam pipe");
546 if(!etp.PropagateToDCA(vtx,b,maxd,dz,covar)) return kFALSE;
548 // update track position and momentum
553 SetPosition(mom,kFALSE);
559 //______________________________________________________________________________
560 Bool_t AliAODTrack::GetPxPyPz(Double_t p[3]) const
562 //---------------------------------------------------------------------
563 // This function returns the global track momentum components
564 //---------------------------------------------------------------------
565 p[0]=Px(); p[1]=Py(); p[2]=Pz();
569 //______________________________________________________________________________
570 Float_t AliAODTrack::GetTPCClusterInfo(Int_t nNeighbours/*=3*/, Int_t type/*=0*/, Int_t row0, Int_t row1) const
573 // TPC cluster information
574 // type 0: get fraction of found/findable clusters with neighbourhood definition
575 // 1: findable clusters with neighbourhood definition
578 // definition of findable clusters:
579 // a cluster is defined as findable if there is another cluster
580 // within +- nNeighbours pad rows. The idea is to overcome threshold
581 // effects with a very simple algorithm.
584 if (type==2) return fTPCClusterMap.CountBits();
588 Int_t last=-nNeighbours;
590 for (Int_t i=row0; i<row1; ++i){
591 //look to current row
592 if (fTPCClusterMap[i]) {
598 //look to nNeighbours before
599 if ((i-last)<=nNeighbours) {
603 //look to nNeighbours after
604 for (Int_t j=i+1; j<i+1+nNeighbours; ++j){
605 if (fTPCClusterMap[j]){
611 if (type==1) return findable;
616 fraction=(Float_t)found/(Float_t)findable;
621 return 0; // undefined type - default value
624 //______________________________________________________________________________
625 Double_t AliAODTrack::GetTRDslice(Int_t plane, Int_t slice) const {
627 // return TRD Pid information
629 if (!fDetPid) return -1;
630 Double32_t *trdSlices=fDetPid->GetTRDsignal();
631 if (!trdSlices) return -1;
632 if ((plane<0) || (plane>=kTRDnPlanes)) {
636 Int_t ns=fDetPid->GetTRDnSlices();
637 if ((slice<-1) || (slice>=ns)) {
641 if(slice>=0) return trdSlices[plane*ns + slice];
643 // return average of the dEdx measurements
644 Double_t q=0.; Double32_t *s = &trdSlices[plane*ns];
645 for (Int_t i=0; i<ns; i++, s++) if((*s)>0.) q+=(*s);
649 //______________________________________________________________________________
650 Double_t AliAODTrack::GetTRDmomentum(Int_t plane, Double_t */*sp*/) const
652 //Returns momentum estimation
653 // in TRD layer "plane".
655 if (!fDetPid) return -1;
656 Float_t *trdMomentum=fDetPid->GetTRDmomentum();
661 if ((plane<0) || (plane>=kTRDnPlanes)) {
665 return trdMomentum[plane];
668 //_______________________________________________________________________
669 Int_t AliAODTrack::GetTOFBunchCrossing(Double_t b) const
671 // Returns the number of bunch crossings after trigger (assuming 25ns spacing)
672 const double kSpacing = 25e3; // min interbanch spacing
673 const double kShift = 0;
674 Int_t bcid = -1; // defualt one
675 if (!IsOn(kTOFout) || !IsOn(kESDpid)) return bcid; // no info
677 double tdif = GetTOFsignal();
678 if (IsOn(kTIME)) { // integrated time info is there
679 int pid = (int)GetMostProbablePID();
681 GetIntegratedTimes(ttimes);
684 else { // assume integrated time info from TOF radius and momentum
685 const double kRTOF = 385.;
686 const double kCSpeed = 3.e-2; // cm/ps
688 if (p<0.001) p = 1.0;
690 double path = kRTOF; // mean TOF radius
691 if (TMath::Abs(b)>kAlmost0) { // account for curvature
692 double curv = Pt()/(b*kB2C);
694 double tgl = Pz()/Pt();
695 path = 2./curv*TMath::ASin(kRTOF*curv/2.)*TMath::Sqrt(1.+tgl*tgl);
698 tdif -= path/kCSpeed*TMath::Sqrt(1.+m*m/(p*p));
700 bcid = TMath::Nint((tdif - kShift)/kSpacing);