Skip to content

Commit 51221cc

Browse files
skundu692skundu692
authored
[PWGLF] Add new event-mixing process function (#12916)
1 parent 116001f commit 51221cc

File tree

1 file changed

+108
-6
lines changed

1 file changed

+108
-6
lines changed

‎PWGLF/Tasks/Strangeness/lambdaspincorrderived.cxx‎

Lines changed: 108 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -138,10 +138,10 @@ struct lambdaspincorrderived {
138138
histos.add("hSparseAntiLambdaLambda", "hSparseAntiLambdLambda", HistType::kTHnSparseF, {configThnAxisInvMass, configThnAxisInvMass, configThnAxisPol, configThnAxisR}, true);
139139
histos.add("hSparseAntiLambdaAntiLambda", "hSparseAntiLambdaAntiLambda", HistType::kTHnSparseF, {configThnAxisInvMass, configThnAxisInvMass, configThnAxisPol, configThnAxisR}, true);
140140

141-
histos.add("hSparseLambdaLambdaMixed", "hSparseLambdaLambdaMixed", HistType::kTHnSparseF, {configThnAxisInvMass, configThnAxisInvMass, configThnAxisPol, configThnAxisRapidity, configThnAxisR}, true);
142-
histos.add("hSparseLambdaAntiLambdaMixed", "hSparseLambdaAntiLambdaMixed", HistType::kTHnSparseF, {configThnAxisInvMass, configThnAxisInvMass, configThnAxisPol, configThnAxisRapidity, configThnAxisR}, true);
143-
histos.add("hSparseAntiLambdaLambdaMixed", "hSparseAntiLambdaLambdaMixed", HistType::kTHnSparseF, {configThnAxisInvMass, configThnAxisInvMass, configThnAxisPol, configThnAxisRapidity, configThnAxisR}, true);
144-
histos.add("hSparseAntiLambdaAntiLambdaMixed", "hSparseAntiLambdaAntiLambdaMixed", HistType::kTHnSparseF, {configThnAxisInvMass, configThnAxisInvMass, configThnAxisPol, configThnAxisRapidity, configThnAxisR}, true);
141+
histos.add("hSparseLambdaLambdaMixed", "hSparseLambdaLambdaMixed", HistType::kTHnSparseF, {configThnAxisInvMass, configThnAxisInvMass, configThnAxisPol, configThnAxisR}, true);
142+
histos.add("hSparseLambdaAntiLambdaMixed", "hSparseLambdaAntiLambdaMixed", HistType::kTHnSparseF, {configThnAxisInvMass, configThnAxisInvMass, configThnAxisPol, configThnAxisR}, true);
143+
histos.add("hSparseAntiLambdaLambdaMixed", "hSparseAntiLambdaLambdaMixed", HistType::kTHnSparseF, {configThnAxisInvMass, configThnAxisInvMass, configThnAxisPol, configThnAxisR}, true);
144+
histos.add("hSparseAntiLambdaAntiLambdaMixed", "hSparseAntiLambdaAntiLambdaMixed", HistType::kTHnSparseF, {configThnAxisInvMass, configThnAxisInvMass, configThnAxisPol, configThnAxisR}, true);
145145

146146
histos.add("hSparseRapLambdaLambda", "hSparseRapLambdaLambda", HistType::kTHnSparseF, {configThnAxisInvMass, configThnAxisInvMass, configThnAxisPol, configThnAxisRapidity}, true);
147147
histos.add("hSparseRapLambdaAntiLambda", "hSparseRapLambdaAntiLambda", HistType::kTHnSparseF, {configThnAxisInvMass, configThnAxisInvMass, configThnAxisPol, configThnAxisRapidity}, true);
@@ -316,24 +316,27 @@ struct lambdaspincorrderived {
316316
weight3 = mixpairweight * hweight3->GetBinContent(hweight3->FindBin(particle1.Pt(), particle1.Eta(), RecoDecay::constrainAngle(particle1.Phi(), 0.0F, harmonic)));
317317
weight4 = mixpairweight * hweight4->GetBinContent(hweight4->FindBin(particle1.Pt(), particle1.Eta(), RecoDecay::constrainAngle(particle1.Phi(), 0.0F, harmonic)));
318318
}
319-
histos.fill(HIST("hPtYMix"), particle1.Pt(), particle1.Rapidity());
319+
320320
if (tag1 == 0 && tag2 == 0) {
321+
histos.fill(HIST("hPtYMix"), particle1.Pt(), particle1.Rapidity(), weight1);
321322
histos.fill(HIST("hSparseLambdaLambdaMixed"), particle1.M(), particle2.M(), cosThetaDiff, deltaR, weight1);
322323
histos.fill(HIST("hSparseRapLambdaLambdaMixed"), particle1.M(), particle2.M(), cosThetaDiff, deltaRap, weight1);
323324
histos.fill(HIST("hSparsePhiLambdaLambdaMixed"), particle1.M(), particle2.M(), costhetaz1costhetaz2, deltaPhi, weight1);
324325
histos.fill(HIST("hLambdaMixForLL"), particle1.Pt(), particle1.Eta(), RecoDecay::constrainAngle(particle1.Phi(), 0.0F, harmonic), weight1);
325326
} else if (tag1 == 0 && tag2 == 1) {
327+
histos.fill(HIST("hPtYMix"), particle1.Pt(), particle1.Rapidity(), weight2);
326328
histos.fill(HIST("hSparseLambdaAntiLambdaMixed"), particle1.M(), particle2.M(), cosThetaDiff, deltaR, weight2);
327329
histos.fill(HIST("hSparseRapLambdaAntiLambdaMixed"), particle1.M(), particle2.M(), cosThetaDiff, deltaRap, weight2);
328330
histos.fill(HIST("hSparsePhiLambdaAntiLambdaMixed"), particle1.M(), particle2.M(), costhetaz1costhetaz2, deltaPhi, weight2);
329331
histos.fill(HIST("hLambdaMixForLAL"), particle1.Pt(), particle1.Eta(), RecoDecay::constrainAngle(particle1.Phi(), 0.0F, harmonic), weight2);
330332
} else if (tag1 == 1 && tag2 == 0) {
333+
histos.fill(HIST("hPtYMix"), particle1.Pt(), particle1.Rapidity(), weight3);
331334
histos.fill(HIST("hSparseAntiLambdaLambdaMixed"), particle1.M(), particle2.M(), cosThetaDiff, deltaR, weight3);
332335
histos.fill(HIST("hSparseRapAntiLambdaLambdaMixed"), particle1.M(), particle2.M(), cosThetaDiff, deltaRap, weight3);
333336
histos.fill(HIST("hSparsePhiAntiLambdaLambdaMixed"), particle1.M(), particle2.M(), costhetaz1costhetaz2, deltaPhi, weight3);
334-
335337
histos.fill(HIST("hLambdaMixForALL"), particle1.Pt(), particle1.Eta(), RecoDecay::constrainAngle(particle1.Phi(), 0.0F, harmonic), weight3);
336338
} else if (tag1 == 1 && tag2 == 1) {
339+
histos.fill(HIST("hPtYMix"), particle1.Pt(), particle1.Rapidity(), weight4);
337340
histos.fill(HIST("hSparseAntiLambdaAntiLambdaMixed"), particle1.M(), particle2.M(), cosThetaDiff, deltaR, weight4);
338341
histos.fill(HIST("hSparseRapAntiLambdaAntiLambdaMixed"), particle1.M(), particle2.M(), cosThetaDiff, deltaRap, weight4);
339342
histos.fill(HIST("hSparsePhiAntiLambdaAntiLambdaMixed"), particle1.M(), particle2.M(), costhetaz1costhetaz2, deltaPhi, weight4);
@@ -554,6 +557,105 @@ struct lambdaspincorrderived {
554557
} // end primary-event loop
555558
}
556559
PROCESS_SWITCH(lambdaspincorrderived, processMEV2, "Process data ME", false);
560+
561+
void processMEV3(EventCandidates const& collisions, AllTrackCandidates const& V0s)
562+
{
563+
// one pool (deque) per mixing bin; each entry holds (collision index, slice of its V0s)
564+
auto nBins = colBinning.getAllBinsCount();
565+
std::vector<std::deque<std::pair<int, AllTrackCandidates>>> eventPools(nBins);
566+
567+
for (auto& collision1 : collisions) {
568+
// select mixing bin for this event
569+
const int bin = colBinning.getBin(std::make_tuple(collision1.posz(), collision1.cent()));
570+
571+
// all V0s from the current event
572+
auto poolA = V0s.sliceBy(tracksPerCollisionV0, collision1.index());
573+
574+
// loop over same-event candidate pairs (t1,t2)
575+
for (auto& [t1, t2] : soa::combinations(o2::soa::CombinationsFullIndexPolicy(poolA, poolA))) {
576+
if (!selectionV0(t1) || !selectionV0(t2))
577+
continue;
578+
if (t2.index() <= t1.index())
579+
continue; // unique unordered pairs
580+
if (t1.protonIndex() == t2.protonIndex())
581+
continue; // no shared daughter
582+
if (t1.pionIndex() == t2.pionIndex())
583+
continue;
584+
585+
// --- First pass over previous events in this bin: count replacements and build a list of usable pools
586+
int mixes = 0;
587+
struct PoolView {
588+
AllTrackCandidates* pool;
589+
int nRepl;
590+
int collIdx;
591+
};
592+
std::vector<PoolView> usable;
593+
int totalRepl = 0;
594+
595+
for (auto it = eventPools[bin].rbegin();
596+
it != eventPools[bin].rend() && mixes < nEvtMixing; ++it, ++mixes) {
597+
const int collision2idx = it->first;
598+
auto& poolB = it->second;
599+
600+
// (defensive; shouldn't happen because we push the current event after mixing)
601+
if (collision2idx == collision1.index())
602+
continue;
603+
604+
int nRepl = 0;
605+
for (auto& t3 : poolB) {
606+
if (selectionV0(t3) && checkKinematics(t1, t3))
607+
++nRepl;
608+
}
609+
if (nRepl > 0) {
610+
usable.push_back(PoolView{&poolB, nRepl, collision2idx});
611+
totalRepl += nRepl;
612+
}
613+
}
614+
615+
if (totalRepl == 0)
616+
continue;
617+
const float w = 1.0f / static_cast<float>(totalRepl); // global normalization: sum of weights over all replacements = 1
618+
619+
// --- Second pass: fill with normalized weight w
620+
for (auto& pv : usable) {
621+
auto& poolB = *pv.pool;
622+
for (auto& t3 : poolB) {
623+
if (!(selectionV0(t3) && checkKinematics(t1, t3)))
624+
continue;
625+
626+
// build 4-vectors for the mixed pair (t3 from prior event replaces t1; t2 stays from current event)
627+
proton = ROOT::Math::PtEtaPhiMVector(t3.protonPt(), t3.protonEta(), t3.protonPhi(), o2::constants::physics::MassProton);
628+
lambda = ROOT::Math::PtEtaPhiMVector(t3.lambdaPt(), t3.lambdaEta(), t3.lambdaPhi(), t3.lambdaMass());
629+
proton2 = ROOT::Math::PtEtaPhiMVector(t2.protonPt(), t2.protonEta(), t2.protonPhi(), o2::constants::physics::MassProton);
630+
lambda2 = ROOT::Math::PtEtaPhiMVector(t2.lambdaPt(), t2.lambdaEta(), t2.lambdaPhi(), t2.lambdaMass());
631+
632+
float dPhi = std::fabs(
633+
RecoDecay::constrainAngle(lambda.Phi(), 0.0F, harmonic) -
634+
RecoDecay::constrainAngle(lambda2.Phi(), 0.0F, harmonic));
635+
histos.fill(HIST("deltaPhiMix"), dPhi, w);
636+
637+
if (t3.v0Status() == 0 && t2.v0Status() == 0) {
638+
fillHistograms(0, 0, lambda, lambda2, proton, proton2, 1, w);
639+
} else if (t3.v0Status() == 0 && t2.v0Status() == 1) {
640+
fillHistograms(0, 1, lambda, lambda2, proton, proton2, 1, w);
641+
} else if (t3.v0Status() == 1 && t2.v0Status() == 0) {
642+
fillHistograms(1, 0, lambda, lambda2, proton, proton2, 1, w);
643+
} else if (t3.v0Status() == 1 && t2.v0Status() == 1) {
644+
fillHistograms(1, 1, lambda, lambda2, proton, proton2, 1, w);
645+
}
646+
}
647+
}
648+
} // end same-event pair loop
649+
650+
// after mixing with prior events, push current event into the pool
651+
auto sliced = V0s.sliceBy(tracksPerCollisionV0, collision1.index());
652+
eventPools[bin].emplace_back(collision1.index(), std::move(sliced));
653+
if (static_cast<int>(eventPools[bin].size()) > nEvtMixing) {
654+
eventPools[bin].pop_front();
655+
}
656+
} // end primary-event loop
657+
}
658+
PROCESS_SWITCH(lambdaspincorrderived, processMEV3, "Process data ME", false);
557659
};
558660
WorkflowSpec defineDataProcessing(ConfigContext const& cfgc)
559661
{

0 commit comments

Comments
 (0)