ClassImp(AliFemtoEventReaderESDChainKine)
#if !(ST_NO_NAMESPACES)
- using namespace units;
+using namespace units;
#endif
using namespace std;
(fEvent->IsTriggerClassFired("CINT1-B-NOPF-FASTNOTRD")))
hbtEvent->SetTriggerCluster(1);
else if ((fEvent->IsTriggerClassFired("CSH1WU-B-NOPF-ALL")) ||
- (fEvent->IsTriggerClassFired("CSH1-B-NOPF-ALLNOTRD")))
+ (fEvent->IsTriggerClassFired("CSH1-B-NOPF-ALLNOTRD")))
hbtEvent->SetTriggerCluster(2);
else
hbtEvent->SetTriggerCluster(0);
if (cdh) {
TList *tGenHeaders = cdh->GetHeaders();
for (int ihead = 0; ihead<tGenHeaders->GetEntries(); ihead++) {
- hdh = dynamic_cast<AliGenHijingEventHeader *> (fGenHeader);
- if (hdh) break;
+ hdh = dynamic_cast<AliGenHijingEventHeader *> (fGenHeader);
+ if (hdh) break;
}
}
}
if (hdh)
- {
- tReactionPlane = hdh->ReactionPlaneAngle();
- cout << "Got reaction plane " << tReactionPlane << endl;
- }
+ {
+ tReactionPlane = hdh->ReactionPlaneAngle();
+ cout << "Got reaction plane " << tReactionPlane << endl;
+ }
hbtEvent->SetReactionPlaneAngle(tReactionPlane);
for (int ip=0; ip<fStack->GetNtrack(); ip++) {
motherpart = fStack->Particle(ip);
if (motherpart->GetDaughter(0) > 0)
- motherids[motherpart->GetDaughter(0)] = ip;
+ motherids[motherpart->GetDaughter(0)] = ip;
if (motherpart->GetDaughter(1) > 0)
- motherids[motherpart->GetDaughter(1)] = ip;
+ motherids[motherpart->GetDaughter(1)] = ip;
// if (motherpart->GetPdgCode() == 211) {
// cout << "Mother " << ip << " has daughters "
}
for (int i=0;i<nofTracks;i++)
- {
- // cout << "Reading track " << i << endl;
- bool tGoodMomentum=true; //flaga to chcek if we can read momentum of this track
-
-
- const AliESDtrack *esdtrack=fEvent->GetTrack(i);//getting next track
- // const AliESDfriendTrack *tESDfriendTrack = esdtrack->GetFriendTrack();
-
-
- if ((esdtrack->GetStatus() & AliESDtrack::kTPCrefit) &&
- (esdtrack->GetStatus() & AliESDtrack::kITSrefit)) {
- if (esdtrack->GetTPCNcls() > 70)
- if (esdtrack->GetTPCchi2()/esdtrack->GetTPCNcls() < 4.0) {
- if (TMath::Abs(esdtrack->Eta()) < 1.2) {
- esdtrack->GetImpactParameters(b,bCov);
- if ((b[0]<0.2) && (b[1] < 0.25)) {
- tNormMult++;
- tTotalPt += esdtrack->Pt();
- }
- }
- }
- }
- else if (esdtrack->GetStatus() & AliESDtrack::kTPCrefit) {
- if (esdtrack->GetTPCNcls() > 100)
- if (esdtrack->GetTPCchi2()/esdtrack->GetTPCNcls() < 4.0) {
- if (TMath::Abs(esdtrack->Eta()) < 1.2) {
- esdtrack->GetImpactParameters(b,bCov);
- if ((b[0]<2.4) && (b[1] < 3.2)) {
- tNormMult++;
- tTotalPt += esdtrack->Pt();
- }
- }
- }
- }
+ {
+ // cout << "Reading track " << i << endl;
+ bool tGoodMomentum=true; //flaga to chcek if we can read momentum of this track
+
+
+ const AliESDtrack *esdtrack=fEvent->GetTrack(i);//getting next track
+ // const AliESDfriendTrack *tESDfriendTrack = esdtrack->GetFriendTrack();
+
+
+ if ((esdtrack->GetStatus() & AliESDtrack::kTPCrefit) &&
+ (esdtrack->GetStatus() & AliESDtrack::kITSrefit)) {
+ if (esdtrack->GetTPCNcls() > 70)
+ if (esdtrack->GetTPCchi2()/esdtrack->GetTPCNcls() < 4.0) {
+ if (TMath::Abs(esdtrack->Eta()) < 1.2) {
+ esdtrack->GetImpactParameters(b,bCov);
+ if ((b[0]<0.2) && (b[1] < 0.25)) {
+ tNormMult++;
+ tTotalPt += esdtrack->Pt();
+ }
+ }
+ }
+ }
+ else if (esdtrack->GetStatus() & AliESDtrack::kTPCrefit) {
+ if (esdtrack->GetTPCNcls() > 100)
+ if (esdtrack->GetTPCchi2()/esdtrack->GetTPCNcls() < 4.0) {
+ if (TMath::Abs(esdtrack->Eta()) < 1.2) {
+ esdtrack->GetImpactParameters(b,bCov);
+ if ((b[0]<2.4) && (b[1] < 3.2)) {
+ tNormMult++;
+ tTotalPt += esdtrack->Pt();
+ }
+ }
+ }
+ }
- hbtEvent->SetZDCEMEnergy(tTotalPt);
- // if (esdtrack->GetStatus() & AliESDtrack::kTPCrefit)
- // if (esdtrack->GetTPCNcls() > 80)
- // if (esdtrack->GetTPCchi2()/esdtrack->GetTPCNcls() < 6.0)
- // if (esdtrack->GetConstrainedParam())
- // if (fabs(esdtrack->GetConstrainedParam()->Eta()) < 0.5)
- // if (esdtrack->GetConstrainedParam()->Pt() < 1.0) {
- // if (esdtrack->GetSign() > 0)
- // tNormMultPos++;
- // else if (esdtrack->GetSign() < 0)
- // tNormMultNeg--;
- // }
-
- // If reading ITS-only tracks, reject all with TPC
- if (fTrackType == kITSOnly) {
- if (esdtrack->GetStatus() & AliESDtrack::kTPCrefit) continue;
- if (!(esdtrack->GetStatus() & AliESDtrack::kITSrefit)) continue;
- if (esdtrack->GetStatus() & AliESDtrack::kTPCin) continue;
- UChar_t iclm = esdtrack->GetITSClusterMap();
- Int_t incls = 0;
- for (int iter=0; iter<6; iter++) if (iclm&(1<<iter)) incls++;
- if (incls<=3) {
- if (Debug()>1) cout << "Rejecting track with " << incls << " clusters" << endl;
- continue;
- }
+ hbtEvent->SetZDCEMEnergy(tTotalPt);
+ // if (esdtrack->GetStatus() & AliESDtrack::kTPCrefit)
+ // if (esdtrack->GetTPCNcls() > 80)
+ // if (esdtrack->GetTPCchi2()/esdtrack->GetTPCNcls() < 6.0)
+ // if (esdtrack->GetConstrainedParam())
+ // if (fabs(esdtrack->GetConstrainedParam()->Eta()) < 0.5)
+ // if (esdtrack->GetConstrainedParam()->Pt() < 1.0) {
+ // if (esdtrack->GetSign() > 0)
+ // tNormMultPos++;
+ // else if (esdtrack->GetSign() < 0)
+ // tNormMultNeg--;
+ // }
+
+ // If reading ITS-only tracks, reject all with TPC
+ if (fTrackType == kITSOnly) {
+ if (esdtrack->GetStatus() & AliESDtrack::kTPCrefit) continue;
+ if (!(esdtrack->GetStatus() & AliESDtrack::kITSrefit)) continue;
+ if (esdtrack->GetStatus() & AliESDtrack::kTPCin) continue;
+ UChar_t iclm = esdtrack->GetITSClusterMap();
+ Int_t incls = 0;
+ for (int iter=0; iter<6; iter++) if (iclm&(1<<iter)) incls++;
+ if (incls<=3) {
+ if (Debug()>1) cout << "Rejecting track with " << incls << " clusters" << endl;
+ continue;
}
+ }
- AliFemtoTrack* trackCopy = new AliFemtoTrack();
- trackCopy->SetCharge((short)esdtrack->GetSign());
+ AliFemtoTrack* trackCopy = new AliFemtoTrack();
+ trackCopy->SetCharge((short)esdtrack->GetSign());
- //in aliroot we have AliPID
- //0-electron 1-muon 2-pion 3-kaon 4-proton 5-photon 6-pi0 7-neutron 8-kaon0 9-eleCon
- //we use only 5 first
- double esdpid[5];
- esdtrack->GetESDpid(esdpid);
- trackCopy->SetPidProbElectron(esdpid[0]);
- trackCopy->SetPidProbMuon(esdpid[1]);
- trackCopy->SetPidProbPion(esdpid[2]);
- trackCopy->SetPidProbKaon(esdpid[3]);
- trackCopy->SetPidProbProton(esdpid[4]);
+ //in aliroot we have AliPID
+ //0-electron 1-muon 2-pion 3-kaon 4-proton 5-photon 6-pi0 7-neutron 8-kaon0 9-eleCon
+ //we use only 5 first
+ double esdpid[5];
+ esdtrack->GetESDpid(esdpid);
+ trackCopy->SetPidProbElectron(esdpid[0]);
+ trackCopy->SetPidProbMuon(esdpid[1]);
+ trackCopy->SetPidProbPion(esdpid[2]);
+ trackCopy->SetPidProbKaon(esdpid[3]);
+ trackCopy->SetPidProbProton(esdpid[4]);
- esdpid[0] = -100000.0;
- esdpid[1] = -100000.0;
- esdpid[2] = -100000.0;
- esdpid[3] = -100000.0;
- esdpid[4] = -100000.0;
+ esdpid[0] = -100000.0;
+ esdpid[1] = -100000.0;
+ esdpid[2] = -100000.0;
+ esdpid[3] = -100000.0;
+ esdpid[4] = -100000.0;
- double tTOF = 0.0;
+ double tTOF = 0.0;
- if (esdtrack->GetStatus()&AliESDtrack::kTOFpid) {
- tTOF = esdtrack->GetTOFsignal();
- esdtrack->GetIntegratedTimes(esdpid);
- }
+ if (esdtrack->GetStatus()&AliESDtrack::kTOFpid) {
+ tTOF = esdtrack->GetTOFsignal();
+ esdtrack->GetIntegratedTimes(esdpid);
+ }
- trackCopy->SetTofExpectedTimes(tTOF-esdpid[2], tTOF-esdpid[3], tTOF-esdpid[4]);
+ trackCopy->SetTofExpectedTimes(tTOF-esdpid[2], tTOF-esdpid[3], tTOF-esdpid[4]);
- ////// TPC ////////////////////////////////////////////
- float nsigmaTPCK=-1000.;
- float nsigmaTPCPi=-1000.;
- float nsigmaTPCP=-1000.;
+ ////// TPC ////////////////////////////////////////////
+ float nsigmaTPCK=-1000.;
+ float nsigmaTPCPi=-1000.;
+ float nsigmaTPCP=-1000.;
- if ((fESDpid) && (esdtrack->IsOn(AliESDtrack::kTPCpid))){
- nsigmaTPCK = fESDpid->NumberOfSigmasTPC(esdtrack,AliPID::kKaon);
- nsigmaTPCPi = fESDpid->NumberOfSigmasTPC(esdtrack,AliPID::kPion);
- nsigmaTPCP = fESDpid->NumberOfSigmasTPC(esdtrack,AliPID::kProton);
+ if ((fESDpid) && (esdtrack->IsOn(AliESDtrack::kTPCpid))){
+ nsigmaTPCK = fESDpid->NumberOfSigmasTPC(esdtrack,AliPID::kKaon);
+ nsigmaTPCPi = fESDpid->NumberOfSigmasTPC(esdtrack,AliPID::kPion);
+ nsigmaTPCP = fESDpid->NumberOfSigmasTPC(esdtrack,AliPID::kProton);
- }
- trackCopy->SetNSigmaTPCPi(nsigmaTPCPi);
- trackCopy->SetNSigmaTPCK(nsigmaTPCK);
- trackCopy->SetNSigmaTPCP(nsigmaTPCP);
+ }
+ trackCopy->SetNSigmaTPCPi(nsigmaTPCPi);
+ trackCopy->SetNSigmaTPCK(nsigmaTPCK);
+ trackCopy->SetNSigmaTPCP(nsigmaTPCP);
- ///// TOF ///////////////////////////////////////////////
+ ///// TOF ///////////////////////////////////////////////
- float vp=-1000.;
- float nsigmaTOFPi=-1000.;
- float nsigmaTOFK=-1000.;
- float nsigmaTOFP=-1000.;
+ float vp=-1000.;
+ float nsigmaTOFPi=-1000.;
+ float nsigmaTOFK=-1000.;
+ float nsigmaTOFP=-1000.;
- if ((esdtrack->GetStatus()&AliESDtrack::kTOFpid) &&
- (esdtrack->GetStatus()&AliESDtrack::kTOFout) &&
- (esdtrack->GetStatus()&AliESDtrack::kTIME))
+ if ((esdtrack->GetStatus()&AliESDtrack::kTOFpid) &&
+ (esdtrack->GetStatus()&AliESDtrack::kTOFout) &&
+ (esdtrack->GetStatus()&AliESDtrack::kTIME))
{
//if ((esdtrack->GetStatus()&AliESDtrack::kTOFpid) &&
//(esdtrack->GetStatus()&AliESDtrack::kTIME)){
// collect info from ESDpid class
- if ((fESDpid) && (esdtrack->IsOn(AliESDtrack::kTOFpid))) {
+ if ((fESDpid) && (esdtrack->IsOn(AliESDtrack::kTOFpid))) {
double tZero = fESDpid->GetTOFResponse().GetStartTime(esdtrack->P());
}
}
- trackCopy->SetVTOF(vp);
- trackCopy->SetNSigmaTOFPi(nsigmaTOFPi);
- trackCopy->SetNSigmaTOFK(nsigmaTOFK);
- trackCopy->SetNSigmaTOFP(nsigmaTOFP);
+ trackCopy->SetVTOF(vp);
+ trackCopy->SetNSigmaTOFPi(nsigmaTOFPi);
+ trackCopy->SetNSigmaTOFK(nsigmaTOFK);
+ trackCopy->SetNSigmaTOFP(nsigmaTOFP);
- double pxyz[3];
- double rxyz[3];
- double impact[2];
- double covimpact[3];
+ double pxyz[3];
+ double rxyz[3];
+ double impact[2];
+ double covimpact[3];
- if (fUseTPCOnly) {
- if (!esdtrack->GetTPCInnerParam()) {
- cout << "No TPC inner param !" << endl;
- delete trackCopy;
- continue;
- }
-
- AliExternalTrackParam *param = new AliExternalTrackParam(*esdtrack->GetTPCInnerParam());
- param->GetXYZ(rxyz);
- param->PropagateToDCA(fEvent->GetPrimaryVertexTPC(), (fEvent->GetMagneticField()), 10000, impact, covimpact);
- param->GetPxPyPz(pxyz);//reading noconstarined momentum
-
- if (fReadInner == true) {
- AliFemtoModelHiddenInfo *tInfo = new AliFemtoModelHiddenInfo();
- tInfo->SetPDGPid(211);
- tInfo->SetTrueMomentum(pxyz[0], pxyz[1], pxyz[2]);
- tInfo->SetMass(0.13957);
- // tInfo->SetEmissionPoint(rxyz[0], rxyz[1], rxyz[2], 0.0);
- // tInfo->SetEmissionPoint(fV1[0], fV1[1], fV1[2], 0.0);
- tInfo->SetEmissionPoint(rxyz[0]-fV1[0], rxyz[1]-fV1[1], rxyz[2]-fV1[2], 0.0);
- trackCopy->SetHiddenInfo(tInfo);
- }
+ if (fUseTPCOnly) {
+ if (!esdtrack->GetTPCInnerParam()) {
+ cout << "No TPC inner param !" << endl;
+ delete trackCopy;
+ continue;
+ }
- if (fRotateToEventPlane) {
- double tPhi = TMath::ATan2(pxyz[1], pxyz[0]);
- double tRad = TMath::Hypot(pxyz[0], pxyz[1]);
-
- pxyz[0] = tRad*TMath::Cos(tPhi - tReactionPlane);
- pxyz[1] = tRad*TMath::Sin(tPhi - tReactionPlane);
- }
-
- AliFemtoThreeVector v(pxyz[0],pxyz[1],pxyz[2]);
- if (v.Mag() < 0.0001) {
- // cout << "Found 0 momentum ???? " << pxyz[0] << " " << pxyz[1] << " " << pxyz[2] << endl;
- delete trackCopy;
- continue;
- }
- trackCopy->SetP(v);//setting momentum
- trackCopy->SetPt(sqrt(pxyz[0]*pxyz[0]+pxyz[1]*pxyz[1]));
-
- const AliFmThreeVectorD kP(pxyz[0],pxyz[1],pxyz[2]);
- const AliFmThreeVectorD kOrigin(fV1[0],fV1[1],fV1[2]);
- //setting helix I do not if it is ok
- AliFmPhysicalHelixD helix(kP,kOrigin,(double)(fEvent->GetMagneticField())*kilogauss,(double)(trackCopy->Charge()));
- trackCopy->SetHelix(helix);
-
- //some stuff which could be useful
- trackCopy->SetImpactD(impact[0]);
- trackCopy->SetImpactZ(impact[1]);
- trackCopy->SetCdd(covimpact[0]);
- trackCopy->SetCdz(covimpact[1]);
- trackCopy->SetCzz(covimpact[2]);
- trackCopy->SetSigmaToVertex(GetSigmaToVertex(impact, covimpact));
-
- delete param;
+ AliExternalTrackParam *param = new AliExternalTrackParam(*esdtrack->GetTPCInnerParam());
+ param->GetXYZ(rxyz);
+ param->PropagateToDCA(fEvent->GetPrimaryVertexTPC(), (fEvent->GetMagneticField()), 10000, impact, covimpact);
+ param->GetPxPyPz(pxyz);//reading noconstarined momentum
+
+ if (fReadInner == true) {
+ AliFemtoModelHiddenInfo *tInfo = new AliFemtoModelHiddenInfo();
+ tInfo->SetPDGPid(211);
+ tInfo->SetTrueMomentum(pxyz[0], pxyz[1], pxyz[2]);
+ tInfo->SetMass(0.13957);
+ // tInfo->SetEmissionPoint(rxyz[0], rxyz[1], rxyz[2], 0.0);
+ // tInfo->SetEmissionPoint(fV1[0], fV1[1], fV1[2], 0.0);
+ tInfo->SetEmissionPoint(rxyz[0]-fV1[0], rxyz[1]-fV1[1], rxyz[2]-fV1[2], 0.0);
+ trackCopy->SetHiddenInfo(tInfo);
}
- else {
- if (fReadInner == true) {
-
- if (esdtrack->GetTPCInnerParam()) {
- AliExternalTrackParam *param = new AliExternalTrackParam(*esdtrack->GetTPCInnerParam());
- trackCopy->SetInnerMomentum(esdtrack->GetTPCmomentum());
- param->GetXYZ(rxyz);
- // param->PropagateToDCA(fEvent->GetPrimaryVertex(), (fEvent->GetMagneticField()), 10000);
- param->GetPxPyPz(pxyz);//reading noconstarined momentum
- delete param;
-
- AliFemtoModelHiddenInfo *tInfo = new AliFemtoModelHiddenInfo();
- tInfo->SetPDGPid(211);
- tInfo->SetTrueMomentum(pxyz[0], pxyz[1], pxyz[2]);
- tInfo->SetMass(0.13957);
- // tInfo->SetEmissionPoint(rxyz[0], rxyz[1], rxyz[2], 0.0);
- //tInfo->SetEmissionPoint(fV1[0], fV1[1], fV1[2], 0.0);
- tInfo->SetEmissionPoint(rxyz[0]-fV1[0], rxyz[1]-fV1[1], rxyz[2]-fV1[2], 0.0);
- trackCopy->SetHiddenInfo(tInfo);
- }
- }
+ if (fRotateToEventPlane) {
+ double tPhi = TMath::ATan2(pxyz[1], pxyz[0]);
+ double tRad = TMath::Hypot(pxyz[0], pxyz[1]);
- if(fTrackType == kGlobal) {
- if (fConstrained==true)
- tGoodMomentum=esdtrack->GetConstrainedPxPyPz(pxyz); //reading constrained momentum
- else
- tGoodMomentum=esdtrack->GetPxPyPz(pxyz);//reading noconstarined momentum
- }
- else if (fTrackType == kTPCOnly) {
- if (esdtrack->GetTPCInnerParam())
- esdtrack->GetTPCInnerParam()->GetPxPyPz(pxyz);
- else {
- delete trackCopy;
- continue;
- }
- }
- else if (fTrackType == kITSOnly) {
- if (fConstrained==true)
- tGoodMomentum=esdtrack->GetConstrainedPxPyPz(pxyz); //reading constrained momentum
- else
- tGoodMomentum=esdtrack->GetPxPyPz(pxyz);//reading noconstarined momentum
- }
+ pxyz[0] = tRad*TMath::Cos(tPhi - tReactionPlane);
+ pxyz[1] = tRad*TMath::Sin(tPhi - tReactionPlane);
+ }
- if (fRotateToEventPlane) {
- double tPhi = TMath::ATan2(pxyz[1], pxyz[0]);
- double tRad = TMath::Hypot(pxyz[0], pxyz[1]);
-
- pxyz[0] = tRad*TMath::Cos(tPhi - tReactionPlane);
- pxyz[1] = tRad*TMath::Sin(tPhi - tReactionPlane);
- }
-
- AliFemtoThreeVector v(pxyz[0],pxyz[1],pxyz[2]);
- if (v.Mag() < 0.0001) {
- // cout << "Found 0 momentum ???? " << pxyz[0] << " " << pxyz[1] << " " << pxyz[2] << endl;
- delete trackCopy;
- continue;
- }
-
- trackCopy->SetP(v);//setting momentum
- trackCopy->SetPt(sqrt(pxyz[0]*pxyz[0]+pxyz[1]*pxyz[1]));
- const AliFmThreeVectorD kP(pxyz[0],pxyz[1],pxyz[2]);
- const AliFmThreeVectorD kOrigin(fV1[0],fV1[1],fV1[2]);
- //setting helix I do not if it is ok
- AliFmPhysicalHelixD helix(kP,kOrigin,(double)(fEvent->GetMagneticField())*kilogauss,(double)(trackCopy->Charge()));
- trackCopy->SetHelix(helix);
-
- //some stuff which could be useful
- float imp[2];
- float cim[3];
- esdtrack->GetImpactParameters(imp,cim);
-
- impact[0] = imp[0];
- impact[1] = imp[1];
- covimpact[0] = cim[0];
- covimpact[1] = cim[1];
- covimpact[2] = cim[2];
-
- trackCopy->SetImpactD(impact[0]);
- trackCopy->SetImpactZ(impact[1]);
- trackCopy->SetCdd(covimpact[0]);
- trackCopy->SetCdz(covimpact[1]);
- trackCopy->SetCzz(covimpact[2]);
- trackCopy->SetSigmaToVertex(GetSigmaToVertex(impact,covimpact));
+ AliFemtoThreeVector v(pxyz[0],pxyz[1],pxyz[2]);
+ if (v.Mag() < 0.0001) {
+ // cout << "Found 0 momentum ???? " << pxyz[0] << " " << pxyz[1] << " " << pxyz[2] << endl;
+ delete trackCopy;
+ continue;
+ }
+ trackCopy->SetP(v);//setting momentum
+ trackCopy->SetPt(sqrt(pxyz[0]*pxyz[0]+pxyz[1]*pxyz[1]));
+
+ const AliFmThreeVectorD kP(pxyz[0],pxyz[1],pxyz[2]);
+ const AliFmThreeVectorD kOrigin(fV1[0],fV1[1],fV1[2]);
+ //setting helix I do not if it is ok
+ AliFmPhysicalHelixD helix(kP,kOrigin,(double)(fEvent->GetMagneticField())*kilogauss,(double)(trackCopy->Charge()));
+ trackCopy->SetHelix(helix);
+
+ //some stuff which could be useful
+ trackCopy->SetImpactD(impact[0]);
+ trackCopy->SetImpactZ(impact[1]);
+ trackCopy->SetCdd(covimpact[0]);
+ trackCopy->SetCdz(covimpact[1]);
+ trackCopy->SetCzz(covimpact[2]);
+ trackCopy->SetSigmaToVertex(GetSigmaToVertex(impact, covimpact));
+
+ delete param;
+ }
+ else {
+ if (fReadInner == true) {
+
+ if (esdtrack->GetTPCInnerParam()) {
+ AliExternalTrackParam *param = new AliExternalTrackParam(*esdtrack->GetTPCInnerParam());
+ trackCopy->SetInnerMomentum(esdtrack->GetTPCmomentum());
+ param->GetXYZ(rxyz);
+ // param->PropagateToDCA(fEvent->GetPrimaryVertex(), (fEvent->GetMagneticField()), 10000);
+ param->GetPxPyPz(pxyz);//reading noconstarined momentum
+ delete param;
+
+ AliFemtoModelHiddenInfo *tInfo = new AliFemtoModelHiddenInfo();
+ tInfo->SetPDGPid(211);
+ tInfo->SetTrueMomentum(pxyz[0], pxyz[1], pxyz[2]);
+ tInfo->SetMass(0.13957);
+ // tInfo->SetEmissionPoint(rxyz[0], rxyz[1], rxyz[2], 0.0);
+ //tInfo->SetEmissionPoint(fV1[0], fV1[1], fV1[2], 0.0);
+ tInfo->SetEmissionPoint(rxyz[0]-fV1[0], rxyz[1]-fV1[1], rxyz[2]-fV1[2], 0.0);
+ trackCopy->SetHiddenInfo(tInfo);
+ }
}
- trackCopy->SetTrackId(esdtrack->GetID());
- trackCopy->SetFlags(esdtrack->GetStatus());
- trackCopy->SetLabel(esdtrack->GetLabel());
- trackCopy->SetITSchi2(esdtrack->GetITSchi2());
- trackCopy->SetITSncls(esdtrack->GetNcls(0));
- trackCopy->SetTPCchi2(esdtrack->GetTPCchi2());
- trackCopy->SetTPCncls(esdtrack->GetTPCNcls());
- trackCopy->SetTPCnclsF(esdtrack->GetTPCNclsF());
- trackCopy->SetTPCsignalN((short)esdtrack->GetTPCsignalN()); //due to bug in aliesdtrack class
- trackCopy->SetTPCsignalS(esdtrack->GetTPCsignalSigma());
+ if(fTrackType == kGlobal) {
+ if (fConstrained==true)
+ tGoodMomentum=esdtrack->GetConstrainedPxPyPz(pxyz); //reading constrained momentum
+ else
+ tGoodMomentum=esdtrack->GetPxPyPz(pxyz);//reading noconstarined momentum
+ }
+ else if (fTrackType == kTPCOnly) {
+ if (esdtrack->GetTPCInnerParam())
+ esdtrack->GetTPCInnerParam()->GetPxPyPz(pxyz);
+ else {
+ delete trackCopy;
+ continue;
+ }
+ }
+ else if (fTrackType == kITSOnly) {
+ if (fConstrained==true)
+ tGoodMomentum=esdtrack->GetConstrainedPxPyPz(pxyz); //reading constrained momentum
+ else
+ tGoodMomentum=esdtrack->GetPxPyPz(pxyz);//reading noconstarined momentum
+ }
+
+ if (fRotateToEventPlane) {
+ double tPhi = TMath::ATan2(pxyz[1], pxyz[0]);
+ double tRad = TMath::Hypot(pxyz[0], pxyz[1]);
+ pxyz[0] = tRad*TMath::Cos(tPhi - tReactionPlane);
+ pxyz[1] = tRad*TMath::Sin(tPhi - tReactionPlane);
+ }
- trackCopy->SetTPCClusterMap(esdtrack->GetTPCClusterMap());
- trackCopy->SetTPCSharedMap(esdtrack->GetTPCSharedMap());
+ AliFemtoThreeVector v(pxyz[0],pxyz[1],pxyz[2]);
+ if (v.Mag() < 0.0001) {
+ // cout << "Found 0 momentum ???? " << pxyz[0] << " " << pxyz[1] << " " << pxyz[2] << endl;
+ delete trackCopy;
+ continue;
+ }
- double xtpc[3];
- esdtrack->GetInnerXYZ(xtpc);
- xtpc[2] -= fV1[2];
- trackCopy->SetNominalTPCEntrancePoint(xtpc);
+ trackCopy->SetP(v);//setting momentum
+ trackCopy->SetPt(sqrt(pxyz[0]*pxyz[0]+pxyz[1]*pxyz[1]));
+ const AliFmThreeVectorD kP(pxyz[0],pxyz[1],pxyz[2]);
+ const AliFmThreeVectorD kOrigin(fV1[0],fV1[1],fV1[2]);
+ //setting helix I do not if it is ok
+ AliFmPhysicalHelixD helix(kP,kOrigin,(double)(fEvent->GetMagneticField())*kilogauss,(double)(trackCopy->Charge()));
+ trackCopy->SetHelix(helix);
+
+ //some stuff which could be useful
+ float imp[2];
+ float cim[3];
+ esdtrack->GetImpactParameters(imp,cim);
+
+ impact[0] = imp[0];
+ impact[1] = imp[1];
+ covimpact[0] = cim[0];
+ covimpact[1] = cim[1];
+ covimpact[2] = cim[2];
+
+ trackCopy->SetImpactD(impact[0]);
+ trackCopy->SetImpactZ(impact[1]);
+ trackCopy->SetCdd(covimpact[0]);
+ trackCopy->SetCdz(covimpact[1]);
+ trackCopy->SetCzz(covimpact[2]);
+ trackCopy->SetSigmaToVertex(GetSigmaToVertex(impact,covimpact));
+ }
- esdtrack->GetOuterXYZ(xtpc);
- xtpc[2] -= fV1[2];
- trackCopy->SetNominalTPCExitPoint(xtpc);
+ trackCopy->SetTrackId(esdtrack->GetID());
+ trackCopy->SetFlags(esdtrack->GetStatus());
+ trackCopy->SetLabel(esdtrack->GetLabel());
- int indexes[3];
- for (int ik=0; ik<3; ik++) {
- indexes[ik] = esdtrack->GetKinkIndex(ik);
- }
- trackCopy->SetKinkIndexes(indexes);
+ trackCopy->SetITSchi2(esdtrack->GetITSchi2());
+ trackCopy->SetITSncls(esdtrack->GetNcls(0));
+ trackCopy->SetTPCchi2(esdtrack->GetTPCchi2());
+ trackCopy->SetTPCncls(esdtrack->GetTPCNcls());
+ trackCopy->SetTPCnclsF(esdtrack->GetTPCNclsF());
+ trackCopy->SetTPCsignalN((short)esdtrack->GetTPCsignalN()); //due to bug in aliesdtrack class
+ trackCopy->SetTPCsignalS(esdtrack->GetTPCsignalSigma());
+ trackCopy->SetTPCsignal(esdtrack->GetTPCsignal());
- for (int ii=0; ii<6; ii++){
- trackCopy->SetITSHitOnLayer(ii,esdtrack->HasPointOnITSLayer(ii));
- }
+ cout << "TPC signal = " << esdtrack->GetTPCsignal() << endl;
+ trackCopy->SetTPCClusterMap(esdtrack->GetTPCClusterMap());
+ trackCopy->SetTPCSharedMap(esdtrack->GetTPCSharedMap());
- // Fill the hidden information with the simulated data
- if (TMath::Abs(esdtrack->GetLabel()) < fStack->GetNtrack()) {
- TParticle *tPart = fStack->Particle(TMath::Abs(esdtrack->GetLabel()));
+ double xtpc[3];
+ esdtrack->GetInnerXYZ(xtpc);
+ xtpc[2] -= fV1[2];
+ trackCopy->SetNominalTPCEntrancePoint(xtpc);
- // Check the mother information
+ esdtrack->GetOuterXYZ(xtpc);
+ xtpc[2] -= fV1[2];
+ trackCopy->SetNominalTPCExitPoint(xtpc);
- // Using the new way of storing the freeze-out information
- // Final state particle is stored twice on the stack
- // one copy (mother) is stored with original freeze-out information
- // and is not tracked
- // the other one (daughter) is stored with primary vertex position
- // and is tracked
+ int indexes[3];
+ for (int ik=0; ik<3; ik++) {
+ indexes[ik] = esdtrack->GetKinkIndex(ik);
+ }
+ trackCopy->SetKinkIndexes(indexes);
- // Freeze-out coordinates
- double fpx=0.0, fpy=0.0, fpz=0.0, fpt=0.0;
- fpx = tPart->Vx() - fV1[0];
- fpy = tPart->Vy() - fV1[1];
- fpz = tPart->Vz() - fV1[2];
- fpt = tPart->T();
+ for (int ii=0; ii<6; ii++){
+ trackCopy->SetITSHitOnLayer(ii,esdtrack->HasPointOnITSLayer(ii));
+ }
- AliFemtoModelGlobalHiddenInfo *tInfo = new AliFemtoModelGlobalHiddenInfo();
- tInfo->SetGlobalEmissionPoint(fpx, fpy, fpz);
- trackCopy->SetGlobalEmissionPoint(fpx, fpy, fpz);
- fpx *= 1e13;
- fpy *= 1e13;
- fpz *= 1e13;
- fpt *= 1e13;
+ // Fill the hidden information with the simulated data
+ if (TMath::Abs(esdtrack->GetLabel()) < fStack->GetNtrack()) {
+ TParticle *tPart = fStack->Particle(TMath::Abs(esdtrack->GetLabel()));
+ // Check the mother information
- if (fOnlyPrimaries && !fStack->IsPhysicalPrimary(TMath::Abs(esdtrack->GetLabel()))){
- delete trackCopy;
- continue;
- }
+ // Using the new way of storing the freeze-out information
+ // Final state particle is stored twice on the stack
+ // one copy (mother) is stored with original freeze-out information
+ // and is not tracked
+ // the other one (daughter) is stored with primary vertex position
+ // and is tracked
- // fillDCA
+ // Freeze-out coordinates
+ double fpx=0.0, fpy=0.0, fpz=0.0, fpt=0.0;
+ fpx = tPart->Vx() - fV1[0];
+ fpy = tPart->Vy() - fV1[1];
+ fpz = tPart->Vz() - fV1[2];
+ fpt = tPart->T();
- if (TMath::Abs(impact[0]) > 0.001) {
- if (fStack->IsPhysicalPrimary(TMath::Abs(esdtrack->GetLabel()))){
- trackCopy->SetImpactDprim(impact[0]);
- //cout << "prim" << endl;
+ AliFemtoModelGlobalHiddenInfo *tInfo = new AliFemtoModelGlobalHiddenInfo();
+ tInfo->SetGlobalEmissionPoint(fpx, fpy, fpz);
+ trackCopy->SetGlobalEmissionPoint(fpx, fpy, fpz);
- }
- else if (fStack->IsSecondaryFromWeakDecay(TMath::Abs(esdtrack->GetLabel()))) {
- trackCopy->SetImpactDweak(impact[0]);
- //cout << "wea" << endl;
- }
- else if (fStack->IsSecondaryFromMaterial(TMath::Abs(esdtrack->GetLabel()))) {
- trackCopy->SetImpactDmat(impact[0]);
- //cout << "mat" << endl;
- }
- }
- // end fillDCA
-
- // cout << "Looking for mother ids " << endl;
- if (motherids[TMath::Abs(esdtrack->GetLabel())]>0) {
- // cout << "Got mother id" << endl;
- TParticle *mother = fStack->Particle(motherids[TMath::Abs(esdtrack->GetLabel())]);
- // Check if this is the same particle stored twice on the stack
- if ((mother->GetPdgCode() == tPart->GetPdgCode() || (mother->Px() == tPart->Px()))) {
- // It is the same particle
- // Read in the original freeze-out information
- // and convert it from to [fm]
-
- // EPOS style
- fpx = mother->Vx()*1e13*0.197327;
- fpy = mother->Vy()*1e13*0.197327;
- fpz = mother->Vz()*1e13*0.197327;
- fpt = mother->T() *1e13*0.197327;
-
-
- // Therminator style
+ fpx *= 1e13;
+ fpy *= 1e13;
+ fpz *= 1e13;
+ fpt *= 1e13;
+
+
+ if (fOnlyPrimaries && !fStack->IsPhysicalPrimary(TMath::Abs(esdtrack->GetLabel()))){
+ delete trackCopy;
+ continue;
+ }
+
+ // fillDCA
+
+ if (TMath::Abs(impact[0]) > 0.001) {
+ if (fStack->IsPhysicalPrimary(TMath::Abs(esdtrack->GetLabel()))){
+ trackCopy->SetImpactDprim(impact[0]);
+ //cout << "prim" << endl;
+
+ }
+ else if (fStack->IsSecondaryFromWeakDecay(TMath::Abs(esdtrack->GetLabel()))) {
+ trackCopy->SetImpactDweak(impact[0]);
+ //cout << "wea" << endl;
+ }
+ else if (fStack->IsSecondaryFromMaterial(TMath::Abs(esdtrack->GetLabel()))) {
+ trackCopy->SetImpactDmat(impact[0]);
+ //cout << "mat" << endl;
+ }
+ }
+ // end fillDCA
+
+ // cout << "Looking for mother ids " << endl;
+ if (motherids[TMath::Abs(esdtrack->GetLabel())]>0) {
+ // cout << "Got mother id" << endl;
+ TParticle *mother = fStack->Particle(motherids[TMath::Abs(esdtrack->GetLabel())]);
+ // Check if this is the same particle stored twice on the stack
+ if ((mother->GetPdgCode() == tPart->GetPdgCode() || (mother->Px() == tPart->Px()))) {
+ // It is the same particle
+ // Read in the original freeze-out information
+ // and convert it from to [fm]
+
+ // EPOS style
+ fpx = mother->Vx()*1e13*0.197327;
+ fpy = mother->Vy()*1e13*0.197327;
+ fpz = mother->Vz()*1e13*0.197327;
+ fpt = mother->T() *1e13*0.197327;
+
+
+ // Therminator style
// fpx = mother->Vx()*1e13;
// fpy = mother->Vy()*1e13;
// fpz = mother->Vz()*1e13;
// fpt = mother->T() *1e13*3e10;
- }
- }
-
- if (fRotateToEventPlane) {
- double tPhi = TMath::ATan2(fpy, fpx);
- double tRad = TMath::Hypot(fpx, fpy);
-
- fpx = tRad*TMath::Cos(tPhi - tReactionPlane);
- fpy = tRad*TMath::Sin(tPhi - tReactionPlane);
- }
-
- tInfo->SetPDGPid(tPart->GetPdgCode());
- trackCopy->SetPDGPid(tPart->GetPdgCode());
-
- if (fRotateToEventPlane) {
- double tPhi = TMath::ATan2(tPart->Py(), tPart->Px());
- double tRad = TMath::Hypot(tPart->Px(), tPart->Py());
-
- tInfo->SetTrueMomentum(tRad*TMath::Cos(tPhi - tReactionPlane),
- tRad*TMath::Sin(tPhi - tReactionPlane),
- tPart->Pz());
- trackCopy->SetTrueMomentum(tRad*TMath::Cos(tPhi - tReactionPlane),
- tRad*TMath::Sin(tPhi - tReactionPlane),
- tPart->Pz());
- }
- else
- {
- tInfo->SetTrueMomentum(tPart->Px(), tPart->Py(), tPart->Pz());
- trackCopy->SetTrueMomentum(tPart->Px(), tPart->Py(), tPart->Pz());
- }
- Double_t mass2 = (tPart->Energy() *tPart->Energy() -
- tPart->Px()*tPart->Px() -
- tPart->Py()*tPart->Py() -
- tPart->Pz()*tPart->Pz());
- if (mass2>0.0)
- {
- tInfo->SetMass(TMath::Sqrt(mass2));
- trackCopy->SetMass(TMath::Sqrt(mass2));
- }
- else
- {
- tInfo->SetMass(0.0);
- trackCopy->SetMass(0.0);
- }
-
- tInfo->SetEmissionPoint(fpx, fpy, fpz, fpt);
- trackCopy->SetEmissionPoint(fpx, fpy, fpz, fpt);
- trackCopy->SetHiddenInfo(tInfo);
- }
- else {
- AliFemtoModelGlobalHiddenInfo *tInfo = new AliFemtoModelGlobalHiddenInfo();
- tInfo->SetMass(0.0);
- double fpx=0.0, fpy=0.0, fpz=0.0, fpt=0.0;
- fpx = fV1[0]*1e13;
- fpy = fV1[1]*1e13;
- fpz = fV1[2]*1e13;
- fpt = 0.0;
- tInfo->SetEmissionPoint(fpx, fpy, fpz, fpt);
-
- tInfo->SetTrueMomentum(pxyz[0],pxyz[1],pxyz[2]);
-
- trackCopy->SetHiddenInfo(tInfo);
+ }
}
- // cout << "Got freeze-out " << fpx << " " << fpy << " " << fpz << " " << fpt << " " << mass2 << " " << tPart->GetPdgCode() << endl;
- //decision if we want this track
- //if we using diffrent labels we want that this label was use for first time
- //if we use hidden info we want to have match between sim data and ESD
+ if (fRotateToEventPlane) {
+ double tPhi = TMath::ATan2(fpy, fpx);
+ double tRad = TMath::Hypot(fpx, fpy);
- bool trackAccept = true;
- if (isKaonAnalysis == true && trackCopy->GetPDGPid() != 321) {
- trackAccept = false;
+ fpx = tRad*TMath::Cos(tPhi - tReactionPlane);
+ fpy = tRad*TMath::Sin(tPhi - tReactionPlane);
}
- if (tGoodMomentum==true && trackAccept == true)
- {
+ tInfo->SetPDGPid(tPart->GetPdgCode());
+ trackCopy->SetPDGPid(tPart->GetPdgCode());
- hbtEvent->TrackCollection()->push_back(trackCopy);//adding track to analysis
- realnofTracks++;//real number of tracks
- // delete trackCopy;
- }
+ if (fRotateToEventPlane) {
+ double tPhi = TMath::ATan2(tPart->Py(), tPart->Px());
+ double tRad = TMath::Hypot(tPart->Px(), tPart->Py());
+
+ tInfo->SetTrueMomentum(tRad*TMath::Cos(tPhi - tReactionPlane),
+ tRad*TMath::Sin(tPhi - tReactionPlane),
+ tPart->Pz());
+ trackCopy->SetTrueMomentum(tRad*TMath::Cos(tPhi - tReactionPlane),
+ tRad*TMath::Sin(tPhi - tReactionPlane),
+ tPart->Pz());
+ }
+ else
+ {
+ tInfo->SetTrueMomentum(tPart->Px(), tPart->Py(), tPart->Pz());
+ trackCopy->SetTrueMomentum(tPart->Px(), tPart->Py(), tPart->Pz());
+ }
+ Double_t mass2 = (tPart->Energy() *tPart->Energy() -
+ tPart->Px()*tPart->Px() -
+ tPart->Py()*tPart->Py() -
+ tPart->Pz()*tPart->Pz());
+ if (mass2>0.0)
+ {
+ tInfo->SetMass(TMath::Sqrt(mass2));
+ trackCopy->SetMass(TMath::Sqrt(mass2));
+ }
else
- {
- delete trackCopy;
- }
+ {
+ tInfo->SetMass(0.0);
+ trackCopy->SetMass(0.0);
+ }
+
+ tInfo->SetEmissionPoint(fpx, fpy, fpz, fpt);
+ trackCopy->SetEmissionPoint(fpx, fpy, fpz, fpt);
+ trackCopy->SetHiddenInfo(tInfo);
+ }
+ else {
+ AliFemtoModelGlobalHiddenInfo *tInfo = new AliFemtoModelGlobalHiddenInfo();
+ tInfo->SetMass(0.0);
+ double fpx=0.0, fpy=0.0, fpz=0.0, fpt=0.0;
+ fpx = fV1[0]*1e13;
+ fpy = fV1[1]*1e13;
+ fpz = fV1[2]*1e13;
+ fpt = 0.0;
+ tInfo->SetEmissionPoint(fpx, fpy, fpz, fpt);
+
+ tInfo->SetTrueMomentum(pxyz[0],pxyz[1],pxyz[2]);
+
+ trackCopy->SetHiddenInfo(tInfo);
+ }
+ // cout << "Got freeze-out " << fpx << " " << fpy << " " << fpz << " " << fpt << " " << mass2 << " " << tPart->GetPdgCode() << endl;
+
+ //decision if we want this track
+ //if we using diffrent labels we want that this label was use for first time
+ //if we use hidden info we want to have match between sim data and ESD
+
+ bool trackAccept = true;
+ if (isKaonAnalysis == true && trackCopy->GetPDGPid() != 321) {
+ trackAccept = false;
+ }
+ if (tGoodMomentum==true && trackAccept == true)
+ {
+
+ hbtEvent->TrackCollection()->push_back(trackCopy);//adding track to analysis
+ realnofTracks++;//real number of tracks
+ // delete trackCopy;
+ }
+ else
+ {
+ delete trackCopy;
}
+ }
+
delete [] motherids;
hbtEvent->SetNumberOfTracks(realnofTracks);//setting number of track which we read in event
else if (fEstEventMult == kSPDLayer1)
hbtEvent->SetNormalizedMult(fEvent->GetMultiplicity()->GetNumberOfITSClusters(1));
else if (fEstEventMult == kVZERO)
- {
- Float_t multV0 = 0;
- for (Int_t i=0; i<64; i++)
- multV0 += fEvent->GetVZEROData()->GetMultiplicity(i);
- hbtEvent->SetNormalizedMult(multV0);
- }
+ {
+ Float_t multV0 = 0;
+ for (Int_t i=0; i<64; i++)
+ multV0 += fEvent->GetVZEROData()->GetMultiplicity(i);
+ hbtEvent->SetNormalizedMult(multV0);
+ }
else if (fEstEventMult == kCentrality) {
// centrality between 0 (central) and 1 (very peripheral)
if (cent) {
if (cent->GetCentralityPercentile("V0M") < 0.00001)
- hbtEvent->SetNormalizedMult(-1);
+ hbtEvent->SetNormalizedMult(-1);
else
- hbtEvent->SetNormalizedMult(lrint(10.0*cent->GetCentralityPercentile("V0M")));
+ hbtEvent->SetNormalizedMult(lrint(10.0*cent->GetCentralityPercentile("V0M")));
if (Debug()>1) printf ("Set Centrality %i %f %li\n", hbtEvent->UncorrectedNumberOfPrimaries(),
- 10.0*cent->GetCentralityPercentile("V0M"), lrint(10.0*cent->GetCentralityPercentile("V0M")));
+ 10.0*cent->GetCentralityPercentile("V0M"), lrint(10.0*cent->GetCentralityPercentile("V0M")));
}
}
else if (fEstEventMult == kCentralityV0A) {
if (cent) {
if (cent->GetCentralityPercentile("V0A") < 0.00001)
- hbtEvent->SetNormalizedMult(-1);
+ hbtEvent->SetNormalizedMult(-1);
else
- hbtEvent->SetNormalizedMult(lrint(10.0*cent->GetCentralityPercentile("V0A")));
+ hbtEvent->SetNormalizedMult(lrint(10.0*cent->GetCentralityPercentile("V0A")));
if (Debug()>1) printf ("Set Centrality %i %f %li\n", hbtEvent->UncorrectedNumberOfPrimaries(),
- 10.0*cent->GetCentralityPercentile("V0A"), lrint(10.0*cent->GetCentralityPercentile("V0A")));
+ 10.0*cent->GetCentralityPercentile("V0A"), lrint(10.0*cent->GetCentralityPercentile("V0A")));
}
}
else if (fEstEventMult == kCentralityV0C) {
if (cent) {
if (cent->GetCentralityPercentile("V0C") < 0.00001)
- hbtEvent->SetNormalizedMult(-1);
+ hbtEvent->SetNormalizedMult(-1);
else
- hbtEvent->SetNormalizedMult(lrint(10.0*cent->GetCentralityPercentile("V0C")));
+ hbtEvent->SetNormalizedMult(lrint(10.0*cent->GetCentralityPercentile("V0C")));
if (Debug()>1) printf ("Set Centrality %i %f %li\n", hbtEvent->UncorrectedNumberOfPrimaries(),
- 10.0*cent->GetCentralityPercentile("V0C"), lrint(10.0*cent->GetCentralityPercentile("V0C")));
+ 10.0*cent->GetCentralityPercentile("V0C"), lrint(10.0*cent->GetCentralityPercentile("V0C")));
}
}
else if (fEstEventMult == kCentralityZNA) {
if (cent) {
if (cent->GetCentralityPercentile("ZNA") < 0.00001)
- hbtEvent->SetNormalizedMult(-1);
+ hbtEvent->SetNormalizedMult(-1);
else
- hbtEvent->SetNormalizedMult(lrint(10.0*cent->GetCentralityPercentile("ZNA")));
+ hbtEvent->SetNormalizedMult(lrint(10.0*cent->GetCentralityPercentile("ZNA")));
if (Debug()>1) printf ("Set Centrality %i %f %li\n", hbtEvent->UncorrectedNumberOfPrimaries(),
- 10.0*cent->GetCentralityPercentile("ZNA"), lrint(10.0*cent->GetCentralityPercentile("ZNA")));
+ 10.0*cent->GetCentralityPercentile("ZNA"), lrint(10.0*cent->GetCentralityPercentile("ZNA")));
}
}
else if (fEstEventMult == kCentralityZNC) {
if (cent) {
if (cent->GetCentralityPercentile("ZNC") < 0.00001)
- hbtEvent->SetNormalizedMult(-1);
+ hbtEvent->SetNormalizedMult(-1);
else
- hbtEvent->SetNormalizedMult(lrint(10.0*cent->GetCentralityPercentile("ZNC")));
+ hbtEvent->SetNormalizedMult(lrint(10.0*cent->GetCentralityPercentile("ZNC")));
if (Debug()>1) printf ("Set Centrality %i %f %li\n", hbtEvent->UncorrectedNumberOfPrimaries(),
- 10.0*cent->GetCentralityPercentile("ZNC"), lrint(10.0*cent->GetCentralityPercentile("ZNC")));
+ 10.0*cent->GetCentralityPercentile("ZNC"), lrint(10.0*cent->GetCentralityPercentile("ZNC")));
}
}
else if (fEstEventMult == kCentralityCL1) {
if (cent) {
if (cent->GetCentralityPercentile("CL1") < 0.00001)
- hbtEvent->SetNormalizedMult(-1);
+ hbtEvent->SetNormalizedMult(-1);
else
- hbtEvent->SetNormalizedMult(lrint(10.0*cent->GetCentralityPercentile("CL1")));
+ hbtEvent->SetNormalizedMult(lrint(10.0*cent->GetCentralityPercentile("CL1")));
if (Debug()>1) printf ("Set Centrality %i %f %li\n", hbtEvent->UncorrectedNumberOfPrimaries(),
- 10.0*cent->GetCentralityPercentile("CL1"), lrint(10.0*cent->GetCentralityPercentile("CL1")));
+ 10.0*cent->GetCentralityPercentile("CL1"), lrint(10.0*cent->GetCentralityPercentile("CL1")));
}
}
else if (fEstEventMult == kCentralityCL0) {
if (cent) {
if (cent->GetCentralityPercentile("CL0") < 0.00001)
- hbtEvent->SetNormalizedMult(-1);
+ hbtEvent->SetNormalizedMult(-1);
else
- hbtEvent->SetNormalizedMult(lrint(10.0*cent->GetCentralityPercentile("CL0")));
+ hbtEvent->SetNormalizedMult(lrint(10.0*cent->GetCentralityPercentile("CL0")));
if (Debug()>1) printf ("Set Centrality %i %f %li\n", hbtEvent->UncorrectedNumberOfPrimaries(),
- 10.0*cent->GetCentralityPercentile("CL0"), lrint(10.0*cent->GetCentralityPercentile("CL0")));
+ 10.0*cent->GetCentralityPercentile("CL0"), lrint(10.0*cent->GetCentralityPercentile("CL0")));
}
}
else if (fEstEventMult == kCentralityTRK) {
if (cent) {
if (cent->GetCentralityPercentile("TRK") < 0.00001)
- hbtEvent->SetNormalizedMult(-1);
+ hbtEvent->SetNormalizedMult(-1);
else
- hbtEvent->SetNormalizedMult(lrint(10.0*cent->GetCentralityPercentile("TRK")));
+ hbtEvent->SetNormalizedMult(lrint(10.0*cent->GetCentralityPercentile("TRK")));
if (Debug()>1) printf ("Set Centrality %i %f %li\n", hbtEvent->UncorrectedNumberOfPrimaries(),
- 10.0*cent->GetCentralityPercentile("TRK"), lrint(10.0*cent->GetCentralityPercentile("TRK")));
+ 10.0*cent->GetCentralityPercentile("TRK"), lrint(10.0*cent->GetCentralityPercentile("TRK")));
}
}
else if (fEstEventMult == kCentralityTKL) {
if (cent) {
if (cent->GetCentralityPercentile("TKL") < 0.00001)
- hbtEvent->SetNormalizedMult(-1);
+ hbtEvent->SetNormalizedMult(-1);
else
- hbtEvent->SetNormalizedMult(lrint(10.0*cent->GetCentralityPercentile("TKL")));
+ hbtEvent->SetNormalizedMult(lrint(10.0*cent->GetCentralityPercentile("TKL")));
if (Debug()>1) printf ("Set Centrality %i %f %li\n", hbtEvent->UncorrectedNumberOfPrimaries(),
- 10.0*cent->GetCentralityPercentile("TKL"), lrint(10.0*cent->GetCentralityPercentile("TKL")));
+ 10.0*cent->GetCentralityPercentile("TKL"), lrint(10.0*cent->GetCentralityPercentile("TKL")));
}
}
else if (fEstEventMult == kCentralityNPA) {
if (cent) {
if (cent->GetCentralityPercentile("NPA") < 0.00001)
- hbtEvent->SetNormalizedMult(-1);
+ hbtEvent->SetNormalizedMult(-1);
else
- hbtEvent->SetNormalizedMult(lrint(10.0*cent->GetCentralityPercentile("NPA")));
+ hbtEvent->SetNormalizedMult(lrint(10.0*cent->GetCentralityPercentile("NPA")));
if (Debug()>1) printf ("Set Centrality %i %f %li\n", hbtEvent->UncorrectedNumberOfPrimaries(),
- 10.0*cent->GetCentralityPercentile("NPA"), lrint(10.0*cent->GetCentralityPercentile("NPA")));
+ 10.0*cent->GetCentralityPercentile("NPA"), lrint(10.0*cent->GetCentralityPercentile("NPA")));
}
}
else if (fEstEventMult == kCentralityCND) {
if (cent) {
if (cent->GetCentralityPercentile("CND") < 0.00001)
- hbtEvent->SetNormalizedMult(-1);
+ hbtEvent->SetNormalizedMult(-1);
else
- hbtEvent->SetNormalizedMult(lrint(10.0*cent->GetCentralityPercentile("CND")));
+ hbtEvent->SetNormalizedMult(lrint(10.0*cent->GetCentralityPercentile("CND")));
if (Debug()>1) printf ("Set Centrality %i %f %li\n", hbtEvent->UncorrectedNumberOfPrimaries(),
- 10.0*cent->GetCentralityPercentile("CND"), lrint(10.0*cent->GetCentralityPercentile("CND")));
+ 10.0*cent->GetCentralityPercentile("CND"), lrint(10.0*cent->GetCentralityPercentile("CND")));
}
}
else if (fEstEventMult == kCentralityFMD) {
if (cent) {
if (cent->GetCentralityPercentile("FMD") < 0.00001)
- hbtEvent->SetNormalizedMult(-1);
+ hbtEvent->SetNormalizedMult(-1);
else
- hbtEvent->SetNormalizedMult(lrint(10.0*cent->GetCentralityPercentile("FMD")));
+ hbtEvent->SetNormalizedMult(lrint(10.0*cent->GetCentralityPercentile("FMD")));
if (Debug()>1) printf ("Set Centrality %i %f %li\n", hbtEvent->UncorrectedNumberOfPrimaries(),
- 10.0*cent->GetCentralityPercentile("FMD"), lrint(10.0*cent->GetCentralityPercentile("FMD")));
+ 10.0*cent->GetCentralityPercentile("FMD"), lrint(10.0*cent->GetCentralityPercentile("FMD")));
}
}
}
if(fReadV0)
- {
- for (Int_t i = 0; i < fEvent->GetNumberOfV0s(); i++) {
- AliESDv0* esdv0 = fEvent->GetV0(i);
- if (!esdv0) continue;
- //if(esdv0->GetNDaughters()>2) continue;
- //if(esdv0->GetNProngs()>2) continue;
- if(esdv0->Charge()!=0) continue;
- AliESDtrack *trackPos = fEvent->GetTrack(esdv0->GetPindex());
- if(!trackPos) continue;
- AliESDtrack *trackNeg = fEvent->GetTrack(esdv0->GetNindex());
- if(!trackNeg) continue;
- if(trackPos->Charge()==trackNeg->Charge()) continue;
- AliFemtoV0* trackCopyV0 = new AliFemtoV0();
- CopyESDtoFemtoV0(esdv0, trackCopyV0, fEvent);
- hbtEvent->V0Collection()->push_back(trackCopyV0);
- //cout<<"Pushback v0 to v0collection"<<endl;
- }
+ {
+ for (Int_t i = 0; i < fEvent->GetNumberOfV0s(); i++) {
+ AliESDv0* esdv0 = fEvent->GetV0(i);
+ if (!esdv0) continue;
+ //if(esdv0->GetNDaughters()>2) continue;
+ //if(esdv0->GetNProngs()>2) continue;
+ if(esdv0->Charge()!=0) continue;
+ AliESDtrack *trackPos = fEvent->GetTrack(esdv0->GetPindex());
+ if(!trackPos) continue;
+ AliESDtrack *trackNeg = fEvent->GetTrack(esdv0->GetNindex());
+ if(!trackNeg) continue;
+ if(trackPos->Charge()==trackNeg->Charge()) continue;
+ AliFemtoV0* trackCopyV0 = new AliFemtoV0();
+ CopyESDtoFemtoV0(esdv0, trackCopyV0, fEvent);
+ hbtEvent->V0Collection()->push_back(trackCopyV0);
+ //cout<<"Pushback v0 to v0collection"<<endl;
}
+ }
fCurEvent++;
AliESDtrack *trackneg = tESDevent->GetTrack(tESDv0->GetNindex()); //AliAODTrack *trackneg = (AliAODTrack*)tESDv0->GetDaughter(1);
if(trackpos && trackneg)
+ {
+ tFemtoV0->SetdcaPosToPrimVertex(TMath::Abs(trackpos->GetD(fPrimaryVtxPosition[0],fPrimaryVtxPosition[1],tESDevent->GetMagneticField())));
+ tFemtoV0->SetdcaNegToPrimVertex(TMath::Abs(trackneg->GetD(fPrimaryVtxPosition[0],fPrimaryVtxPosition[1],tESDevent->GetMagneticField())));
+ double MomPos[3];
+ trackpos->PxPyPz(MomPos);
+ tFemtoV0->SetmomPosX(MomPos[0]);
+ tFemtoV0->SetmomPosY(MomPos[1]);
+ tFemtoV0->SetmomPosZ(MomPos[2]);
+ AliFemtoThreeVector mompos(MomPos[0],MomPos[1],MomPos[2]);
+ tFemtoV0->SetmomPos(mompos);
+
+ double MomNeg[3];
+ trackneg->PxPyPz(MomNeg);
+ tFemtoV0->SetmomNegX(MomNeg[0]);
+ tFemtoV0->SetmomNegY(MomNeg[1]);
+ tFemtoV0->SetmomNegZ(MomNeg[2]);
+ AliFemtoThreeVector momneg(MomNeg[0],MomNeg[1],MomNeg[2]);
+ tFemtoV0->SetmomNeg(momneg);
+
+ tFemtoV0->SetptPos(trackpos->Pt());
+ tFemtoV0->SetptotPos(trackpos->P());
+ tFemtoV0->SetptNeg(trackneg->Pt());
+ tFemtoV0->SetptotNeg(trackneg->P());
+
+ tFemtoV0->SetidNeg(trackneg->GetID());
+ //cout<<"tESDv0->GetNegID(): "<<tESDv0->GetNegID()<<endl;
+ //cout<<"tFemtoV0->IdNeg(): "<<tFemtoV0->IdNeg()<<endl;
+ tFemtoV0->SetidPos(trackpos->GetID());
+
+ tFemtoV0->SetEtaPos(trackpos->Eta());
+ tFemtoV0->SetEtaNeg(trackneg->Eta());
+
+ tFemtoV0->SetEtaPos(trackpos->Eta()); //tESDv0->PseudoRapPos()
+ tFemtoV0->SetEtaNeg(trackneg->Eta()); //tESDv0->PseudoRapNeg()
+ tFemtoV0->SetTPCNclsPos(trackpos->GetTPCNcls());
+ tFemtoV0->SetTPCNclsNeg(trackneg->GetTPCNcls());
+ tFemtoV0->SetTPCclustersPos(trackpos->GetTPCClusterMap());
+ tFemtoV0->SetTPCclustersNeg(trackneg->GetTPCClusterMap());
+ tFemtoV0->SetTPCsharingPos(trackpos->GetTPCSharedMap());
+ tFemtoV0->SetTPCsharingNeg(trackneg->GetTPCSharedMap());
+ tFemtoV0->SetNdofPos(trackpos->GetTPCchi2()/trackpos->GetTPCNcls());
+ tFemtoV0->SetNdofNeg(trackneg->GetTPCchi2()/trackneg->GetTPCNcls());
+ tFemtoV0->SetStatusPos(trackpos->GetStatus());
+ tFemtoV0->SetStatusNeg(trackneg->GetStatus());
+
+ float bfield = 5*fMagFieldSign;
+ float globalPositionsAtRadiiPos[9][3];
+ GetGlobalPositionAtGlobalRadiiThroughTPC(trackpos,bfield,globalPositionsAtRadiiPos);
+ double tpcEntrancePos[3]={globalPositionsAtRadiiPos[0][0],globalPositionsAtRadiiPos[0][1],globalPositionsAtRadiiPos[0][2]};
+ double tpcExitPos[3]={globalPositionsAtRadiiPos[8][0],globalPositionsAtRadiiPos[8][1],globalPositionsAtRadiiPos[8][2]};
+
+ float globalPositionsAtRadiiNeg[9][3];
+ GetGlobalPositionAtGlobalRadiiThroughTPC(trackneg,bfield,globalPositionsAtRadiiNeg);
+ double tpcEntranceNeg[3]={globalPositionsAtRadiiNeg[0][0],globalPositionsAtRadiiNeg[0][1],globalPositionsAtRadiiNeg[0][2]};
+ double tpcExitNeg[3]={globalPositionsAtRadiiNeg[8][0],globalPositionsAtRadiiNeg[8][1],globalPositionsAtRadiiNeg[8][2]};
+
+ AliFemtoThreeVector tmpVec;
+ tmpVec.SetX(tpcEntrancePos[0]); tmpVec.SetX(tpcEntrancePos[1]); tmpVec.SetX(tpcEntrancePos[2]);
+ tFemtoV0->SetNominalTpcEntrancePointPos(tmpVec);
+
+ tmpVec.SetX(tpcExitPos[0]); tmpVec.SetX(tpcExitPos[1]); tmpVec.SetX(tpcExitPos[2]);
+ tFemtoV0->SetNominalTpcExitPointPos(tmpVec);
+
+ tmpVec.SetX(tpcEntranceNeg[0]); tmpVec.SetX(tpcEntranceNeg[1]); tmpVec.SetX(tpcEntranceNeg[2]);
+ tFemtoV0->SetNominalTpcEntrancePointNeg(tmpVec);
+
+ tmpVec.SetX(tpcExitNeg[0]); tmpVec.SetX(tpcExitNeg[1]); tmpVec.SetX(tpcExitNeg[2]);
+ tFemtoV0->SetNominalTpcExitPointNeg(tmpVec);
+
+ AliFemtoThreeVector vecTpcPos[9];
+ AliFemtoThreeVector vecTpcNeg[9];
+ for(int i=0;i<9;i++)
{
- tFemtoV0->SetdcaPosToPrimVertex(TMath::Abs(trackpos->GetD(fPrimaryVtxPosition[0],fPrimaryVtxPosition[1],tESDevent->GetMagneticField())));
- tFemtoV0->SetdcaNegToPrimVertex(TMath::Abs(trackneg->GetD(fPrimaryVtxPosition[0],fPrimaryVtxPosition[1],tESDevent->GetMagneticField())));
- double MomPos[3];
- trackpos->PxPyPz(MomPos);
- tFemtoV0->SetmomPosX(MomPos[0]);
- tFemtoV0->SetmomPosY(MomPos[1]);
- tFemtoV0->SetmomPosZ(MomPos[2]);
- AliFemtoThreeVector mompos(MomPos[0],MomPos[1],MomPos[2]);
- tFemtoV0->SetmomPos(mompos);
-
- double MomNeg[3];
- trackneg->PxPyPz(MomNeg);
- tFemtoV0->SetmomNegX(MomNeg[0]);
- tFemtoV0->SetmomNegY(MomNeg[1]);
- tFemtoV0->SetmomNegZ(MomNeg[2]);
- AliFemtoThreeVector momneg(MomNeg[0],MomNeg[1],MomNeg[2]);
- tFemtoV0->SetmomNeg(momneg);
-
- tFemtoV0->SetptPos(trackpos->Pt());
- tFemtoV0->SetptotPos(trackpos->P());
- tFemtoV0->SetptNeg(trackneg->Pt());
- tFemtoV0->SetptotNeg(trackneg->P());
-
- tFemtoV0->SetidNeg(trackneg->GetID());
- //cout<<"tESDv0->GetNegID(): "<<tESDv0->GetNegID()<<endl;
- //cout<<"tFemtoV0->IdNeg(): "<<tFemtoV0->IdNeg()<<endl;
- tFemtoV0->SetidPos(trackpos->GetID());
-
- tFemtoV0->SetEtaPos(trackpos->Eta());
- tFemtoV0->SetEtaNeg(trackneg->Eta());
-
- tFemtoV0->SetEtaPos(trackpos->Eta()); //tESDv0->PseudoRapPos()
- tFemtoV0->SetEtaNeg(trackneg->Eta()); //tESDv0->PseudoRapNeg()
- tFemtoV0->SetTPCNclsPos(trackpos->GetTPCNcls());
- tFemtoV0->SetTPCNclsNeg(trackneg->GetTPCNcls());
- tFemtoV0->SetTPCclustersPos(trackpos->GetTPCClusterMap());
- tFemtoV0->SetTPCclustersNeg(trackneg->GetTPCClusterMap());
- tFemtoV0->SetTPCsharingPos(trackpos->GetTPCSharedMap());
- tFemtoV0->SetTPCsharingNeg(trackneg->GetTPCSharedMap());
- tFemtoV0->SetNdofPos(trackpos->GetTPCchi2()/trackpos->GetTPCNcls());
- tFemtoV0->SetNdofNeg(trackneg->GetTPCchi2()/trackneg->GetTPCNcls());
- tFemtoV0->SetStatusPos(trackpos->GetStatus());
- tFemtoV0->SetStatusNeg(trackneg->GetStatus());
-
- float bfield = 5*fMagFieldSign;
- float globalPositionsAtRadiiPos[9][3];
- GetGlobalPositionAtGlobalRadiiThroughTPC(trackpos,bfield,globalPositionsAtRadiiPos);
- double tpcEntrancePos[3]={globalPositionsAtRadiiPos[0][0],globalPositionsAtRadiiPos[0][1],globalPositionsAtRadiiPos[0][2]};
- double tpcExitPos[3]={globalPositionsAtRadiiPos[8][0],globalPositionsAtRadiiPos[8][1],globalPositionsAtRadiiPos[8][2]};
-
- float globalPositionsAtRadiiNeg[9][3];
- GetGlobalPositionAtGlobalRadiiThroughTPC(trackneg,bfield,globalPositionsAtRadiiNeg);
- double tpcEntranceNeg[3]={globalPositionsAtRadiiNeg[0][0],globalPositionsAtRadiiNeg[0][1],globalPositionsAtRadiiNeg[0][2]};
- double tpcExitNeg[3]={globalPositionsAtRadiiNeg[8][0],globalPositionsAtRadiiNeg[8][1],globalPositionsAtRadiiNeg[8][2]};
-
- AliFemtoThreeVector tmpVec;
- tmpVec.SetX(tpcEntrancePos[0]); tmpVec.SetX(tpcEntrancePos[1]); tmpVec.SetX(tpcEntrancePos[2]);
- tFemtoV0->SetNominalTpcEntrancePointPos(tmpVec);
-
- tmpVec.SetX(tpcExitPos[0]); tmpVec.SetX(tpcExitPos[1]); tmpVec.SetX(tpcExitPos[2]);
- tFemtoV0->SetNominalTpcExitPointPos(tmpVec);
-
- tmpVec.SetX(tpcEntranceNeg[0]); tmpVec.SetX(tpcEntranceNeg[1]); tmpVec.SetX(tpcEntranceNeg[2]);
- tFemtoV0->SetNominalTpcEntrancePointNeg(tmpVec);
-
- tmpVec.SetX(tpcExitNeg[0]); tmpVec.SetX(tpcExitNeg[1]); tmpVec.SetX(tpcExitNeg[2]);
- tFemtoV0->SetNominalTpcExitPointNeg(tmpVec);
-
- AliFemtoThreeVector vecTpcPos[9];
- AliFemtoThreeVector vecTpcNeg[9];
- for(int i=0;i<9;i++)
- {
- vecTpcPos[i].SetX(globalPositionsAtRadiiPos[i][0]); vecTpcPos[i].SetY(globalPositionsAtRadiiPos[i][1]); vecTpcPos[i].SetZ(globalPositionsAtRadiiPos[i][2]);
- vecTpcNeg[i].SetX(globalPositionsAtRadiiNeg[i][0]); vecTpcNeg[i].SetY(globalPositionsAtRadiiNeg[i][1]); vecTpcNeg[i].SetZ(globalPositionsAtRadiiNeg[i][2]);
- }
- tFemtoV0->SetNominalTpcPointPos(vecTpcPos);
- tFemtoV0->SetNominalTpcPointNeg(vecTpcNeg);
-
- tFemtoV0->SetTPCMomentumPos(trackpos->GetTPCInnerParam()->P()); //trackpos->GetTPCmomentum();
- tFemtoV0->SetTPCMomentumNeg(trackneg->GetTPCInnerParam()->P()); //trackneg->GetTPCmomentum();
-
- tFemtoV0->SetdedxPos(trackpos->GetTPCsignal());
- tFemtoV0->SetdedxNeg(trackneg->GetTPCsignal());
+ vecTpcPos[i].SetX(globalPositionsAtRadiiPos[i][0]); vecTpcPos[i].SetY(globalPositionsAtRadiiPos[i][1]); vecTpcPos[i].SetZ(globalPositionsAtRadiiPos[i][2]);
+ vecTpcNeg[i].SetX(globalPositionsAtRadiiNeg[i][0]); vecTpcNeg[i].SetY(globalPositionsAtRadiiNeg[i][1]); vecTpcNeg[i].SetZ(globalPositionsAtRadiiNeg[i][2]);
+ }
+ tFemtoV0->SetNominalTpcPointPos(vecTpcPos);
+ tFemtoV0->SetNominalTpcPointNeg(vecTpcNeg);
+ tFemtoV0->SetTPCMomentumPos(trackpos->GetTPCInnerParam()->P()); //trackpos->GetTPCmomentum();
+ tFemtoV0->SetTPCMomentumNeg(trackneg->GetTPCInnerParam()->P()); //trackneg->GetTPCmomentum();
- if (fESDpid) {
- tFemtoV0->SetPosNSigmaTPCK(fESDpid->NumberOfSigmasTPC(trackpos,AliPID::kKaon));
- tFemtoV0->SetNegNSigmaTPCK(fESDpid->NumberOfSigmasTPC(trackneg,AliPID::kKaon));
- tFemtoV0->SetPosNSigmaTPCP(fESDpid->NumberOfSigmasTPC(trackpos,AliPID::kProton));
- tFemtoV0->SetNegNSigmaTPCP(fESDpid->NumberOfSigmasTPC(trackneg,AliPID::kProton));
- tFemtoV0->SetPosNSigmaTPCPi(fESDpid->NumberOfSigmasTPC(trackpos,AliPID::kPion));
- tFemtoV0->SetNegNSigmaTPCPi(fESDpid->NumberOfSigmasTPC(trackneg,AliPID::kPion));
- }
- else {
- tFemtoV0->SetPosNSigmaTPCK(-1000);
- tFemtoV0->SetNegNSigmaTPCK(-1000);
- tFemtoV0->SetPosNSigmaTPCP(-1000);
- tFemtoV0->SetNegNSigmaTPCP(-1000);
- tFemtoV0->SetPosNSigmaTPCPi(-1000);
- tFemtoV0->SetNegNSigmaTPCPi(-1000);
- }
+ tFemtoV0->SetdedxPos(trackpos->GetTPCsignal());
+ tFemtoV0->SetdedxNeg(trackneg->GetTPCsignal());
+
+
+ if (fESDpid) {
+ tFemtoV0->SetPosNSigmaTPCK(fESDpid->NumberOfSigmasTPC(trackpos,AliPID::kKaon));
+ tFemtoV0->SetNegNSigmaTPCK(fESDpid->NumberOfSigmasTPC(trackneg,AliPID::kKaon));
+ tFemtoV0->SetPosNSigmaTPCP(fESDpid->NumberOfSigmasTPC(trackpos,AliPID::kProton));
+ tFemtoV0->SetNegNSigmaTPCP(fESDpid->NumberOfSigmasTPC(trackneg,AliPID::kProton));
+ tFemtoV0->SetPosNSigmaTPCPi(fESDpid->NumberOfSigmasTPC(trackpos,AliPID::kPion));
+ tFemtoV0->SetNegNSigmaTPCPi(fESDpid->NumberOfSigmasTPC(trackneg,AliPID::kPion));
+ }
+ else {
+ tFemtoV0->SetPosNSigmaTPCK(-1000);
+ tFemtoV0->SetNegNSigmaTPCK(-1000);
+ tFemtoV0->SetPosNSigmaTPCP(-1000);
+ tFemtoV0->SetNegNSigmaTPCP(-1000);
+ tFemtoV0->SetPosNSigmaTPCPi(-1000);
+ tFemtoV0->SetNegNSigmaTPCPi(-1000);
+ }
- if((tFemtoV0->StatusPos()&AliESDtrack::kTOFpid)==0 || (tFemtoV0->StatusPos()&AliESDtrack::kTIME)==0 || (tFemtoV0->StatusPos()&AliESDtrack::kTOFout)==0)
- {
- if((tFemtoV0->StatusNeg()&AliESDtrack::kTOFpid)==0 || (tFemtoV0->StatusNeg()&AliESDtrack::kTIME)==0 || (tFemtoV0->StatusNeg()&AliESDtrack::kTOFout)==0)
+ if((tFemtoV0->StatusPos()&AliESDtrack::kTOFpid)==0 || (tFemtoV0->StatusPos()&AliESDtrack::kTIME)==0 || (tFemtoV0->StatusPos()&AliESDtrack::kTOFout)==0)
+ {
+ if((tFemtoV0->StatusNeg()&AliESDtrack::kTOFpid)==0 || (tFemtoV0->StatusNeg()&AliESDtrack::kTIME)==0 || (tFemtoV0->StatusNeg()&AliESDtrack::kTOFout)==0)
{
tFemtoV0->SetPosNSigmaTOFK(-1000);
tFemtoV0->SetNegNSigmaTOFK(-1000);
tFemtoV0->SetPosNSigmaTOFPi(-1000);
tFemtoV0->SetNegNSigmaTOFPi(-1000);
}
- }
- else
- {
- if (fESDpid) {
- tFemtoV0->SetPosNSigmaTOFK(fESDpid->NumberOfSigmasTOF(trackpos,AliPID::kKaon));
- tFemtoV0->SetNegNSigmaTOFK(fESDpid->NumberOfSigmasTOF(trackneg,AliPID::kKaon));
- tFemtoV0->SetPosNSigmaTOFP(fESDpid->NumberOfSigmasTOF(trackpos,AliPID::kProton));
- tFemtoV0->SetNegNSigmaTOFP(fESDpid->NumberOfSigmasTOF(trackneg,AliPID::kProton));
- tFemtoV0->SetPosNSigmaTOFPi(fESDpid->NumberOfSigmasTOF(trackpos,AliPID::kPion));
- tFemtoV0->SetNegNSigmaTOFPi(fESDpid->NumberOfSigmasTOF(trackneg,AliPID::kPion));
- }
- else {
- tFemtoV0->SetPosNSigmaTOFK(-1000);
- tFemtoV0->SetNegNSigmaTOFK(-1000);
- tFemtoV0->SetPosNSigmaTOFP(-1000);
- tFemtoV0->SetNegNSigmaTOFP(-1000);
- tFemtoV0->SetPosNSigmaTOFPi(-1000);
- tFemtoV0->SetNegNSigmaTOFPi(-1000);
- }
- }
-
- if ((TMath::Abs(trackpos->GetLabel()) < fStack->GetNtrack()) && (TMath::Abs(trackneg->GetLabel()) < fStack->GetNtrack())) {
-
- //cout<<"tESDv0->GetPdgCode(): "<<tESDv0->GetPdgCode()<<endl;
- // cout<<"Labels: "<<trackpos->GetLabel()<<" "<<trackneg->GetLabel()<<endl;
- AliFemtoModelHiddenInfo *tInfo = new AliFemtoModelHiddenInfo();
- //TParticle *tPart = fStack->Particle(tESDv0->GetLabel()); //zle
+ }
+ else
+ {
+ if (fESDpid) {
+ tFemtoV0->SetPosNSigmaTOFK(fESDpid->NumberOfSigmasTOF(trackpos,AliPID::kKaon));
+ tFemtoV0->SetNegNSigmaTOFK(fESDpid->NumberOfSigmasTOF(trackneg,AliPID::kKaon));
+ tFemtoV0->SetPosNSigmaTOFP(fESDpid->NumberOfSigmasTOF(trackpos,AliPID::kProton));
+ tFemtoV0->SetNegNSigmaTOFP(fESDpid->NumberOfSigmasTOF(trackneg,AliPID::kProton));
+ tFemtoV0->SetPosNSigmaTOFPi(fESDpid->NumberOfSigmasTOF(trackpos,AliPID::kPion));
+ tFemtoV0->SetNegNSigmaTOFPi(fESDpid->NumberOfSigmasTOF(trackneg,AliPID::kPion));
+ }
+ else {
+ tFemtoV0->SetPosNSigmaTOFK(-1000);
+ tFemtoV0->SetNegNSigmaTOFK(-1000);
+ tFemtoV0->SetPosNSigmaTOFP(-1000);
+ tFemtoV0->SetNegNSigmaTOFP(-1000);
+ tFemtoV0->SetPosNSigmaTOFPi(-1000);
+ tFemtoV0->SetNegNSigmaTOFPi(-1000);
+ }
+ }
- int labelpos = TMath::Abs(trackpos->GetLabel());
- int labelneg = TMath::Abs(trackneg->GetLabel());
- TParticle *tPartPos = fStack->Particle(labelpos);
- TParticle *tPartNeg = fStack->Particle(labelneg);
+ if ((TMath::Abs(trackpos->GetLabel()) < fStack->GetNtrack()) && (TMath::Abs(trackneg->GetLabel()) < fStack->GetNtrack())) {
+ //cout<<"tESDv0->GetPdgCode(): "<<tESDv0->GetPdgCode()<<endl;
+ // cout<<"Labels: "<<trackpos->GetLabel()<<" "<<trackneg->GetLabel()<<endl;
+ AliFemtoModelHiddenInfo *tInfo = new AliFemtoModelHiddenInfo();
+ //TParticle *tPart = fStack->Particle(tESDv0->GetLabel()); //zle
- double impactpos[2];
- double impactneg[2];
- double covimpact[3];
+ int labelpos = TMath::Abs(trackpos->GetLabel());
+ int labelneg = TMath::Abs(trackneg->GetLabel());
+ TParticle *tPartPos = fStack->Particle(labelpos);
+ TParticle *tPartNeg = fStack->Particle(labelneg);
- AliExternalTrackParam *parampos = new AliExternalTrackParam(*trackpos->GetTPCInnerParam());
- parampos->PropagateToDCA(fEvent->GetPrimaryVertexTPC(), (fEvent->GetMagneticField()), 10000, impactpos, covimpact);
- AliExternalTrackParam *paramneg = new AliExternalTrackParam(*trackneg->GetTPCInnerParam());
- paramneg->PropagateToDCA(fEvent->GetPrimaryVertexTPC(), (fEvent->GetMagneticField()), 10000, impactneg, covimpact);
+ double impactpos[2];
+ double impactneg[2];
+ double covimpact[3];
- // fillDCA
- if (TMath::Abs(impactpos[0]) > 0.001) {
- if (fStack->IsPhysicalPrimary(labelpos)){
- tFemtoV0->SetImpactDprimPos(impactpos[0]);
- }
- else if (fStack->IsSecondaryFromWeakDecay(labelpos)) {
- tFemtoV0->SetImpactDweakPos(impactpos[0]);
- //cout << "wea" << endl;
- }
- else if (fStack->IsSecondaryFromMaterial(labelpos)) {
- tFemtoV0->SetImpactDmatPos(impactpos[0]);
- //cout << "mat" << endl;
- }
- }
- if (TMath::Abs(impactneg[0]) > 0.001) {
- if (fStack->IsPhysicalPrimary(labelneg)){
- tFemtoV0->SetImpactDprimNeg(impactneg[0]);
- //cout << "prim" << endl;
- }
- else if (fStack->IsSecondaryFromWeakDecay(labelneg)) {
- tFemtoV0->SetImpactDweakNeg(impactneg[0]);
- //cout << "wea" << endl;
- }
- else if (fStack->IsSecondaryFromMaterial(labelneg)) {
- tFemtoV0->SetImpactDmatNeg(impactneg[0]);
- //cout << "mat" << endl;
- }
+ AliExternalTrackParam *parampos = new AliExternalTrackParam(*trackpos->GetTPCInnerParam());
+ parampos->PropagateToDCA(fEvent->GetPrimaryVertexTPC(), (fEvent->GetMagneticField()), 10000, impactpos, covimpact);
+ AliExternalTrackParam *paramneg = new AliExternalTrackParam(*trackneg->GetTPCInnerParam());
+ paramneg->PropagateToDCA(fEvent->GetPrimaryVertexTPC(), (fEvent->GetMagneticField()), 10000, impactneg, covimpact);
+
+
+ // fillDCA
+ if (TMath::Abs(impactpos[0]) > 0.001) {
+ if (fStack->IsPhysicalPrimary(labelpos)){
+ tFemtoV0->SetImpactDprimPos(impactpos[0]);
+ }
+ else if (fStack->IsSecondaryFromWeakDecay(labelpos)) {
+ tFemtoV0->SetImpactDweakPos(impactpos[0]);
+ //cout << "wea" << endl;
+ }
+ else if (fStack->IsSecondaryFromMaterial(labelpos)) {
+ tFemtoV0->SetImpactDmatPos(impactpos[0]);
+ //cout << "mat" << endl;
+ }
+ }
+ if (TMath::Abs(impactneg[0]) > 0.001) {
+ if (fStack->IsPhysicalPrimary(labelneg)){
+ tFemtoV0->SetImpactDprimNeg(impactneg[0]);
+ //cout << "prim" << endl;
+ }
+ else if (fStack->IsSecondaryFromWeakDecay(labelneg)) {
+ tFemtoV0->SetImpactDweakNeg(impactneg[0]);
+ //cout << "wea" << endl;
+ }
+ else if (fStack->IsSecondaryFromMaterial(labelneg)) {
+ tFemtoV0->SetImpactDmatNeg(impactneg[0]);
+ //cout << "mat" << endl;
+ }
- }
- // end fillDCA
+ }
+ // end fillDCA
- //tInfo->SetPDGPid();
- //tInfo->SetMass();
- //tInfo->SetTrueMomentum();
- //tInfo->SetEmissionPoint();
+ //tInfo->SetPDGPid();
+ //tInfo->SetMass();
+ //tInfo->SetTrueMomentum();
+ //tInfo->SetEmissionPoint();
- // Freeze-out coordinates
- double fpx=0.0, fpy=0.0, fpz=0.0, fpt=0.0;
+ // Freeze-out coordinates
+ double fpx=0.0, fpy=0.0, fpz=0.0, fpt=0.0;
- fpx = tPartPos->Vx() - fPrimaryVtxPosition[0];
- fpy = tPartPos->Vy() - fPrimaryVtxPosition[1];
- fpz = tPartPos->Vz() - fPrimaryVtxPosition[2];
- fpt = tPartPos->T();
+ fpx = tPartPos->Vx() - fPrimaryVtxPosition[0];
+ fpy = tPartPos->Vy() - fPrimaryVtxPosition[1];
+ fpz = tPartPos->Vz() - fPrimaryVtxPosition[2];
+ fpt = tPartPos->T();
- fpx *= 1e13;
- fpy *= 1e13;
- fpz *= 1e13;
- fpt *= 1e13;
+ fpx *= 1e13;
+ fpy *= 1e13;
+ fpz *= 1e13;
+ fpt *= 1e13;
- tInfo->SetPDGPidPos(tPartPos->GetPdgCode());
- tInfo->SetMassPos(tPartPos->GetMass());
- tInfo->SetTrueMomentumPos(tPartPos->Px(),tPartPos->Py(),tPartPos->Pz());
- tInfo->SetEmissionPointPos(fpx,fpy,fpz,fpt);
+ tInfo->SetPDGPidPos(tPartPos->GetPdgCode());
+ tInfo->SetMassPos(tPartPos->GetMass());
+ tInfo->SetTrueMomentumPos(tPartPos->Px(),tPartPos->Py(),tPartPos->Pz());
+ tInfo->SetEmissionPointPos(fpx,fpy,fpz,fpt);
- fpx = tPartNeg->Vx() - fPrimaryVtxPosition[0];
- fpy = tPartNeg->Vy() - fPrimaryVtxPosition[1];
- fpz = tPartNeg->Vz() - fPrimaryVtxPosition[2];
- fpt = tPartNeg->T();
+ fpx = tPartNeg->Vx() - fPrimaryVtxPosition[0];
+ fpy = tPartNeg->Vy() - fPrimaryVtxPosition[1];
+ fpz = tPartNeg->Vz() - fPrimaryVtxPosition[2];
+ fpt = tPartNeg->T();
- fpx *= 1e13;
- fpy *= 1e13;
- fpz *= 1e13;
- fpt *= 1e13;
+ fpx *= 1e13;
+ fpy *= 1e13;
+ fpz *= 1e13;
+ fpt *= 1e13;
- tInfo->SetPDGPidNeg(tPartNeg->GetPdgCode());
- tInfo->SetMassNeg(tPartNeg->GetMass());
- tInfo->SetTrueMomentumNeg(tPartNeg->Px(),tPartNeg->Py(),tPartNeg->Pz());
- tInfo->SetEmissionPointNeg(fpx,fpy,fpz,fpt);
+ tInfo->SetPDGPidNeg(tPartNeg->GetPdgCode());
+ tInfo->SetMassNeg(tPartNeg->GetMass());
+ tInfo->SetTrueMomentumNeg(tPartNeg->Px(),tPartNeg->Py(),tPartNeg->Pz());
+ tInfo->SetEmissionPointNeg(fpx,fpy,fpz,fpt);
- tFemtoV0->SetHiddenInfo(tInfo);
- }
+ tFemtoV0->SetHiddenInfo(tInfo);
}
+ }
else
- {
- tFemtoV0->SetStatusPos(999);
- tFemtoV0->SetStatusNeg(999);
- }
+ {
+ tFemtoV0->SetStatusPos(999);
+ tFemtoV0->SetStatusNeg(999);
+ }
tFemtoV0->SetOnFlyStatusV0(tESDv0->GetOnFlyStatus());
}
// Bigger loop has bad precision, we're nearly one centimeter too far, go back in small steps.
while (globalRadius>Rwanted[iR]){
- x-=.1;
- // printf("propagating to x %5.2f\n",x);
- if(!etp.PropagateTo(x,bfield))break;
- etp.GetXYZ(xyz); // GetXYZ returns global coordinates
- globalRadius = TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1]); //Idea to speed up: compare squared radii
+ x-=.1;
+ // printf("propagating to x %5.2f\n",x);
+ if(!etp.PropagateTo(x,bfield))break;
+ etp.GetXYZ(xyz); // GetXYZ returns global coordinates
+ globalRadius = TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1]); //Idea to speed up: compare squared radii
}
//printf("At Radius:%05.2f (local x %5.2f). Setting position to x %4.1f y %4.1f z %4.1f\n",globalRadius,x,xyz[0],xyz[1],xyz[2]);
globalPositionsAtRadii[iR][0]=xyz[0];