Skip to content

Add antilambda inv mass, noo CPR for opposite charge particles, fix f… #5736

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
zchochul merged 1 commit into from
Apr 20, 2024
Merged

Add antilambda inv mass, noo CPR for opposite charge particles, fix f… #5736

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
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
75 changes: 56 additions & 19 deletions PWGCF/FemtoDream/Core/femtoDreamDetaDphiStar.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -41,11 +41,13 @@ class FemtoDreamDetaDphiStar
/// Destructor
virtual ~FemtoDreamDetaDphiStar() = default;
/// Initialization of the histograms and setting required values
void init(HistogramRegistry* registry, HistogramRegistry* registryQA, float ldeltaPhiMax, float ldeltaEtaMax, bool lplotForEveryRadii, int meORse = 0, bool oldversion = true)
void init(HistogramRegistry* registry, HistogramRegistry* registryQA, float ldeltaPhiMax, float ldeltaEtaMax, bool lplotForEveryRadii, int meORse = 0, bool oldversion = true, float Q3Limit = 8., bool isMELambda = false)
{
deltaPhiMax = ldeltaPhiMax;
deltaEtaMax = ldeltaEtaMax;
plotForEveryRadii = lplotForEveryRadii;
upperQ3LimitForPlotting = Q3Limit;
isMixedEventLambda = isMELambda;
runOldVersion = oldversion;
mHistogramRegistry = registry;
mHistogramRegistryQA = registryQA;
Expand Down Expand Up @@ -75,7 +77,7 @@ class FemtoDreamDetaDphiStar
}
/// Check if pair is close or not
template <typename Part, typename Parts>
bool isClosePair(Part const& part1, Part const& part2, Parts const& particles, float lmagfield)
bool isClosePair(Part const& part1, Part const& part2, Parts const& particles, float lmagfield, float Q3 = 999.)
{
magfield = lmagfield;

Expand All @@ -87,12 +89,25 @@ class FemtoDreamDetaDphiStar
return false;
}
auto deta = part1.eta() - part2.eta();
auto dphiAvg = AveragePhiStar(part1, part2, 0);
histdetadpi[0][0]->Fill(deta, dphiAvg);
if (pow(dphiAvg, 2) / pow(deltaPhiMax, 2) + pow(deta, 2) / pow(deltaEtaMax, 2) < 1.) {
return true;
bool sameCharge = false;
auto dphiAvg = AveragePhiStar(part1, part2, 0, &sameCharge);
if (Q3 == 999) {
histdetadpi[0][0]->Fill(deta, dphiAvg);
} else if (Q3 < upperQ3LimitForPlotting) {
histdetadpi[0][0]->Fill(deta, dphiAvg);
}
if (sameCharge) {
if (pow(dphiAvg, 2) / pow(deltaPhiMax, 2) + pow(deta, 2) / pow(deltaEtaMax, 2) < 1.) {
return true;
} else {
if (Q3 == 999) {
histdetadpi[0][1]->Fill(deta, dphiAvg);
} else if (Q3 < upperQ3LimitForPlotting) {
histdetadpi[0][1]->Fill(deta, dphiAvg);
}
return false;
}
} else {
histdetadpi[0][1]->Fill(deta, dphiAvg);
return false;
}

Expand All @@ -106,15 +121,31 @@ class FemtoDreamDetaDphiStar

bool pass = false;
for (int i = 0; i < 2; i++) {
auto indexOfDaughter = part2.index() - 2 + i;
int indexOfDaughter;
if (isMixedEventLambda) {
indexOfDaughter = part2.globalIndex() - 2 + i;
} else {
indexOfDaughter = part2.index() - 2 + i;
}
auto daughter = particles.begin() + indexOfDaughter;
auto deta = part1.eta() - daughter.eta();
auto dphiAvg = AveragePhiStar(part1, *daughter, i);
histdetadpi[i][0]->Fill(deta, dphiAvg);
if (pow(dphiAvg, 2) / pow(deltaPhiMax, 2) + pow(deta, 2) / pow(deltaEtaMax, 2) < 1.) {
pass = true;
} else {
histdetadpi[i][1]->Fill(deta, dphiAvg);
bool sameCharge = false;
auto dphiAvg = AveragePhiStar(part1, *daughter, i, &sameCharge);
if (Q3 == 999) {
histdetadpi[i][0]->Fill(deta, dphiAvg);
} else if (Q3 < upperQ3LimitForPlotting) {
histdetadpi[i][0]->Fill(deta, dphiAvg);
}
if (sameCharge) {
if (pow(dphiAvg, 2) / pow(deltaPhiMax, 2) + pow(deta, 2) / pow(deltaEtaMax, 2) < 1.) {
pass = true;
} else {
if (Q3 == 999) {
histdetadpi[i][1]->Fill(deta, dphiAvg);
} else if (Q3 < upperQ3LimitForPlotting) {
histdetadpi[i][1]->Fill(deta, dphiAvg);
}
}
}
}
return pass;
Expand Down Expand Up @@ -153,6 +184,8 @@ class FemtoDreamDetaDphiStar
float deltaEtaMax;
float magfield;
bool plotForEveryRadii = false;
bool isMixedEventLambda = false;
float upperQ3LimitForPlotting = 8.;
// a possible bug was found, but this must be tested on hyperloop with larger statistics
// possiboility to run old code is turned on so a proper comparison of both code versions can be done
bool runOldVersion = true;
Expand All @@ -163,12 +196,12 @@ class FemtoDreamDetaDphiStar
/// Calculate phi at all required radii stored in tmpRadiiTPC
/// Magnetic field to be provided in Tesla
template <typename T>
void PhiAtRadiiTPC(const T& part, std::vector<float>& tmpVec)
int PhiAtRadiiTPC(const T& part, std::vector<float>& tmpVec)
{

float phi0 = part.phi();
// Start: Get the charge from cutcontainer using masks
float charge = 0.;
int charge = 0.;
if ((part.cut() & kSignMinusMask) == kValue0 && (part.cut() & kSignPlusMask) == kValue0) {
charge = 0;
} else if ((part.cut() & kSignPlusMask) == kSignPlusMask) {
Expand All @@ -194,16 +227,20 @@ class FemtoDreamDetaDphiStar
}
}
}
return charge;
}

/// Calculate average phi
template <typename T1, typename T2>
float AveragePhiStar(const T1& part1, const T2& part2, int iHist)
float AveragePhiStar(const T1& part1, const T2& part2, int iHist, bool* sameCharge)
{
std::vector<float> tmpVec1;
std::vector<float> tmpVec2;
PhiAtRadiiTPC(part1, tmpVec1);
PhiAtRadiiTPC(part2, tmpVec2);
auto charge1 = PhiAtRadiiTPC(part1, tmpVec1);
auto charge2 = PhiAtRadiiTPC(part2, tmpVec2);
if (charge1 == charge2) {
*sameCharge = true;
}
int num = tmpVec1.size();
int meaningfulEntries = num;
float dPhiAvg = 0;
Expand Down
11 changes: 10 additions & 1 deletion PWGCF/FemtoDream/Tasks/femtoDreamCollisionMasker.cxx
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -254,6 +254,10 @@ struct femoDreamCollisionMasker {
FilterPtMax.at(CollisionMasks::kPartOne).push_back(option.defaultValue.get<float>());
} else if (option.name.compare(std::string("ConfMinpT")) == 0) {
FilterPtMin.at(CollisionMasks::kPartOne).push_back(option.defaultValue.get<float>());
} else if (option.name.compare(std::string("ConfMaxDCAxy")) == 0) {
FilterTempFitVarMax.at(CollisionMasks::kPartOne).push_back(option.defaultValue.get<float>());
} else if (option.name.compare(std::string("ConfMinDCAxy")) == 0) {
FilterTempFitVarMin.at(CollisionMasks::kPartOne).push_back(option.defaultValue.get<float>());
}
}
} else if (device.name.find("femto-dream-triplet-task-track-track-v0") != std::string::npos) {
Expand All @@ -272,6 +276,10 @@ struct femoDreamCollisionMasker {
FilterPtMax.at(CollisionMasks::kPartOne).push_back(option.defaultValue.get<float>());
} else if (option.name.compare(std::string("ConfMinpT")) == 0) {
FilterPtMin.at(CollisionMasks::kPartOne).push_back(option.defaultValue.get<float>());
} else if (option.name.compare(std::string("ConfMaxDCAxy")) == 0) {
FilterTempFitVarMax.at(CollisionMasks::kPartOne).push_back(option.defaultValue.get<float>());
} else if (option.name.compare(std::string("ConfMinDCAxy")) == 0) {
FilterTempFitVarMin.at(CollisionMasks::kPartOne).push_back(option.defaultValue.get<float>());
} else if (option.name.compare(std::string("ConfCutV0")) == 0) {
V0CutBits.at(CollisionMasks::kPartThree).push_back(option.defaultValue.get<femtodreamparticle::cutContainerType>());
} else if (option.name.compare(std::string("Conf_ChildPos_CutV0")) == 0) {
Expand Down Expand Up @@ -362,7 +370,8 @@ struct femoDreamCollisionMasker {
continue;
}
// check filter cuts
if (track.pt() < FilterPtMin.at(P).at(index) || track.pt() > FilterPtMax.at(P).at(index)) {
if (track.pt() < FilterPtMin.at(P).at(index) || track.pt() > FilterPtMax.at(P).at(index) ||
track.tempFitVar() > FilterTempFitVarMax.at(P).at(index) || track.tempFitVar() < FilterTempFitVarMin.at(P).at(index)) {
// if they are not passed, skip the particle
continue;
}
Expand Down
39 changes: 23 additions & 16 deletions PWGCF/FemtoDream/Tasks/femtoDreamTripletTaskTrackTrackTrack.cxx
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -11,9 +11,6 @@

/// \file femtoDreamTripletTaskTrackTrackTrack.cxx
/// \brief Tasks that reads the track tables and creates track triplets; only three identical particles can be used
/// \author Andi Mathis, TU München, andreas.mathis@ph.tum.de
/// \author Georgios Mantzaridis, TU München, georgios.mantzaridis@tum.de
/// \author Anton Riedel, TU München, anton.riedel@tum.de
/// \author Laura Serksnyte, TU München, laura.serksnyte@tum.de

#include <vector>
Expand Down Expand Up @@ -58,6 +55,8 @@ struct femtoDreamTripletTaskTrackTrackTrack {
Configurable<float> ConfTracksInMixedEvent{"ConfTracksInMixedEvent", 1, "Number of tracks of interest, contained in the mixed event sample: 1 - only events with at least one track of interest are used in mixing; ...; 3 - only events with at least three track of interest are used in mixing. Max value is 3"};
Configurable<float> ConfMaxpT{"ConfMaxpT", 4.05f, "Maximum transverse momentum of the particles"};
Configurable<float> ConfMinpT{"ConfMinpT", 0.3f, "Minimum transverse momentum of the particles"};
Configurable<float> ConfMaxDCAxy{"ConfMaxDCAxy", -0.1f, "Maximum DCAxy of the particles"};
Configurable<float> ConfMinDCAxy{"ConfMinDCAxy", 0.1f, "Minimum DCAxy of the particles"};
Configurable<float> ConfPIDthrMom{"ConfPIDthrMom", 1.f, "Momentum threshold from which TPC and TOF are required for PID"};
Configurable<o2::aod::femtodreamparticle::cutContainerType> ConfTPCPIDBit{"ConfTPCPIDBit", 16, "PID TPC bit from cutCulator "};
Configurable<o2::aod::femtodreamparticle::cutContainerType> ConfTPCTOFPIDBit{"ConfTPCTOFPIDBit", 8, "PID TPCTOF bit from cutCulator"};
Expand All @@ -73,12 +72,16 @@ struct femtoDreamTripletTaskTrackTrackTrack {
ifnode(aod::femtodreamparticle::pt * (nexp(aod::femtodreamparticle::eta) + nexp(-1.f * aod::femtodreamparticle::eta)) / 2.f <= ConfPIDthrMom, ncheckbit(aod::femtodreamparticle::pidcut, ConfTPCPIDBit), ncheckbit(aod::femtodreamparticle::pidcut, ConfTPCTOFPIDBit)) &&
(ncheckbit(aod::femtodreamparticle::cut, ConfCutPart)) &&
(aod::femtodreamparticle::pt < ConfMaxpT) &&
(aod::femtodreamparticle::pt > ConfMinpT);
(aod::femtodreamparticle::pt > ConfMinpT) &&
(aod::femtodreamparticle::tempFitVar < ConfMaxDCAxy) &&
(aod::femtodreamparticle::tempFitVar > ConfMinDCAxy);
Partition<soa::Join<aod::FDParticles, aod::FDMCLabels>> SelectedPartsMC = (aod::femtodreamparticle::partType == uint8_t(aod::femtodreamparticle::ParticleType::kTrack)) &&
ifnode(aod::femtodreamparticle::pt * (nexp(aod::femtodreamparticle::eta) + nexp(-1.f * aod::femtodreamparticle::eta)) / 2.f <= ConfPIDthrMom, ncheckbit(aod::femtodreamparticle::pidcut, ConfTPCPIDBit), ncheckbit(aod::femtodreamparticle::pidcut, ConfTPCTOFPIDBit)) &&
(ncheckbit(aod::femtodreamparticle::cut, ConfCutPart)) &&
(aod::femtodreamparticle::pt < ConfMaxpT) &&
(aod::femtodreamparticle::pt > ConfMinpT);
(aod::femtodreamparticle::pt > ConfMinpT) &&
(aod::femtodreamparticle::tempFitVar < ConfMaxDCAxy) &&
(aod::femtodreamparticle::tempFitVar > ConfMinDCAxy);

/// Histogramming of Selected Particles
FemtoDreamParticleHisto<aod::femtodreamparticle::ParticleType::kTrack, 1> trackHistoSelectedParts;
Expand All @@ -105,6 +108,7 @@ struct femtoDreamTripletTaskTrackTrackTrack {
Configurable<bool> ConfCPRPlotPerRadii{"ConfCPRPlotPerRadii", false, "Plot CPR per radii"};
Configurable<float> ConfCPRdeltaPhiMax{"ConfCPRdeltaPhiMax", 0.01, "Max. Delta Phi for Close Pair Rejection"};
Configurable<float> ConfCPRdeltaEtaMax{"ConfCPRdeltaEtaMax", 0.01, "Max. Delta Eta for Close Pair Rejection"};
Configurable<float> ConfMaxQ3IncludedInCPRPlots{"ConfMaxQ3IncludedInCPRPlots", 8., "Maximum Q3, for which the pair CPR is included in plots"};
ConfigurableAxis ConfDummy{"ConfDummy", {1, 0, 1}, "Dummy axis"};

FemtoDreamContainerThreeBody<femtoDreamContainerThreeBody::EventType::same, femtoDreamContainerThreeBody::Observable::Q3> sameEventCont;
Expand Down Expand Up @@ -138,8 +142,8 @@ struct femtoDreamTripletTaskTrackTrackTrack {
mixedEventCont.setPDGCodes(ConfPDGCodePart, ConfPDGCodePart, ConfPDGCodePart);
pairCleaner.init(&qaRegistry); // SERKSNYTE : later check if init should be updated to have 3 separate histos
if (ConfIsCPR.value) {
pairCloseRejectionSE.init(&resultRegistry, &qaRegistry, ConfCPRdeltaPhiMax.value, ConfCPRdeltaEtaMax.value, ConfCPRPlotPerRadii.value, 1, ConfUseOLD_possiblyWrong_CPR);
pairCloseRejectionME.init(&resultRegistry, &qaRegistry, ConfCPRdeltaPhiMax.value, ConfCPRdeltaEtaMax.value, ConfCPRPlotPerRadii.value, 2, ConfUseOLD_possiblyWrong_CPR);
pairCloseRejectionSE.init(&resultRegistry, &qaRegistry, ConfCPRdeltaPhiMax.value, ConfCPRdeltaEtaMax.value, ConfCPRPlotPerRadii.value, 1, ConfUseOLD_possiblyWrong_CPR, ConfMaxQ3IncludedInCPRPlots);
pairCloseRejectionME.init(&resultRegistry, &qaRegistry, ConfCPRdeltaPhiMax.value, ConfCPRdeltaEtaMax.value, ConfCPRPlotPerRadii.value, 2, ConfUseOLD_possiblyWrong_CPR, ConfMaxQ3IncludedInCPRPlots);
}

// get masses
Expand All @@ -158,7 +162,9 @@ struct femtoDreamTripletTaskTrackTrackTrack {
containsNameValuePair(device.options, "ConfTPCTOFPIDBit", ConfTPCTOFPIDBit.value) &&
containsNameValuePair(device.options, "ConfPIDthrMom", ConfPIDthrMom.value) &&
containsNameValuePair(device.options, "ConfMaxpT", ConfMaxpT.value) &&
containsNameValuePair(device.options, "ConfMinpT", ConfMinpT.value)) {
containsNameValuePair(device.options, "ConfMinpT", ConfMinpT.value) &&
containsNameValuePair(device.options, "ConfMaxDCAxy", ConfMaxDCAxy.value) &&
containsNameValuePair(device.options, "ConfMinDCAxy", ConfMinDCAxy.value)) {
mask.set(index);
MaskBit = static_cast<aod::femtodreamcollision::BitMaskType>(mask.to_ulong());
LOG(info) << "Device name matched: " << device.name;
Expand Down Expand Up @@ -203,14 +209,16 @@ struct femtoDreamTripletTaskTrackTrackTrack {

/// Now build the combinations
for (auto& [p1, p2, p3] : combinations(CombinationsStrictlyUpperIndexPolicy(groupSelectedParts, groupSelectedParts, groupSelectedParts))) {
auto Q3 = FemtoDreamMath::getQ3(p1, mMassOne, p2, mMassTwo, p3, mMassThree);

if (ConfIsCPR.value) {
if (pairCloseRejectionSE.isClosePair(p1, p2, parts, magFieldTesla)) {
if (pairCloseRejectionSE.isClosePair(p1, p2, parts, magFieldTesla, Q3)) {
continue;
}
if (pairCloseRejectionSE.isClosePair(p2, p3, parts, magFieldTesla)) {
if (pairCloseRejectionSE.isClosePair(p2, p3, parts, magFieldTesla, Q3)) {
continue;
}
if (pairCloseRejectionSE.isClosePair(p1, p3, parts, magFieldTesla)) {
if (pairCloseRejectionSE.isClosePair(p1, p3, parts, magFieldTesla, Q3)) {
continue;
}
}
Expand All @@ -226,7 +234,6 @@ struct femtoDreamTripletTaskTrackTrackTrack {
continue;
}
// fill pT of all three particles as a function of Q3 for lambda calculations
auto Q3 = FemtoDreamMath::getQ3(p1, mMassOne, p2, mMassTwo, p3, mMassThree);
ThreeBodyQARegistry.fill(HIST("TripletTaskQA/particle_pT_in_Triplet_SE"), p1.pt(), p2.pt(), p3.pt(), Q3);
sameEventCont.setTriplet<isMC>(p1, p2, p3, multCol, Q3);
}
Expand Down Expand Up @@ -322,20 +329,20 @@ struct femtoDreamTripletTaskTrackTrackTrack {
void doMixedEvent(PartitionType groupPartsOne, PartitionType groupPartsTwo, PartitionType groupPartsThree, PartType parts, float magFieldTesla, int multCol)
{
for (auto& [p1, p2, p3] : combinations(CombinationsFullIndexPolicy(groupPartsOne, groupPartsTwo, groupPartsThree))) {
auto Q3 = FemtoDreamMath::getQ3(p1, mMassOne, p2, mMassTwo, p3, mMassThree);
if (ConfIsCPR.value) {
if (pairCloseRejectionME.isClosePair(p1, p2, parts, magFieldTesla)) {
if (pairCloseRejectionME.isClosePair(p1, p2, parts, magFieldTesla, Q3)) {
continue;
}
if (pairCloseRejectionME.isClosePair(p2, p3, parts, magFieldTesla)) {
if (pairCloseRejectionME.isClosePair(p2, p3, parts, magFieldTesla, Q3)) {
continue;
}

if (pairCloseRejectionME.isClosePair(p1, p3, parts, magFieldTesla)) {
if (pairCloseRejectionME.isClosePair(p1, p3, parts, magFieldTesla, Q3)) {
continue;
}
}
// fill pT of all three particles as a function of Q3 for lambda calculations
auto Q3 = FemtoDreamMath::getQ3(p1, mMassOne, p2, mMassTwo, p3, mMassThree);
ThreeBodyQARegistry.fill(HIST("TripletTaskQA/particle_pT_in_Triplet_ME"), p1.pt(), p2.pt(), p3.pt(), Q3);
mixedEventCont.setTriplet<isMC>(p1, p2, p3, multCol, Q3);
}
Expand Down
Loading