Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
34 changes: 17 additions & 17 deletions PWGCF/FemtoDream/Core/femtoDreamCollisionSelection.h
Original file line number Diff line number Diff line change
Expand Up @@ -223,10 +223,10 @@ class FemtoDreamCollisionSelection
mHistogramQn->add("Event/centVsqnVsSpher", "; cent; qn; Sphericity", kTH3F, {{10, 0, 100}, {100, 0, 1000}, {100, 0, 1}});
mHistogramQn->add("Event/qnBin", "; qnBin; entries", kTH1F, {{20, 0, 20}});
mHistogramQn->add("Event/psiEP", "; #Psi_{EP} (deg); entries", kTH1F, {{100, 0, 180}});
mHistogramQn->add("Event/epReso_FT0CTPC", "; cent; qnBin; reso_ft0c_tpc", kTH2F, {{10, 0, 100},{10,0,10}});
mHistogramQn->add("Event/epReso_FT0ATPC", "; cent; qnBin; reso_ft0a_tpc", kTH2F, {{10, 0, 100},{10,0,10}});
mHistogramQn->add("Event/epReso_FT0CFT0A", "; cent; qnBin; reso_ft0c_ft0a", kTH2F, {{10, 0, 100},{10,0,10}});
mHistogramQn->add("Event/epReso_count", "; cent; qnBin; count", kTH2F, {{10, 0, 100},{10,0,10}});
mHistogramQn->add("Event/epReso_FT0CTPC", "; cent; qnBin; reso_ft0c_tpc", kTH2F, {{10, 0, 100}, {10, 0, 10}});
mHistogramQn->add("Event/epReso_FT0ATPC", "; cent; qnBin; reso_ft0a_tpc", kTH2F, {{10, 0, 100}, {10, 0, 10}});
mHistogramQn->add("Event/epReso_FT0CFT0A", "; cent; qnBin; reso_ft0c_ft0a", kTH2F, {{10, 0, 100}, {10, 0, 10}});
mHistogramQn->add("Event/epReso_count", "; cent; qnBin; count", kTH2F, {{10, 0, 100}, {10, 0, 10}});

return;
}
Expand Down Expand Up @@ -337,13 +337,13 @@ class FemtoDreamCollisionSelection
float computeqnVec(T const& col, int qvecMod = 0)
{
double qn = -999.f;
if (qvecMod == 0){
if (qvecMod == 0) {
qn = std::sqrt(col.qvecFT0CReVec()[0] * col.qvecFT0CReVec()[0] + col.qvecFT0CImVec()[0] * col.qvecFT0CImVec()[0]) * std::sqrt(col.sumAmplFT0C());
} else if (qvecMod == 1){
} else if (qvecMod == 1) {
qn = std::sqrt(col.qvecFT0AReVec()[0] * col.qvecFT0AReVec()[0] + col.qvecFT0AImVec()[0] * col.qvecFT0AImVec()[0]) * std::sqrt(col.sumAmplFT0A());
} else {
LOGP(error, "no selected detector of Qvec for ESE ");
return qn;
LOGP(error, "no selected detector of Qvec for ESE ");
return qn;
}
return qn;
}
Expand All @@ -357,20 +357,20 @@ class FemtoDreamCollisionSelection
float computeEP(T const& col, int nmode, int qvecMod)
{
double EP = -999.f;
if (qvecMod == 0){
if (qvecMod == 0) {
EP = ((1. / nmode) * (TMath::ATan2(col.qvecFT0CImVec()[0], col.qvecFT0CReVec()[0])));
} else if (qvecMod == 1) {
EP = ((1. / nmode) * (TMath::ATan2(col.qvecFT0AImVec()[0], col.qvecFT0AReVec()[0])));
EP = ((1. / nmode) * (TMath::ATan2(col.qvecFT0AImVec()[0], col.qvecFT0AReVec()[0])));
} else if (qvecMod == 2) {
EP = ((1. / nmode) * (TMath::ATan2(col.qvecTPCallImVec()[0], col.qvecTPCallReVec()[0])));
EP = ((1. / nmode) * (TMath::ATan2(col.qvecTPCallImVec()[0], col.qvecTPCallReVec()[0])));
} else {
LOGP(error, "no selected detector of Qvec for EP");
return EP;
LOGP(error, "no selected detector of Qvec for EP");
return EP;
}
if (EP < 0){

if (EP < 0) {
EP += TMath::Pi();
} // atan2 return in rad -pi/2-pi/2, then make it 0-pi
} // atan2 return in rad -pi/2-pi/2, then make it 0-pi
return EP;
}

Expand All @@ -383,7 +383,7 @@ class FemtoDreamCollisionSelection
{
const float psi_ft0c = ((1. / nmode) * (TMath::ATan2(col.qvecFT0CImVec()[0], col.qvecFT0CReVec()[0])));
const float psi_ft0a = ((1. / nmode) * (TMath::ATan2(col.qvecFT0AImVec()[0], col.qvecFT0AReVec()[0])));
const float psi_tpc = ((1. / nmode) * (TMath::ATan2(col.qvecTPCallImVec()[0], col.qvecTPCallReVec()[0])));
const float psi_tpc = ((1. / nmode) * (TMath::ATan2(col.qvecTPCallImVec()[0], col.qvecTPCallReVec()[0])));

mHistogramQn->fill(HIST("Event/epReso_FT0CTPC"), centrality, mQnBin + 0.f, std::cos((psi_ft0c - psi_tpc) * nmode));
mHistogramQn->fill(HIST("Event/epReso_FT0ATPC"), centrality, mQnBin + 0.f, std::cos((psi_ft0a - psi_tpc) * nmode));
Expand Down
35 changes: 17 additions & 18 deletions PWGCF/FemtoDream/Core/femtoDreamContainer.h
Original file line number Diff line number Diff line change
Expand Up @@ -225,24 +225,24 @@ class FemtoDreamContainer

template <typename T>
void init_3Dqn_MC(std::string folderName, std::string femtoDKout, std::string femtoDKside, std::string femtoDKlong,
T& femtoDKoutAxis, T& femtoDKsideAxis, T& femtoDKlongAxis, bool smearingByOrigin=false)
T& femtoDKoutAxis, T& femtoDKsideAxis, T& femtoDKlongAxis, bool smearingByOrigin = false)
{
mHistogramRegistry->add((folderName + "/hNoMCtruthPairsCounter").c_str(), "; Counter; Entries", kTH1I, {{1, 0, 1}});
mHistogramRegistry->add((folderName + "/hFakePairsCounter").c_str(), "; Counter; Entries", kTH1I, {{1, 0, 1}});
if (smearingByOrigin) {
const int nOriginBins = o2::aod::femtodreamMCparticle::ParticleOriginMCTruth::kNOriginMCTruthTypes;
mHistogramRegistry->add((folderName + "/relPair3dresolution").c_str(), (";" + femtoDKout + "mctruth;" + femtoDKout + "_reco;" + femtoDKside + "mctruth;" + femtoDKside + "_reco;" +femtoDKlong + "mctruth;" +femtoDKlong + "_reco;" + "MC origin particle 1; MC origin particle 2;").c_str(),
kTHnSparseF, {femtoDKoutAxis, femtoDKoutAxis, femtoDKsideAxis, femtoDKsideAxis, femtoDKlongAxis, femtoDKlongAxis,
{nOriginBins, 0, nOriginBins}, {nOriginBins, 0, nOriginBins}});
mHistogramRegistry->add((folderName + "/relPair3dresolution").c_str(), (";" + femtoDKout + "mctruth;" + femtoDKout + "_reco;" + femtoDKside + "mctruth;" + femtoDKside + "_reco;" + femtoDKlong + "mctruth;" + femtoDKlong + "_reco;" + "MC origin particle 1; MC origin particle 2;").c_str(),
kTHnSparseF, {femtoDKoutAxis, femtoDKoutAxis, femtoDKsideAxis, femtoDKsideAxis, femtoDKlongAxis, femtoDKlongAxis, {nOriginBins, 0, nOriginBins}, {nOriginBins, 0, nOriginBins}});
} else {
mHistogramRegistry->add((folderName + "/relPair3dresolution").c_str(), (";" + femtoDKout + "mctruth;" + femtoDKside + "mctruth;" +femtoDKlong + "mctruth;" + femtoDKout + "_reco;" + femtoDKside + "_reco;" +femtoDKlong + "_reco;").c_str(),
kTHnSparseF, {femtoDKoutAxis, femtoDKoutAxis, femtoDKsideAxis, femtoDKsideAxis, femtoDKlongAxis, femtoDKlongAxis}); }
mHistogramRegistry->add((folderName + "/relPair3dresolution").c_str(), (";" + femtoDKout + "mctruth;" + femtoDKside + "mctruth;" + femtoDKlong + "mctruth;" + femtoDKout + "_reco;" + femtoDKside + "_reco;" + femtoDKlong + "_reco;").c_str(),
kTHnSparseF, {femtoDKoutAxis, femtoDKoutAxis, femtoDKsideAxis, femtoDKsideAxis, femtoDKlongAxis, femtoDKlongAxis});
}
}

template <typename T>
void init_3Dqn(HistogramRegistry* registry,
T& DKoutBins, T& DKsideBins, T& DKlongBins, T& mTBins4D, T& multPercentileBins4D,
bool isMC, ConfigurableAxis qnBins = {"qnBins", {10, 0, 10}, "qn binning"}, ConfigurableAxis pairPhiBins = {"phiBins", {10, 0 - 0.05, TMath::Pi() + 0.05}, "pair phi binning"}, bool smearingByOrigin=false)
bool isMC, ConfigurableAxis qnBins = {"qnBins", {10, 0, 10}, "qn binning"}, ConfigurableAxis pairPhiBins = {"phiBins", {10, 0 - 0.05, TMath::Pi() + 0.05}, "pair phi binning"}, bool smearingByOrigin = false)
{
mHistogramRegistry = registry;

Expand All @@ -268,7 +268,7 @@ class FemtoDreamContainer
init_base_3Dqn(folderName, femtoObsDKout, femtoObsDKside, femtoObsDKlong,
DKoutAxis, DKsideAxis, DKlongAxis, mTAxis4D, multPercentileAxis4D, qnAxis, pairPhiAxis);
init_3Dqn_MC(folderName, femtoObsDKout, femtoObsDKside, femtoObsDKlong,
DKoutAxis, DKsideAxis, DKlongAxis, smearingByOrigin);
DKoutAxis, DKsideAxis, DKlongAxis, smearingByOrigin);
}
}

Expand Down Expand Up @@ -511,16 +511,15 @@ class FemtoDreamContainer
} else {
mHistogramRegistry->fill(HIST(mFolderSuffix[mEventType]) + HIST(o2::aod::femtodreamMCparticle::MCTypeName[o2::aod::femtodreamMCparticle::MCType::kTruth]) + HIST("_3Dqn") + HIST("/relPair3dresolution"), k3dMC[1], k3d[1], k3dMC[2], k3d[2], k3dMC[3], k3d[3]);
}

}

template <bool isMC, typename T1, typename T2>
void setPair_3Dqn(T1 const& part1, T2 const& part2, const float multPercentile, bool IsSameSpecies, const float myQnBin, const float eventPlane, bool smearingByOrigin = false)
{

std::vector<double> k3d = FemtoDreamMath::newpairfunc(part1, mMassOne, part2, mMassTwo, IsSameSpecies);
if (k3d.size()<4){
LOG(error)<<"newpairfunc returned size="<< k3d.size();
if (k3d.size() < 4) {
LOG(error) << "newpairfunc returned size=" << k3d.size();
return;
}
float DKout = k3d[1];
Expand All @@ -538,8 +537,8 @@ class FemtoDreamContainer
if (part1.has_fdMCParticle() && part2.has_fdMCParticle()) {

std::vector<double> k3dMC = FemtoDreamMath::newpairfuncMC(part1.fdMCParticle(), mMassOne, part2.fdMCParticle(), mMassTwo, IsSameSpecies);
if (k3dMC.size()<4){
LOG(error)<<"newpairfunc returned size="<< k3d.size();
if (k3dMC.size() < 4) {
LOG(error) << "newpairfunc returned size=" << k3d.size();
return;
}
const float mTMC = FemtoDreamMath::getmT(part1.fdMCParticle(), mMassOne, part2.fdMCParticle(), mMassTwo);
Expand All @@ -564,8 +563,8 @@ class FemtoDreamContainer
{

std::vector<double> k3d = FemtoDreamMath::newpairfunc(part1, mMassOne, part2, mMassTwo, IsSameSpecies);
if (k3d.size()<4){
LOG(error)<<"newpairfunc returned size="<< k3d.size();
if (k3d.size() < 4) {
LOG(error) << "newpairfunc returned size=" << k3d.size();
return;
}
float DKout = k3d[1];
Expand All @@ -583,10 +582,10 @@ class FemtoDreamContainer
if (part1.has_fdMCParticle() && part2.has_fdMCParticle()) {

std::vector<double> k3dMC = FemtoDreamMath::newpairfuncMC(part1.fdMCParticle(), mMassOne, part2.fdMCParticle(), mMassTwo, IsSameSpecies);
if (k3dMC.size()<4){
LOG(error)<<"newpairfunc returned size="<< k3d.size();
if (k3dMC.size() < 4) {
LOG(error) << "newpairfunc returned size=" << k3d.size();
return;
}
}
const float mTMC = FemtoDreamMath::getmT(part1.fdMCParticle(), mMassOne, part2.fdMCParticle(), mMassTwo);
const float pairPhiEPMC = FemtoDreamMath::getPairPhiEP(part1.fdMCParticle(), mMassOne, part2.fdMCParticle(), mMassTwo, EP1, EP2);

Expand Down
20 changes: 9 additions & 11 deletions PWGCF/FemtoDream/TableProducer/femtoDreamProducerTask.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -314,9 +314,7 @@ struct femtoDreamProducerTask {
if (doprocessData == false && doprocessData_noCentrality == false && doprocessData_CentPbPb == false && doprocessData_CentPbPb_EP == false && doprocessMC == false && doprocessMC_noCentrality == false && doprocessMC_CentPbPb == false && doprocessMC_CentPbPb_EP == false) {
LOGF(fatal, "Neither processData nor processMC enabled. Please choose one.");
}
if ((doprocessData == true && doprocessMC == true) || (doprocessData == true && doprocessMC_noCentrality == true) || (doprocessMC == true && doprocessMC_noCentrality == true) || (doprocessData_noCentrality == true && doprocessData == true) || (doprocessData_noCentrality == true && doprocessMC == true) || (doprocessData_noCentrality == true && doprocessMC_noCentrality == true) || (doprocessData_CentPbPb == true && doprocessData == true) || (doprocessData_CentPbPb == true && doprocessData_noCentrality == true) || (doprocessData_CentPbPb == true && doprocessMC == true) || (doprocessData_CentPbPb == true && doprocessMC_noCentrality == true) || (doprocessData_CentPbPb == true && doprocessMC_CentPbPb == true)
|| (doprocessData_CentPbPb_EP == true && doprocessData == true) || (doprocessData_CentPbPb_EP == true && doprocessData_noCentrality == true) || (doprocessData_CentPbPb_EP == true && doprocessMC == true) || (doprocessData_CentPbPb_EP == true && doprocessMC_noCentrality == true) || (doprocessData_CentPbPb_EP == true && doprocessMC_CentPbPb == true) || (doprocessData_CentPbPb_EP == true && doprocessData_CentPbPb == true) || (doprocessData_CentPbPb_EP == true && doprocessMC_CentPbPb_EP == true)
|| (doprocessMC_CentPbPb_EP == true && doprocessData == true) || (doprocessMC_CentPbPb_EP == true && doprocessData_noCentrality == true) || (doprocessMC_CentPbPb_EP == true && doprocessMC == true) || (doprocessMC_CentPbPb_EP == true && doprocessMC_noCentrality == true) || (doprocessMC_CentPbPb_EP == true && doprocessMC_CentPbPb == true) || (doprocessMC_CentPbPb_EP == true && doprocessData_CentPbPb == true)) {
if ((doprocessData == true && doprocessMC == true) || (doprocessData == true && doprocessMC_noCentrality == true) || (doprocessMC == true && doprocessMC_noCentrality == true) || (doprocessData_noCentrality == true && doprocessData == true) || (doprocessData_noCentrality == true && doprocessMC == true) || (doprocessData_noCentrality == true && doprocessMC_noCentrality == true) || (doprocessData_CentPbPb == true && doprocessData == true) || (doprocessData_CentPbPb == true && doprocessData_noCentrality == true) || (doprocessData_CentPbPb == true && doprocessMC == true) || (doprocessData_CentPbPb == true && doprocessMC_noCentrality == true) || (doprocessData_CentPbPb == true && doprocessMC_CentPbPb == true) || (doprocessData_CentPbPb_EP == true && doprocessData == true) || (doprocessData_CentPbPb_EP == true && doprocessData_noCentrality == true) || (doprocessData_CentPbPb_EP == true && doprocessMC == true) || (doprocessData_CentPbPb_EP == true && doprocessMC_noCentrality == true) || (doprocessData_CentPbPb_EP == true && doprocessMC_CentPbPb == true) || (doprocessData_CentPbPb_EP == true && doprocessData_CentPbPb == true) || (doprocessData_CentPbPb_EP == true && doprocessMC_CentPbPb_EP == true) || (doprocessMC_CentPbPb_EP == true && doprocessData == true) || (doprocessMC_CentPbPb_EP == true && doprocessData_noCentrality == true) || (doprocessMC_CentPbPb_EP == true && doprocessMC == true) || (doprocessMC_CentPbPb_EP == true && doprocessMC_noCentrality == true) || (doprocessMC_CentPbPb_EP == true && doprocessMC_CentPbPb == true) || (doprocessMC_CentPbPb_EP == true && doprocessData_CentPbPb == true)) {
LOGF(fatal,
"Cannot enable more than one process switch at the same time. "
"Please choose one.");
Expand Down Expand Up @@ -1138,7 +1136,7 @@ struct femtoDreamProducerTask {
float myqn = 0.f;
float myEP = 0.f;

if (!isMC || epCal.ConfMCQvec){
if (!isMC || epCal.ConfMCQvec) {
myqn = colCuts.computeqnVec(col, epCal.ConfQvecDetector);
myEP = TMath::RadToDeg() * colCuts.computeEP(col, EPHarmonic, epCal.ConfEPDetector); // psi from rad(0-pi) to deg, in table psi would be in deg,from 0-180
}
Expand Down Expand Up @@ -1288,16 +1286,16 @@ struct femtoDreamProducerTask {
PROCESS_SWITCH(femtoDreamProducerTask, processMC_CentPbPb, "Provide MC data with centrality information for PbPb collisions", false);

void processMC_CentPbPb_EP(aod::FemtoFullCollisionMC_CentPbPb_qvec const& col,
aod::BCsWithTimestamps const&,
soa::Join<aod::FemtoFullTracks, aod::McTrackLabels> const& tracks,
aod::FemtoFullMCgenCollisions const&,
aod::McParticles const&,
soa::Join<o2::aod::V0Datas, aod::McV0Labels> const& fullV0s, /// \todo with FilteredFullV0s
soa::Join<o2::aod::CascDatas, aod::McCascLabels> const& fullCascades)
aod::BCsWithTimestamps const&,
soa::Join<aod::FemtoFullTracks, aod::McTrackLabels> const& tracks,
aod::FemtoFullMCgenCollisions const&,
aod::McParticles const&,
soa::Join<o2::aod::V0Datas, aod::McV0Labels> const& fullV0s, /// \todo with FilteredFullV0s
soa::Join<o2::aod::CascDatas, aod::McCascLabels> const& fullCascades)
{
// get magnetic field for run
initCCDB_Mag_Trig(col.bc_as<aod::BCsWithTimestamps>());
// fill the tables
// fill the tables
fillCollisionsAndTracksAndV0AndCascade<true, false, true, true, true>(col, tracks, tracks, fullV0s, fullCascades);
}
PROCESS_SWITCH(femtoDreamProducerTask, processMC_CentPbPb_EP, "Provide MC data with centrality and q-vector table for PbPb collisions", false);
Expand Down
14 changes: 6 additions & 8 deletions PWGCF/FemtoDream/Tasks/femtoDreamPairTaskTrackTrack.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -713,7 +713,7 @@ struct femtoDreamPairTaskTrackTrack {
float myEP = -999.f;
int myqnBin = -999;

if (!isMC || EPCal.mcQvec){
if (!isMC || EPCal.mcQvec) {
myEP = TMath::DegToRad() * col.eventPlane();
if (EPCal.doQnSeparation || EPCal.do3DFemto) {
myqnBin = epCalculator.myqnBin(col.multV0M(), EPCal.centMax, EPCal.fillFlowQA, EPCal.qnBinSeparator, col.qnVal(), EPCal.numQnBins, EPCal.centBinWidth);
Expand Down Expand Up @@ -748,7 +748,7 @@ struct femtoDreamPairTaskTrackTrack {
}
if (EPCal.do3DFemto) {
sameEventQnCont.setPair_3Dqn<isMC>(p1, p2, col.multV0M(), Option.SameSpecies.value, myqnBin + 0.f, myEP, Option.SmearingByOrigin);
}
}
} else {
if (EPCal.do1DFemto) {
sameEventQnCont.setPair_EP<isMC>(p2, p1, col.multV0M(), EPCal.doQnSeparation, EPCal.doQnSeparation ? myqnBin + 0.f : myEP);
Expand Down Expand Up @@ -801,14 +801,14 @@ struct femtoDreamPairTaskTrackTrack {
/// process function for to call doSameEventEP with MC Data
/// \param col subscribe to the collision table (Data)
/// \param parts subscribe to the femtoDreamParticleTable
void processSameEventEPMC(FilteredMCQnCollision& col,
void processSameEventEPMC(FilteredMCQnCollision& col,
o2::aod::FDMCCollisions&,
soa::Join<o2::aod::FDParticles, o2::aod::FDMCLabels>& parts,
o2::aod::FDMCParticles&)
{
if (EPCal.storeEvtTrkInfo) {
fillCollision<true>(col);
}
}
auto SliceTrk1 = PartitionMCTrk1->sliceByCached(aod::femtodreamparticle::fdCollisionId, col.globalIndex(), cache);
auto SliceTrk2 = PartitionMCTrk2->sliceByCached(aod::femtodreamparticle::fdCollisionId, col.globalIndex(), cache);
if (SliceTrk1.size() == 0 && SliceTrk2.size() == 0) {
Expand All @@ -833,16 +833,14 @@ struct femtoDreamPairTaskTrackTrack {
auto myEP_event1 = -999.f;
auto myEP_event2 = -999.f;

if (!isMC || EPCal.mcQvec){
if (!isMC || EPCal.mcQvec) {
myEP_event1 = TMath::DegToRad() * collision1.eventPlane();
myEP_event2 = TMath::DegToRad() * collision2.eventPlane();
}
else{
} else {
myEP_event1 = 0.f;
myEP_event2 = 0.f;
}


for (auto& [p1, p2] : combinations(CombinationsFullIndexPolicy(SliceTrk1, SliceTrk2))) {
if (Option.CPROn.value) {
if (pairCloseRejectionME.isClosePair(p1, p2, parts, collision1.magField())) {
Expand Down