#include "AliIDFFUtils.h"
#include "AliAODEvent.h"
#include "AliAODMCParticle.h"
#include "AliESDtrack.h"
#include "TClonesArray.h"
ClassImp(AliIDFFUtils);
AliPIDResponse * AliIDFFUtils::fPid=0x0;
Int_t AliIDFFUtils::PDG2Type(const Int_t pdg)
{
Int_t itype = kNOTSELECTED;
switch(pdg){
case 11:
itype = kELECTRON;
break;
case 211:
itype = kPION;
break;
case 2212:
itype = kPROTON;
break;
case 321:
itype = kKAON;
break;
default:
itype = kNOTSELECTED;
break;
}
return itype;
}
THnSparseD *AliIDFFUtils::GetTHn(const TString name)
{
const Int_t nvar = 11;
const TString atitle[nvar]={"TrackEta","JetPt", "TrackTPCsig", "Log10TrackP", "Log10TrackPt", "z", "xi", "pdg", "comb", "tof", "tpc"};
const Int_t nbins[nvar] ={ 4, 15, 1200, Nx(), 50, 30, 60, 7, 7, 7, 7};
const Double_t xmins[nvar]={ 0, 5, 0, Xmin(), Xmin(), 0, 0, -3.5, -3.5, -3.5, -3.5};
const Double_t xmaxs[nvar]={ 0.9, 20, 200, Xmax(), Xmax(), 1.2, 6, 3.5, 3.5, 3.5, 3.5};
THnSparseD * hh = new THnSparseD(name,"", nvar, nbins, xmins, xmaxs);
for(Int_t ia=0; ia<nvar; ia++){
hh->GetAxis(ia)->SetTitle(atitle[ia]);
}
TAxis * ax = 0x0;
Int_t nb = 0;
ax = hh->GetAxis(0);
nb = ax->GetNbins();
const Double_t etas[]={0, 0.2, 0.4, 0.6, 0.9};
ax->Set(nb, etas);
return hh;
}
Int_t AliIDFFUtils::TPCType(const AliAODTrack * trackptr)
{
const AliPIDResponse::EDetPidStatus tpcstatus = fPid->CheckPIDStatus(AliPIDResponse::kTPC, trackptr);
if(tpcstatus != AliPIDResponse::kDetPidOk)
return kNOINFO;
Double_t nsigma[]={-999,-999,-999, -999};
nsigma[kPION] = fPid->NumberOfSigmasTPC( trackptr, AliPID::kPion);
nsigma[kKAON] = fPid->NumberOfSigmasTPC( trackptr, AliPID::kKaon);
nsigma[kPROTON] = fPid->NumberOfSigmasTPC( trackptr, AliPID::kProton);
nsigma[kELECTRON] = fPid->NumberOfSigmasTPC( trackptr, AliPID::kElectron);
const Double_t inclusion=5;
const Double_t exclusion=5;
const Double_t maxsig = 150;
if(trackptr->GetTPCsignal()> maxsig){
if(TMath::Abs(nsigma[kPION])<inclusion && TMath::Abs(nsigma[kKAON])>exclusion && TMath::Abs(nsigma[kPROTON])>exclusion && TMath::Abs(nsigma[kELECTRON])>exclusion) return kPION;
if(TMath::Abs(nsigma[kPION])>exclusion && TMath::Abs(nsigma[kKAON])<inclusion && TMath::Abs(nsigma[kPROTON])>exclusion && TMath::Abs(nsigma[kELECTRON])>exclusion) return kKAON;
if(TMath::Abs(nsigma[kPION])>exclusion && TMath::Abs(nsigma[kKAON])>exclusion && TMath::Abs(nsigma[kPROTON])<inclusion && TMath::Abs(nsigma[kELECTRON])>exclusion) return kPROTON;
if(TMath::Abs(nsigma[kPION])>exclusion && TMath::Abs(nsigma[kKAON])>exclusion && TMath::Abs(nsigma[kPROTON])>exclusion && TMath::Abs(nsigma[kELECTRON])<inclusion) return kELECTRON;
}
return kNOTSELECTED;
}
Int_t AliIDFFUtils::TOFType(const AliAODTrack * trackptr, const Int_t tofmode)
{
const AliPIDResponse::EDetPidStatus tofstatus = fPid->CheckPIDStatus(AliPIDResponse::kTOF, trackptr);
if(tofstatus != AliPIDResponse::kDetPidOk)
return kNOINFO;
Double_t nsigma[]={-999,-999,-999, -999};
nsigma[kPION] = fPid->NumberOfSigmasTOF( trackptr, AliPID::kPion);
nsigma[kKAON] = fPid->NumberOfSigmasTOF( trackptr, AliPID::kKaon);
nsigma[kPROTON] = fPid->NumberOfSigmasTOF( trackptr, AliPID::kProton);
nsigma[kELECTRON] = fPid->NumberOfSigmasTOF( trackptr, AliPID::kElectron);
Double_t inclusion=-999;
Double_t exclusion=-999;
if(tofmode == 1){
inclusion = 2;
exclusion = 2;
}
else if(tofmode == 2){
inclusion = 2;
exclusion = 3;
}
else if(tofmode == 3){
inclusion = 3;
exclusion = 3;
}
else if(tofmode == 4){
inclusion = 3;
exclusion = 4;
}
else{
printf("AliIDFFUtils::TOFType bad tofmode ! %d\n", tofmode); exit(1);
}
const Bool_t kpassEle = kTRUE;
if(TMath::Abs(nsigma[kPION])<inclusion && TMath::Abs(nsigma[kKAON])>exclusion && TMath::Abs(nsigma[kPROTON])>exclusion && kpassEle) return kPION;
if(TMath::Abs(nsigma[kPION])>exclusion && TMath::Abs(nsigma[kKAON])<inclusion && TMath::Abs(nsigma[kPROTON])>exclusion && kpassEle) return kKAON;
if(TMath::Abs(nsigma[kPION])>exclusion && TMath::Abs(nsigma[kKAON])>exclusion && TMath::Abs(nsigma[kPROTON])<inclusion && kpassEle) return kPROTON;
if(TMath::Abs(nsigma[kPION])>exclusion+0.5 && TMath::Abs(nsigma[kKAON])>exclusion && TMath::Abs(nsigma[kPROTON])>exclusion && TMath::Abs(nsigma[kELECTRON])<inclusion) return kELECTRON;
return kNOTSELECTED;
}
Int_t AliIDFFUtils::CombineTPCTOF(const Int_t ktpc, const Int_t ktof)
{
if(ktpc == ktof)
return ktpc;
else if(ktpc < 0 && ktof >= 0 )
return ktof;
else if(ktof < 0 && ktpc >= 0)
return ktpc;
else
return kNOTACCEPTED;
}
void AliIDFFUtils::FillTHn(THnSparseD * hh, Double_t jetpt, const AliAODTrack * track, AliAODEvent *aodevent, const Int_t tofmode)
{
Int_t mcpdg = kNOINFO;
TClonesArray *tca = dynamic_cast<TClonesArray*>(aodevent->FindListObject(AliAODMCParticle::StdBranchName()));
if(tca){
const Int_t mclabel = TMath::Abs(track->GetLabel());
AliAODMCParticle* gentrack = dynamic_cast<AliAODMCParticle*> (tca->At(mclabel));
if(gentrack){
mcpdg = PDG2Type(TMath::Abs(gentrack->GetPdgCode()));
}
}
const Int_t ktpc = TPCType(track);
const Int_t ktof = TOFType(track, tofmode);
const Int_t kcomb = CombineTPCTOF(ktpc, ktof);
const Double_t eps = 1e-6;
const Double_t tracketa = TMath::Abs(track->Eta());
Double_t tracktpc = track->GetTPCsignal();
if(tracktpc>=200)
tracktpc = 200 - eps;
const Double_t tracklogmom = TMath::Log10(track->P());
const Double_t tracklogpt = TMath::Log10(track->Pt());
Double_t zz = -999;
Double_t xi = -999;
if(jetpt<-990){
jetpt = zz = xi = eps;
}
else if(track->Pt()>(1-eps)*jetpt && track->Pt()<(1+eps)*jetpt){
zz = 1-eps;
xi = eps;
}
else if(jetpt>0){
zz = track->Pt()/jetpt;
xi = TMath::Log(1/zz);
}
const Double_t vars[]={tracketa, jetpt, tracktpc, tracklogmom, tracklogpt, zz, xi, (Double_t)mcpdg, (Double_t)kcomb, (Double_t)ktof, (Double_t)ktpc};
hh->Fill(vars);
}
Bool_t AliIDFFUtils::TPCCutPIDN(const AliAODTrack * track)
{
if(track->GetTPCsignalN()<60){
return kFALSE;
}
return kTRUE;
}
Bool_t AliIDFFUtils::TPCCutMIGeo(const AliAODTrack * track, const AliVEvent* evt, TTreeStream * streamer)
{
Short_t sign = track->Charge();
Double_t xyz[50];
Double_t pxpypz[50];
Double_t cv[100];
track->GetXYZ(xyz);
track->GetPxPyPz(pxpypz);
AliExternalTrackParam * par = new AliExternalTrackParam(xyz, pxpypz, cv, sign);
AliESDtrack dummy;
Double_t varGeom = dummy.GetLengthInActiveZone(par,3,236, evt->GetMagneticField(), 0,0);
Double_t varNcr = track->GetTPCClusterInfo(3,1);
Double_t varNcls = track->GetTPCsignalN();
Bool_t cutGeom = varGeom > 1.*(130-5*TMath::Abs(1./track->Pt()));
Bool_t cutNcr = varNcr > 0.85*(130-5*TMath::Abs(1./track->Pt()));
Bool_t cutNcls = varNcls > 0.7*(130-5*TMath::Abs(1./track->Pt()));
Bool_t kout = cutGeom && cutNcr && cutNcls;
if(streamer){
Double_t dedx = track->GetTPCsignal();
(*streamer)<<"tree"<<
"param.="<< par<<
"varGeom="<<varGeom<<
"varNcr="<<varNcr<<
"varNcls="<<varNcls<<
"cutGeom="<<cutGeom<<
"cutNcr="<<cutNcr<<
"cutNcls="<<cutNcls<<
"kout="<<kout<<
"dedx="<<dedx<<
"\n";
}
delete par;
return kout;
}