#include "TDatabasePDG.h"
#include "TPDGCode.h"
#include "TVectorT.h"
#include "AliTrackReference.h"
#include "AliTrackPointArray.h"
#include "AliExternalTrackParam.h"
#include "AliTRDseedV1.h"
#include "AliTRDtrackV1.h"
#include "AliTRDgeometry.h"
#include "AliTRDtrackerV1.h"
#include "AliTRDtrackInfo.h"
ClassImp(AliTRDtrackInfo)
ClassImp(AliTRDtrackInfo::AliMCinfo)
ClassImp(AliTRDtrackInfo::AliESDinfo)
Double_t AliTRDtrackInfo::AliMCinfo::fgKalmanStep = 2.;
Bool_t AliTRDtrackInfo::AliMCinfo::fgKalmanUpdate = kTRUE;
AliTRDtrackInfo::AliTRDtrackInfo():
TObject()
,fNClusters(0)
,fTRDtrack(NULL)
,fMC(NULL)
,fESD()
{
}
AliTRDtrackInfo::AliTRDtrackInfo(const AliTRDtrackInfo &trdInfo):
TObject((const TObject&)trdInfo)
,fNClusters(trdInfo.fNClusters)
,fTRDtrack(NULL)
,fMC(NULL)
,fESD(trdInfo.fESD)
{
if(trdInfo.fMC) fMC = new AliMCinfo(*trdInfo.fMC);
SetTrack(trdInfo.fTRDtrack);
}
AliTRDtrackInfo::AliMCinfo::AliMCinfo()
:fLabel(0)
,fTRDlabel(0)
,fPDG(0)
,fNTrackRefs(0)
,fEta(-999.)
,fPhi(-999.)
,fPt(-1.)
{
memset(fTrackRefs, 0, sizeof(AliTrackReference *) * 12);
}
AliTRDtrackInfo::AliMCinfo::AliMCinfo(const AliMCinfo &mc)
:fLabel(mc.fLabel)
,fTRDlabel(mc.fTRDlabel)
,fPDG(mc.fPDG)
,fNTrackRefs(mc.fNTrackRefs)
,fEta(mc.fEta)
,fPhi(mc.fPhi)
,fPt(mc.fPt)
{
memset(fTrackRefs, 0, sizeof(AliTrackReference *) * 12);
for(Int_t ien = 0; ien < 12; ien++){
if(mc.fTrackRefs[ien])
fTrackRefs[ien] = new AliTrackReference(*(mc.fTrackRefs[ien]));
}
}
Int_t AliTRDtrackInfo::AliMCinfo::GetPID() const
{
switch(fPDG){
case kElectron:
case kPositron: return AliPID::kElectron;
case kMuonPlus:
case kMuonMinus: return AliPID::kMuon;
case kPiPlus:
case kPiMinus: return AliPID::kPion;
case kKPlus:
case kKMinus: return AliPID::kKaon;
case kProton:
case kProtonBar: return AliPID::kProton;
}
return -1;
}
AliTRDtrackInfo::AliESDinfo::AliESDinfo()
:fSteer(0)
,fId(-1)
,fStatus(0)
,fKinkIndex(0)
,fTPCncls(0)
,fTPCdedx(0.)
,fTOFbeta(0.)
,fTOFbc(0)
,fTRDpidQuality(0)
,fTRDnSlices(0)
,fPt(0.)
,fPhi(-999.)
,fEta(-999.)
,fTRDslices(NULL)
,fOP(NULL)
,fTPCout(NULL)
,fITSout(NULL)
,fTPArray(NULL)
{
memset(fTRDr, 0, AliPID::kSPECIES*sizeof(Double32_t));
memset(fTRDv0pid, 0, AliPID::kSPECIES*sizeof(Int_t));
}
AliTRDtrackInfo::AliESDinfo::AliESDinfo(const AliESDinfo &esd)
:fSteer(esd.fSteer)
,fId(esd.fId)
,fStatus(esd.fStatus)
,fKinkIndex(esd.fKinkIndex)
,fTPCncls(esd.fTPCncls)
,fTPCdedx(esd.fTPCdedx)
,fTOFbeta(esd.fTOFbeta)
,fTOFbc(esd.fTOFbc)
,fTRDpidQuality(esd.fTRDpidQuality)
,fTRDnSlices(esd.fTRDnSlices)
,fPt(esd.fPt)
,fPhi(esd.fPhi)
,fEta(esd.fEta)
,fTRDslices(NULL)
,fOP(NULL)
,fTPCout(NULL)
,fITSout(NULL)
,fTPArray(NULL)
{
memcpy(fTRDr, esd.fTRDr, AliPID::kSPECIES*sizeof(Double32_t));
memcpy(fTRDv0pid, esd.fTRDv0pid, AliPID::kSPECIES*sizeof(Int_t));
if(fTRDnSlices){
fTRDslices = new Double32_t[fTRDnSlices];
memcpy(fTRDslices, esd.fTRDslices, fTRDnSlices*sizeof(Double32_t));
}
SetOuterParam(esd.fOP);
SetTPCoutParam(esd.fTPCout);
SetITSoutParam(esd.fITSout);
SetTrackPointArray(esd.fTPArray);
}
AliTRDtrackInfo::~AliTRDtrackInfo()
{
if(fMC) delete fMC;
if(fTRDtrack) delete fTRDtrack;
}
AliTRDtrackInfo::AliMCinfo::~AliMCinfo()
{
fNTrackRefs = 0;
for(Int_t ien = 0; ien < 12; ien++){
if(fTrackRefs[ien]) delete fTrackRefs[ien];
fTrackRefs[ien] = NULL;
}
}
AliTRDtrackInfo::AliESDinfo::~AliESDinfo()
{
if(fTRDnSlices){
delete [] fTRDslices;
fTRDslices = NULL;
fTRDnSlices = 0;
}
if(fOP) delete fOP; fOP = NULL;
if(fTPCout) delete fTPCout; fTPCout = NULL;
if(fITSout) delete fITSout; fITSout = NULL;
if(fTPArray) delete fTPArray;
}
void AliTRDtrackInfo::AliESDinfo::Delete(const Option_t *){
if(fTRDnSlices){
delete [] fTRDslices;
fTRDslices = NULL;
fTRDnSlices = 0;
}
if(fOP) delete fOP; fOP = NULL;
if(fTPCout) delete fTPCout; fTPCout = NULL;
if(fITSout) delete fITSout; fITSout = NULL;
if(fTPArray) delete fTPArray; fTPArray = NULL;
}
AliTRDtrackInfo& AliTRDtrackInfo::operator=(const AliTRDtrackInfo &trdInfo)
{
if(this == &trdInfo) return *this;
fNClusters = trdInfo.fNClusters;
fESD = trdInfo.fESD;
if(trdInfo.fMC){
if(!fMC) fMC = new AliMCinfo(*trdInfo.fMC);
else{
fMC->~AliMCinfo();
new(fMC) AliMCinfo(*trdInfo.fMC);
}
} else {if(fMC) delete fMC; fMC = NULL;}
SetTrack(trdInfo.fTRDtrack);
return *this;
}
AliTRDtrackInfo::AliMCinfo& AliTRDtrackInfo::AliMCinfo::operator=(const AliMCinfo &mc)
{
if(this == &mc) return *this;
fLabel = mc.fLabel;
fTRDlabel = mc.fTRDlabel;
fPDG = mc.fPDG;
fNTrackRefs = mc.fNTrackRefs;
AliTrackReference **itr = &fTrackRefs[0];
AliTrackReference* const *jtr = &mc.fTrackRefs[0];
for(Int_t ien = 0; ien < 2*AliTRDgeometry::kNlayer; ien++, itr++, jtr++){
if((*jtr)){
if(!(*itr)) (*itr) = new AliTrackReference(*(*jtr));
else{
(*itr)->~AliTrackReference();
new(&(*itr)) AliTrackReference(*(*jtr));
}
} else {if((*itr)) delete (*itr); (*itr) = NULL;}
}
return *this;
}
AliTRDtrackInfo::AliESDinfo& AliTRDtrackInfo::AliESDinfo::operator=(const AliESDinfo &esd)
{
if(this == &esd) return *this;
fSteer = esd.fSteer;
fId = esd.fId;
fStatus = esd.fStatus;
fKinkIndex = esd.fKinkIndex;
fTPCncls = esd.fTPCncls;
fTPCdedx = esd.fTPCdedx;
fTOFbeta = esd.fTOFbeta;
fTOFbc = esd.fTOFbc;
fTRDpidQuality= esd.fTRDpidQuality;
fTRDnSlices = esd.fTRDnSlices;
fPt = esd.fPt;
fPhi = esd.fPhi;
fEta = esd.fEta;
memcpy(fTRDr, esd.fTRDr, AliPID::kSPECIES*sizeof(Double32_t));
memcpy(fTRDv0pid, esd.fTRDv0pid, AliPID::kSPECIES*sizeof(Int_t));
if(fTRDnSlices){
if(!fTRDslices) fTRDslices = new Double32_t[fTRDnSlices];
memcpy(fTRDslices, esd.fTRDslices, fTRDnSlices*sizeof(Double32_t));
}
SetOuterParam(esd.fOP);
SetTPCoutParam(esd.fTPCout);
SetITSoutParam(esd.fITSout);
SetTrackPointArray(esd.fTPArray);
return *this;
}
void AliTRDtrackInfo::Delete(const Option_t *)
{
AliDebug(2, Form("track[%p] mc[%p]", (void*)fTRDtrack, (void*)fMC));
fNClusters = 0;
if(fMC) delete fMC; fMC = NULL;
fESD.Delete(NULL);
if(fTRDtrack) delete fTRDtrack; fTRDtrack = NULL;
}
void AliTRDtrackInfo::SetTrack(const AliTRDtrackV1 *track)
{
if(track){
if(!fTRDtrack) fTRDtrack = new AliTRDtrackV1(*track);
else{
fTRDtrack->~AliTRDtrackV1();
new(fTRDtrack) AliTRDtrackV1(*track);
}
if(track->IsOwner()) fTRDtrack->SetOwner();
} else {
if(fTRDtrack) delete fTRDtrack; fTRDtrack = NULL;
}
}
void AliTRDtrackInfo::AliESDinfo::SetTrackPointArray(const AliTrackPointArray *tps)
{
if(tps){
if(!fTPArray) fTPArray = new AliTrackPointArray(*tps);
else{
fTPArray->~AliTrackPointArray();
new(fTPArray) AliTrackPointArray(*tps);
}
} else {
if(fTPArray) delete fTPArray; fTPArray = NULL;
}
}
void AliTRDtrackInfo::AddTrackRef(const AliTrackReference *tref)
{
if(fMC->fNTrackRefs >= 2*AliTRDgeometry::kNlayer){
SetCurved();
return;
}
fMC->fTrackRefs[fMC->fNTrackRefs++] = new AliTrackReference(*tref);
}
AliTrackReference* AliTRDtrackInfo::GetTrackRef(Int_t idx) const
{
if(!fMC) return NULL;
return (idx>=0 && idx < 12) ? fMC->fTrackRefs[idx] : NULL;
}
AliTrackReference* AliTRDtrackInfo::GetTrackRef(const AliTRDseedV1* const tracklet) const
{
if(!fMC) return NULL;
Double_t cw = AliTRDgeometry::CamHght() + AliTRDgeometry::CdrHght();
AliTrackReference * const* jtr = &(fMC->fTrackRefs[0]);
for(Int_t itr = 0; itr < fMC->fNTrackRefs; itr++, ++jtr){
if(!(*jtr)) break;
if(TMath::Abs(tracklet->GetX0() - (*jtr)->LocalX()) < cw) return (*jtr);
}
return NULL;
}
Int_t AliTRDtrackInfo::GetNumberOfClusters() const
{
Int_t n = 0;
if(!fTRDtrack) return 0;
if(fTRDtrack->GetNumberOfTracklets() == 0) return n;
AliTRDseedV1 *tracklet = NULL;
for(Int_t ip=0; ip<AliTRDgeometry::kNlayer; ip++){
if(!(tracklet = const_cast<AliTRDseedV1 *>(fTRDtrack->GetTracklet(ip)))) continue;
n+=tracklet->GetN();
}
return n;
}
void AliTRDtrackInfo::AliESDinfo::SetOuterParam(const AliExternalTrackParam *op)
{
if(op){
if(fOP){
fOP->~AliExternalTrackParam();
new(fOP) AliExternalTrackParam(*op);
} else fOP = new AliExternalTrackParam(*op);
} else {
if(fOP) delete fOP; fOP = NULL;
}
}
void AliTRDtrackInfo::AliESDinfo::SetITSoutParam(const AliExternalTrackParam *op)
{
if(op){
if(fITSout){
fITSout->~AliExternalTrackParam();
new(fITSout) AliExternalTrackParam(*op);
} else fITSout = new AliExternalTrackParam(*op);
} else {
if(fITSout) delete fITSout; fITSout = NULL;
}
}
void AliTRDtrackInfo::AliESDinfo::SetTPCoutParam(const AliExternalTrackParam *op)
{
if(op){
if(fTPCout){
fTPCout->~AliExternalTrackParam();
new(fTPCout) AliExternalTrackParam(*op);
} else fTPCout = new AliExternalTrackParam(*op);
} else {
if(fTPCout) delete fTPCout; fTPCout = NULL;
}
}
Int_t AliTRDtrackInfo::GetNTracklets() const
{
if(!fTRDtrack) return 0;
return fTRDtrack->GetNumberOfTracklets();
}
void AliTRDtrackInfo::SetSlices(Int_t n, Double32_t *s)
{
if(fESD.fTRDnSlices != n){
fESD.fTRDnSlices = 0;
delete [] fESD.fTRDslices;
fESD.fTRDslices = NULL;
}
if(!fESD.fTRDnSlices){
fESD.fTRDnSlices = n;
fESD.fTRDslices = new Double32_t[fESD.fTRDnSlices];
}
memcpy(fESD.fTRDslices, s, n*sizeof(Double32_t));
}
Bool_t AliTRDtrackInfo::AliMCinfo::GetDirections(Float_t &x0, Float_t &y0, Float_t &z0, Float_t &dydx, Float_t &dzdx, Float_t &pt, Float_t &p, Float_t &eta, Float_t &phi, UChar_t &status) const
{
status = 0;
Double_t cw = AliTRDgeometry::CamHght() + AliTRDgeometry::CdrHght();
Int_t nFound = 0;
AliTrackReference *tr[2] = {NULL, NULL};
AliTrackReference * const* jtr = &fTrackRefs[0];
for(Int_t itr = 0; itr < fNTrackRefs; itr++, ++jtr){
if(!(*jtr)) break;
if(TMath::Abs(x0 - (*jtr)->LocalX()) > cw) continue;
tr[nFound++] = (*jtr);
if(nFound == 2) break;
}
if(nFound < 2){
if(!nFound) SETBIT(status, 0);
else SETBIT(status, 1);
return kFALSE;
}
pt=tr[1]->Pt();
p =tr[1]->P();
if(pt < 1.e-3) return kFALSE;
Double_t dx = tr[1]->LocalX() - tr[0]->LocalX();
if(dx <= 0.){
AliWarningGeneral("AliTRDtrackInfo::AliMCinfo::GetDirections()", Form("Track ref with wrong radial distances refX0[%6.3f] refX1[%6.3f]", tr[0]->LocalX(), tr[1]->LocalX()));
SETBIT(status, 2);
return kFALSE;
}
if(TMath::Abs(dx-AliTRDgeometry::CamHght()-AliTRDgeometry::CdrHght())>1.E-3) SETBIT(status, 3);
dydx = (tr[1]->LocalY() - tr[0]->LocalY()) / dx;
dzdx = (tr[1]->Z() - tr[0]->Z()) / dx;
y0 = tr[1]->LocalY();
z0 = tr[1]->Z();
x0 = tr[1]->LocalX();
eta = -TMath::Log(TMath::Tan(0.5 * tr[1]->Theta()));
phi = TMath::ATan2(tr[1]->Y(), tr[1]->X());
return kTRUE;
}
Bool_t AliTRDtrackInfo::AliMCinfo::PropagateKalman(
TVectorD *x, TVectorD *y, TVectorD *z,
TVectorD *dx, TVectorD *dy, TVectorD *dz,
TVectorD *pt, TVectorD *dpt, TVectorD *budget, TVectorD *c, Double_t mass) const
{
if(!fNTrackRefs) return kFALSE;
for(Int_t itr=kNTrackRefs; itr--;){
(*x)[itr] = 0.;(*y)[itr] = 0.;(*z)[itr] = 0.;
(*dx)[itr] = -1.; (*dy)[itr] = 100.; (*dz)[itr] = 100.; (*dpt)[itr] = 100.;
}
AliTrackReference *tr(NULL);
Int_t itr(0); while(!(tr = fTrackRefs[itr])) itr++;
if(tr->Pt()<1.e-3) return kFALSE;
AliTRDtrackV1 tt;
Double_t xyz[3]={tr->X(),tr->Y(),tr->Z()};
Double_t pxyz[3]={tr->Px(),tr->Py(),tr->Pz()};
Double_t var[6] = {1.e-4, 1.e-4, 1.e-4, 1.e-4, 1.e-4, 1.e-4};
Double_t cov[21]={
var[0], 0., 0., 0., 0., 0.,
var[1], 0., 0., 0., 0.,
var[2], 0., 0., 0.,
var[3], 0., 0.,
var[4], 0.,
var[5]
};
TDatabasePDG db;
const TParticlePDG *pdg=db.GetParticle(fPDG);
if(!pdg){
AliWarningGeneral("AliTRDtrackInfo::AliMCinfo::PropagateKalman()", Form("PDG entry missing for code %d. References for track %d", fPDG, fNTrackRefs));
return kFALSE;
}
tt.Set(xyz, pxyz, cov, Short_t(pdg->Charge()));
if(mass<0){
tt.SetMass(pdg->Mass());
} else {
tt.SetMass(mass);
}
Double_t x0(tr->LocalX());
const Double_t *cc(NULL);
for(Int_t ip=0; itr<fNTrackRefs; itr++){
if(!(tr = fTrackRefs[itr])) continue;
if(!AliTRDtrackerV1::PropagateToX(tt, tr->LocalX(), fgKalmanStep)) continue;
tt.GetXYZ(xyz);
(*x)[ip] = xyz[0];
(*y)[ip] = xyz[1];
(*z)[ip] = xyz[2];
(*dx)[ip] = tt.GetX() - x0;
(*dy)[ip] = tt.GetY() - tr->LocalY();
(*dz)[ip] = tt.GetZ() - tr->Z();
(*pt)[ip] = tr->Pt();
(*dpt)[ip] = tt.Pt()- tr->Pt();
(*budget)[ip] = tt.GetBudget(0);
cc = tt.GetCovariance();
for(Int_t ic(0), jp(ip*15); ic<15; ic++, jp++) (*c)[jp]=cc[ic];
ip++;
}
return kTRUE;
}