// Markus.Oldenburg@cern.ch
//-------------------------------------------------------------------------
+#include <TVector3.h>
#include "AliLog.h"
#include "AliExternalTrackParam.h"
#include "AliVVertex.h"
+#include "AliDetectorPID.h"
+#include "AliAODEvent.h"
+#include "AliAODHMPIDrings.h"
+
#include "AliAODTrack.h"
ClassImp(AliAODTrack)
fFlags(0),
fLabel(-999),
fITSMuonClusterMap(0),
+ fMUONtrigHitsMapTrg(0),
+ fMUONtrigHitsMapTrk(0),
fFilterMap(0),
+ fTPCFitMap(),
fTPCClusterMap(),
fTPCSharedMap(),
fTPCnclsF(0),
+ fTPCNCrossedRows(0),
fID(-999),
fCharge(-99),
fType(kUndef),
fCaloIndex(kEMCALNoMatch),
fCovMatrix(NULL),
fDetPid(NULL),
- fProdVertex(NULL)
+ fDetectorPID(NULL),
+ fProdVertex(NULL),
+ fTrackPhiOnEMCal(-999),
+ fTrackEtaOnEMCal(-999),
+ fTPCsignalTuned(0),
+ fAODEvent(NULL)
{
// default constructor
fFlags(0),
fLabel(label),
fITSMuonClusterMap(0),
+ fMUONtrigHitsMapTrg(0),
+ fMUONtrigHitsMapTrk(0),
fFilterMap(selectInfo),
+ fTPCFitMap(),
fTPCClusterMap(),
fTPCSharedMap(),
fTPCnclsF(0),
+ fTPCNCrossedRows(0),
fID(id),
fCharge(charge),
fType(ttype),
fCaloIndex(kEMCALNoMatch),
fCovMatrix(NULL),
fDetPid(NULL),
- fProdVertex(prodVertex)
+ fDetectorPID(NULL),
+ fProdVertex(prodVertex),
+ fTrackPhiOnEMCal(-999),
+ fTrackEtaOnEMCal(-999),
+ fTPCsignalTuned(0),
+ fAODEvent(NULL)
{
// constructor
fFlags(0),
fLabel(label),
fITSMuonClusterMap(0),
+ fMUONtrigHitsMapTrg(0),
+ fMUONtrigHitsMapTrk(0),
fFilterMap(selectInfo),
+ fTPCFitMap(),
fTPCClusterMap(),
fTPCSharedMap(),
fTPCnclsF(0),
+ fTPCNCrossedRows(0),
fID(id),
fCharge(charge),
fType(ttype),
fCaloIndex(kEMCALNoMatch),
fCovMatrix(NULL),
fDetPid(NULL),
- fProdVertex(prodVertex)
+ fDetectorPID(NULL),
+ fProdVertex(prodVertex),
+ fTrackPhiOnEMCal(-999),
+ fTrackEtaOnEMCal(-999),
+ fTPCsignalTuned(0),
+ fAODEvent(NULL)
{
// constructor
// destructor
delete fCovMatrix;
delete fDetPid;
+ delete fDetectorPID;
}
fFlags(trk.fFlags),
fLabel(trk.fLabel),
fITSMuonClusterMap(trk.fITSMuonClusterMap),
+ fMUONtrigHitsMapTrg(trk.fMUONtrigHitsMapTrg),
+ fMUONtrigHitsMapTrk(trk.fMUONtrigHitsMapTrk),
fFilterMap(trk.fFilterMap),
+ fTPCFitMap(trk.fTPCFitMap),
fTPCClusterMap(trk.fTPCClusterMap),
fTPCSharedMap(trk.fTPCSharedMap),
fTPCnclsF(trk.fTPCnclsF),
+ fTPCNCrossedRows(trk.fTPCNCrossedRows),
fID(trk.fID),
fCharge(trk.fCharge),
fType(trk.fType),
fCaloIndex(trk.fCaloIndex),
fCovMatrix(NULL),
fDetPid(NULL),
- fProdVertex(trk.fProdVertex)
+ fDetectorPID(NULL),
+ fProdVertex(trk.fProdVertex),
+ fTrackPhiOnEMCal(trk.fTrackPhiOnEMCal),
+ fTrackEtaOnEMCal(trk.fTrackEtaOnEMCal),
+ fTPCsignalTuned(trk.fTPCsignalTuned),
+ fAODEvent(trk.fAODEvent)
{
// Copy constructor
if(trk.fCovMatrix) fCovMatrix=new AliAODRedCov<6>(*trk.fCovMatrix);
if(trk.fDetPid) fDetPid=new AliAODPid(*trk.fDetPid);
SetPID(trk.fPID);
+ if (trk.fDetectorPID) fDetectorPID = new AliDetectorPID(*trk.fDetectorPID);
}
//______________________________________________________________________________
trk.GetP(fMomentum);
trk.GetPosition(fPosition);
- trk.GetPID(fPID);
-
SetXYAtDCA(trk.XAtDCA(), trk.YAtDCA());
SetPxPyPzAtDCA(trk.PxAtDCA(), trk.PyAtDCA(), trk.PzAtDCA());
-
- fRAtAbsorberEnd = trk.fRAtAbsorberEnd;
-
- fChi2perNDF = trk.fChi2perNDF;
- fChi2MatchTrigger = trk.fChi2MatchTrigger;
-
- fFlags = trk.fFlags;
- fLabel = trk.fLabel;
-
+ fRAtAbsorberEnd = trk.fRAtAbsorberEnd;
+ fChi2perNDF = trk.fChi2perNDF;
+ fChi2MatchTrigger = trk.fChi2MatchTrigger;
+ trk.GetPID(fPID);
+ fFlags = trk.fFlags;
+ fLabel = trk.fLabel;
fITSMuonClusterMap = trk.fITSMuonClusterMap;
- fFilterMap = trk.fFilterMap;
-
- fID = trk.fID;
-
- fCharge = trk.fCharge;
- fType = trk.fType;
-
- fCaloIndex = trk.fCaloIndex;
+ fMUONtrigHitsMapTrg = trk.fMUONtrigHitsMapTrg;
+ fMUONtrigHitsMapTrk = trk.fMUONtrigHitsMapTrk;
+ fFilterMap = trk.fFilterMap;
+ fTPCFitMap = trk.fTPCFitMap;
+ fTPCClusterMap = trk.fTPCClusterMap;
+ fTPCSharedMap = trk.fTPCSharedMap;
+ fTPCnclsF = trk.fTPCnclsF;
+ fTPCNCrossedRows = trk.fTPCNCrossedRows;
+ fID = trk.fID;
+ fCharge = trk.fCharge;
+ fType = trk.fType;
+ fCaloIndex = trk.fCaloIndex;
+ fTrackPhiOnEMCal = trk.fTrackPhiOnEMCal;
+ fTrackEtaOnEMCal = trk.fTrackEtaOnEMCal;
+ fTPCsignalTuned = trk.fTPCsignalTuned;
delete fCovMatrix;
if(trk.fCovMatrix) fCovMatrix=new AliAODRedCov<6>(*trk.fCovMatrix);
else fCovMatrix=NULL;
- fProdVertex = trk.fProdVertex;
+
+ fProdVertex = trk.fProdVertex;
SetUsedForVtxFit(trk.GetUsedForVtxFit());
SetUsedForPrimVtxFit(trk.GetUsedForPrimVtxFit());
+ //detector raw signals
delete fDetPid;
if(trk.fDetPid) fDetPid=new AliAODPid(*trk.fDetPid);
else fDetPid=NULL;
+
+ //calibrated PID cache
+ delete fDetectorPID;
+ fDetectorPID=0x0;
+ if (trk.fDetectorPID) fDetectorPID = new AliDetectorPID(*trk.fDetectorPID);
}
return *this;
//______________________________________________________________________________
-template <class T> void AliAODTrack::SetP(const T *p, const Bool_t cartesian)
+template <typename T> void AliAODTrack::SetP(const T *p, const Bool_t cartesian)
{
// Set the momentum
}
}
+/*
//______________________________________________________________________________
-template <class T> void AliAODTrack::SetPosition(const T *x, const Bool_t dca)
+template <typename T> void AliAODTrack::SetPosition(const T *x, const Bool_t dca)
{
// set the position
fPosition[2] = -999.;
}
}
-
+*/
//______________________________________________________________________________
void AliAODTrack::SetDCA(Double_t d, Double_t z)
{
// return kFALSE is something went wrong
// convert to AliExternalTrackParam
- AliExternalTrackParam etp(this);
+ AliExternalTrackParam etp; etp.CopyFromVTrack(this);
Float_t xstart = etp.GetX();
if(xstart>3.) {
return kTRUE;
}
-//______________________________________________________________________________
-Float_t AliAODTrack::GetTPCClusterInfo(Int_t nNeighbours/*=3*/, Int_t type/*=0*/, Int_t row0, Int_t row1) const
+
+//_______________________________________________________________________
+Float_t AliAODTrack::GetTPCClusterInfo(Int_t nNeighbours/*=3*/, Int_t type/*=0*/, Int_t row0, Int_t row1, Int_t bitType ) const
{
//
- // TPC cluster information
+ // TPC cluster information
// type 0: get fraction of found/findable clusters with neighbourhood definition
// 1: findable clusters with neighbourhood definition
// 2: found clusters
- //
+ // bitType:
+ // 0 - all cluster used
+ // 1 - clusters used for the kalman update
// definition of findable clusters:
// a cluster is defined as findable if there is another cluster
// within +- nNeighbours pad rows. The idea is to overcome threshold
// effects with a very simple algorithm.
//
-
- if (type==2) return fTPCClusterMap.CountBits();
+
Int_t found=0;
Int_t findable=0;
Int_t last=-nNeighbours;
+ const TBits & clusterMap = (bitType%2==0) ? fTPCClusterMap : fTPCFitMap;
- for (Int_t i=row0; i<row1; ++i){
+ Int_t upperBound=clusterMap.GetNbits();
+ if (upperBound>row1) upperBound=row1;
+ for (Int_t i=row0; i<upperBound; ++i){
//look to current row
- if (fTPCClusterMap[i]) {
+ if (clusterMap[i]) {
last=i;
++found;
++findable;
}
//look to nNeighbours after
for (Int_t j=i+1; j<i+1+nNeighbours; ++j){
- if (fTPCClusterMap[j]){
+ if (clusterMap[j]){
++findable;
break;
}
}
}
+ if (type==2) return found;
if (type==1) return findable;
if (type==0){
Float_t fraction=0;
- if (findable>0)
+ if (findable>0)
fraction=(Float_t)found/(Float_t)findable;
- else
+ else
fraction=0;
return fraction;
- }
+ }
return 0; // undefined type - default value
}
+
//______________________________________________________________________________
Double_t AliAODTrack::GetTRDslice(Int_t plane, Int_t slice) const {
//
// return TRD Pid information
//
if (!fDetPid) return -1;
- Double32_t *trdSlices=fDetPid->GetTRDsignal();
+ Double32_t *trdSlices=fDetPid->GetTRDslices();
if (!trdSlices) return -1;
if ((plane<0) || (plane>=kTRDnPlanes)) {
return -1.;
// return number of tracklets calculated from the slices
//
if(!fDetPid) return -1;
+ return fDetPid->GetTRDntrackletsPID();
+}
- Int_t ntracklets = 0, // Number of tracklets / track
- nSlicesTracklet = fDetPid->GetTRDnSlices()/kTRDnPlanes, // Number of slices per tracklet
- nSlicesNonZero = 0; // Number of slices containing a dE/dx measurement
- for(Int_t ily = 0; ily < kTRDnPlanes; ily++){
- // a tracklet is found if it has at least one slice containing a dE/dx measurement
- nSlicesNonZero = 0;
- for(Int_t islice = 0; islice < nSlicesTracklet; islice++){
- if(fDetPid->GetTRDsignal()[nSlicesTracklet * ily + islice] > 0.01) nSlicesNonZero++;
- }
- if(nSlicesNonZero) ntracklets++;
- }
- return ntracklets;
+//______________________________________________________________________________
+UChar_t AliAODTrack::GetTRDncls(Int_t layer) const {
+ //
+ // return number of TRD clusters
+ //
+ if(!fDetPid || layer > 5) return -1;
+ if(layer < 0) return fDetPid->GetTRDncls();
+ else return fDetPid->GetTRDncls(layer);
}
//______________________________________________________________________________
// in TRD layer "plane".
if (!fDetPid) return -1;
- Float_t *trdMomentum=fDetPid->GetTRDmomentum();
+ const Double_t *trdMomentum=fDetPid->GetTRDmomentum();
if (!trdMomentum) {
return -1.;
}
//_______________________________________________________________________
-Int_t AliAODTrack::GetTOFBunchCrossing(Double_t b) const
+Int_t AliAODTrack::GetTOFBunchCrossing(Double_t b, Bool_t) const
{
// Returns the number of bunch crossings after trigger (assuming 25ns spacing)
const double kSpacing = 25e3; // min interbanch spacing
bcid = TMath::Nint((tdif - kShift)/kSpacing);
return bcid;
}
+
+void AliAODTrack::SetDetectorPID(const AliDetectorPID *pid)
+{
+ //
+ // Set the detector PID
+ //
+ if (fDetectorPID) delete fDetectorPID;
+ fDetectorPID=pid;
+
+}
+
+//_____________________________________________________________________________
+Double_t AliAODTrack::GetHMPIDsignal() const
+{
+ if(fAODEvent->GetHMPIDringForTrackID(fID)) return fAODEvent->GetHMPIDringForTrackID(fID)->GetHmpSignal();
+ else return -999.;
+}
+
+//_____________________________________________________________________________
+Double_t AliAODTrack::GetHMPIDoccupancy() const
+{
+ if(fAODEvent->GetHMPIDringForTrackID(fID)) return fAODEvent->GetHMPIDringForTrackID(fID)->GetHmpOccupancy();
+ else return -999.;
+}
+
+//_____________________________________________________________________________
+Int_t AliAODTrack::GetHMPIDcluIdx() const
+{
+ if(fAODEvent->GetHMPIDringForTrackID(fID)) return fAODEvent->GetHMPIDringForTrackID(fID)->GetHmpCluIdx();
+ else return -999;
+}
+
+//_____________________________________________________________________________
+void AliAODTrack::GetHMPIDtrk(Float_t &x, Float_t &y, Float_t &th, Float_t &ph) const
+{
+ x = -999; y = -999.; th = -999.; ph = -999.;
+
+ const AliAODHMPIDrings *ring=fAODEvent->GetHMPIDringForTrackID(fID);
+ if(ring){
+ x = ring->GetHmpTrackX();
+ y = ring->GetHmpTrackY();
+ th = ring->GetHmpTrackTheta();
+ ph = ring->GetHmpTrackPhi();
+ }
+}
+
+//_____________________________________________________________________________
+void AliAODTrack::GetHMPIDmip(Float_t &x,Float_t &y,Int_t &q, Int_t &nph) const
+{
+ x = -999; y = -999.; q = -999; nph = -999;
+
+ const AliAODHMPIDrings *ring=fAODEvent->GetHMPIDringForTrackID(fID);
+ if(ring){
+ x = ring->GetHmpMipX();
+ y = ring->GetHmpMipY();
+ q = (Int_t)ring->GetHmpMipCharge();
+ nph = (Int_t)ring->GetHmpNumOfPhotonClusters();
+ }
+}
+
+//_____________________________________________________________________________
+Bool_t AliAODTrack::GetOuterHmpPxPyPz(Double_t *p) const
+{
+ if(fAODEvent->GetHMPIDringForTrackID(fID)) {fAODEvent->GetHMPIDringForTrackID(fID)->GetHmpMom(p); return kTRUE;}
+
+ else return kFALSE;
+}
+//_____________________________________________________________________________
+Bool_t AliAODTrack::GetXYZAt(Double_t x, Double_t b, Double_t *r) const
+{
+ //---------------------------------------------------------------------
+ // This function returns the global track position extrapolated to
+ // the radial position "x" (cm) in the magnetic field "b" (kG)
+ //---------------------------------------------------------------------
+
+ //conversion of track parameter representation is
+ //based on the implementation of AliExternalTrackParam::Set(...)
+ //maybe some of this code can be moved to AliVTrack to avoid code duplication
+ const double kSafe = 1e-5;
+ Double_t alpha=0.0;
+ Double_t radPos2 = fPosition[0]*fPosition[0]+fPosition[1]*fPosition[1];
+ Double_t radMax = 45.; // approximately ITS outer radius
+ if (radPos2 < radMax*radMax) { // inside the ITS
+ alpha = TMath::ATan2(fMomentum[1],fMomentum[0]);
+ } else { // outside the ITS
+ Float_t phiPos = TMath::Pi()+TMath::ATan2(-fPosition[1], -fPosition[0]);
+ alpha =
+ TMath::DegToRad()*(20*((((Int_t)(phiPos*TMath::RadToDeg()))/20))+10);
+ }
+ //
+ Double_t cs=TMath::Cos(alpha), sn=TMath::Sin(alpha);
+ // protection: avoid alpha being too close to 0 or +-pi/2
+ if (TMath::Abs(sn)<kSafe) {
+ alpha = kSafe;
+ cs=TMath::Cos(alpha);
+ sn=TMath::Sin(alpha);
+ }
+ else if (cs<kSafe) {
+ alpha -= TMath::Sign(kSafe, alpha);
+ cs=TMath::Cos(alpha);
+ sn=TMath::Sin(alpha);
+ }
+
+ // Get the vertex of origin and the momentum
+ TVector3 ver(fPosition[0],fPosition[1],fPosition[2]);
+ TVector3 mom(fMomentum[0],fMomentum[1],fMomentum[2]);
+ //
+ // avoid momenta along axis
+ if (TMath::Abs(mom[0])<kSafe) mom[0] = TMath::Sign(kSafe*TMath::Abs(mom[1]), mom[0]);
+ if (TMath::Abs(mom[1])<kSafe) mom[1] = TMath::Sign(kSafe*TMath::Abs(mom[0]), mom[1]);
+
+ // Rotate to the local coordinate system
+ ver.RotateZ(-alpha);
+ mom.RotateZ(-alpha);
+
+ Double_t param0 = ver.Y();
+ Double_t param1 = ver.Z();
+ Double_t param2 = TMath::Sin(mom.Phi());
+ Double_t param3 = mom.Pz()/mom.Pt();
+ Double_t param4 = TMath::Sign(1/mom.Pt(),(Double_t)fCharge);
+
+ //calculate the propagated coordinates
+ //this is based on AliExternalTrackParam::GetXYZAt(Double_t x, Double_t b, Double_t *r)
+ Double_t dx=x-ver.X();
+ if(TMath::Abs(dx)<=kAlmost0) return GetXYZ(r);
+
+ Double_t f1=param2;
+ Double_t f2=f1 + dx*param4*b*kB2C;
+
+ if (TMath::Abs(f1) >= kAlmost1) return kFALSE;
+ if (TMath::Abs(f2) >= kAlmost1) return kFALSE;
+
+ Double_t r1=TMath::Sqrt((1.-f1)*(1.+f1)), r2=TMath::Sqrt((1.-f2)*(1.+f2));
+ r[0] = x;
+ r[1] = param0 + dx*(f1+f2)/(r1+r2);
+ r[2] = param1 + dx*(r2 + f2*(f1+f2)/(r1+r2))*param3;//Thanks to Andrea & Peter
+
+ return Local2GlobalPosition(r,alpha);
+}
+
+
+//_______________________________________________________
+void AliAODTrack::GetITSdEdxSamples(Double_t s[4]) const
+{
+ // get ITS dedx samples
+ if (!fDetPid) for (int i=4;i--;) s[i]=0;
+ else for (int i=4;i--;) s[i] = fDetPid->GetITSdEdxSample(i);
+}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff