Skip to content
54 changes: 31 additions & 23 deletions PWGHF/HFC/Tasks/taskCharmHadronsFemtoDream.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -125,8 +125,8 @@ struct HfTaskCharmHadronsFemtoDream {
FemtoDreamContainer<femtoDreamContainer::EventType::same, femtoDreamContainer::Observable::kstar> sameEventCont;
FemtoDreamContainer<femtoDreamContainer::EventType::mixed, femtoDreamContainer::Observable::kstar> mixedEventCont;
FemtoDreamPairCleaner<aod::femtodreamparticle::ParticleType::kTrack, aod::femtodreamparticle::ParticleType::kCharmHadron> pairCleaner;
FemtoDreamDetaDphiStar<aod::femtodreamparticle::ParticleType::kTrack, aod::femtodreamparticle::ParticleType::kCharmHadron> pairCloseRejection;

FemtoDreamDetaDphiStar<aod::femtodreamparticle::ParticleType::kTrack, aod::femtodreamparticle::ParticleType::kCharmHadron> pairCloseRejectionSE;
FemtoDreamDetaDphiStar<aod::femtodreamparticle::ParticleType::kTrack, aod::femtodreamparticle::ParticleType::kCharmHadron> pairCloseRejectionME;
Filter eventMultiplicity = aod::femtodreamcollision::multNtr >= eventSel.multMin && aod::femtodreamcollision::multNtr <= eventSel.multMax;
Filter eventMultiplicityPercentile = aod::femtodreamcollision::multV0M >= eventSel.multPercentileMin && aod::femtodreamcollision::multV0M <= eventSel.multPercentileMax;
Filter hfCandSelFilter = aod::fdhf::candidateSelFlag >= charmHadCandSel.value;
Expand Down Expand Up @@ -205,7 +205,8 @@ struct HfTaskCharmHadronsFemtoDream {
mixedEventCont.setPDGCodes(pdgCodeTrack1, charmHadPDGCode);
pairCleaner.init(&registry);
if (useCPR.value) {
pairCloseRejection.init(&registry, &registry, cprDeltaPhiMax.value, cprDeltaEtaMax.value, cprPlotPerRadii.value);
pairCloseRejectionSE.init(&registry, &registry, cprDeltaPhiMax.value, cprDeltaEtaMax.value, cprPlotPerRadii.value, 1);
pairCloseRejectionME.init(&registry, &registry, cprDeltaPhiMax.value, cprDeltaEtaMax.value, cprPlotPerRadii.value, 2);
}
}

Expand All @@ -221,6 +222,17 @@ struct HfTaskCharmHadronsFemtoDream {
}

for (auto const& [p1, p2] : combinations(CombinationsFullIndexPolicy(sliceTrk1, sliceCharmHad))) {

if (useCPR.value) {
if (pairCloseRejectionSE.isClosePair(p1, p2, parts, col.magField())) {
continue;
}
}

if (!pairCleaner.isCleanPair(p1, p2, parts)) {
continue;
}

// proton track charge
float chargeTrack = 0.;
if ((p1.cut() & 2) == 2) {
Expand All @@ -247,16 +259,6 @@ struct HfTaskCharmHadronsFemtoDream {
// partSign = 1 << 1;
// }

if (useCPR.value) {
if (pairCloseRejection.isClosePair(p1, p2, parts, col.magField())) {
continue;
}
}

if (!pairCleaner.isCleanPair(p1, p2, parts)) {
continue;
}

float invMass;
if (p2.candidateSelFlag() == 1) {
invMass = p2.m(std::array{o2::constants::physics::MassProton, o2::constants::physics::MassKPlus, o2::constants::physics::MassPiPlus});
Expand Down Expand Up @@ -303,14 +305,28 @@ struct HfTaskCharmHadronsFemtoDream {
template <bool isMc, typename CollisionType, typename PartType, typename PartitionType1, typename PartitionType2, typename BinningType>
void doMixedEvent(CollisionType const& cols, PartType const& parts, PartitionType1& part1, PartitionType2& part2, BinningType policy)
{
processType = 1 << 1; // for mixed event

for (auto const& [collision1, collision2] : soa::selfCombinations(policy, mixingDepth.value, -1, cols, cols)) {
// Mixed events that contain the pair of interest
processType = 2; // for mixed event

for (auto const& [collision1, collision2] : combinations(soa::CombinationsBlockFullSameIndexPolicy(policy, mixingDepth.value, -1, cols, cols))) {
// make sure that tracks in the same events are not mixed
if (collision1.globalIndex() == collision2.globalIndex()) {
continue;
}
auto sliceTrk1 = part1->sliceByCached(aod::femtodreamparticle::fdCollisionId, collision1.globalIndex(), cache);
auto sliceCharmHad = part2->sliceByCached(aod::femtodreamparticle::fdCollisionId, collision2.globalIndex(), cache);
for (auto& [p1, p2] : combinations(CombinationsFullIndexPolicy(sliceTrk1, sliceCharmHad))) {

if (useCPR.value) {
if (pairCloseRejectionME.isClosePair(p1, p2, parts, collision1.magField())) {
continue;
}
}
if (!pairCleaner.isCleanPair(p1, p2, parts)) {
continue;
}

float chargeTrack = 0.;
if ((p1.cut() & 2) == 2) {
chargeTrack = PositiveCharge;
Expand Down Expand Up @@ -361,14 +377,6 @@ struct HfTaskCharmHadronsFemtoDream {
charmHadMc,
originType);

if (useCPR.value) {
if (pairCloseRejection.isClosePair(p1, p2, parts, collision1.magField())) {
continue;
}
}
if (!pairCleaner.isCleanPair(p1, p2, parts)) {
continue;
}
// if constexpr (!isMc) mixedEventCont.setPair<isMc, true>(p1, p2, collision1.multNtr(), collision1.multV0M(), use4D, extendedPlots, smearingByOrigin);
mixedEventCont.setPair<isMc, true>(p1, p2, collision1.multNtr(), collision1.multV0M(), use4D, extendedPlots, smearingByOrigin);
}
Expand Down