Skip to content

Commit 23d64fa

Browse files
fuchuncuifuchuncui
authored
[PWGCF] Add buffer for mixevent in cascDiHadronCorr.cxx (#16772)
1 parent bcff833 commit 23d64fa

1 file changed

Lines changed: 260 additions & 6 deletions

File tree

‎PWGCF/TwoParticleCorrelations/Tasks/cascDiHadronCorr.cxx‎

Lines changed: 260 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -237,6 +237,35 @@ struct CascDiHadronCorr {
237237
MixedEvent = 3
238238
};
239239

240+
// buffer for mixevents
241+
struct ValidCollision {
242+
struct ValidParticle {
243+
float eta;
244+
float phi;
245+
float pt;
246+
int region;
247+
float efficiency;
248+
float efficiencyError;
249+
int type;
250+
};
251+
float pvz;
252+
float mult;
253+
std::vector<ValidParticle> trigParticles;
254+
std::vector<ValidParticle> assocParticles;
255+
void addValidParticle(float eta, float phi, float pt, int region, float efficiency, float efficiencyError, int type)
256+
{
257+
ValidParticle particle{eta, phi, pt, region, efficiency, efficiencyError, type};
258+
259+
if (type == -1) {
260+
trigParticles.push_back(particle);
261+
} else {
262+
assocParticles.push_back(particle);
263+
}
264+
}
265+
};
266+
using ValidCollisions = std::vector<std::vector<ValidCollision>>;
267+
ValidCollisions validCollisions;
268+
240269
// persistent caches
241270
std::vector<float> efficiencyAssociatedCache;
242271

@@ -332,6 +361,10 @@ struct CascDiHadronCorr {
332361
registry.add("Solo_tracks_assoc", "pT", {HistType::kTH1D, {axisPtAssoc}});
333362
}
334363
}
364+
if (cfgOutputXi || cfgOutputOmega) {
365+
if (!cfgCentTableUnavailable)
366+
registry.add("Invmass", "", {HistType::kTHnSparseF, {{axisPtTrigger, axisInvMass, axisEta, axisCentrality}}});
367+
}
335368

336369
registry.add("eventcount", "bin", {HistType::kTH1F, {{4, 0, 4, "bin"}}}); // histogram to see how many events are in the same and mixed event
337370
if (doprocessMCSame && doprocessOntheflySame) {
@@ -378,11 +411,14 @@ struct CascDiHadronCorr {
378411
{axisVertexEfficiency, "z-vtx (cm)"},
379412
};
380413
std::vector<AxisSpec> userAxis;
381-
userAxis.emplace_back(axisInvMass, "m (GeV/c^2)");
414+
if (cfgOutputXi || cfgOutputOmega)
415+
userAxis.emplace_back(axisInvMass, "m (GeV/c^2)");
382416

383417
same.setObject(new CorrelationContainer("sameEvent", "sameEvent", corrAxis, effAxis, userAxis));
384418
mixed.setObject(new CorrelationContainer("mixedEvent", "mixedEvent", corrAxis, effAxis, userAxis));
385419

420+
validCollisions.resize(registry.get<TH1>(HIST("Nch"))->GetNbinsX() * registry.get<TH1>(HIST("zVtx"))->GetNbinsX());
421+
386422
LOGF(info, "End of init");
387423
}
388424

@@ -427,12 +463,28 @@ struct CascDiHadronCorr {
427463
template <typename TTrack>
428464
bool trackSelected(TTrack track)
429465
{
430-
return ((track.tpcNClsFound() >= cfgCutTPCclu) && (track.tpcNClsCrossedRows() >= cfgCutTPCCrossedRows) && (track.itsNCls() >= cfgCutITSclu));
466+
if (std::abs(track.eta()) > cfgCutEta) {
467+
return false;
468+
}
469+
if (std::abs(track.pt()) < cfgCutPtMin || std::abs(track.pt()) > cfgCutPtMax) {
470+
return false;
471+
}
472+
if ((track.tpcNClsFound() < cfgCutTPCclu) || (track.tpcNClsCrossedRows() < cfgCutTPCCrossedRows) || (track.itsNCls() < cfgCutITSclu)) {
473+
return false;
474+
}
475+
return true;
431476
}
432477

433478
template <typename TTrackCasc>
434479
bool cascSelected(TTrackCasc casc, float posX, float posY, float posZ)
435480
{
481+
if (std::abs(casc.eta()) > cfgCutEta) {
482+
return false;
483+
}
484+
if (std::abs(casc.pt()) < cfgCutPtMin || std::abs(casc.pt()) > cfgCutPtMax) {
485+
return false;
486+
}
487+
436488
auto bachelor = casc.template bachelor_as<DaughterTracks>();
437489
auto posdau = casc.template posTrack_as<DaughterTracks>();
438490
auto negdau = casc.template negTrack_as<DaughterTracks>();
@@ -680,6 +732,98 @@ struct CascDiHadronCorr {
680732
return dPhiStar;
681733
}
682734

735+
template <CorrelationContainer::CFStep step, typename TTracks, typename TCollision>
736+
void fillCorrelations(TTracks tracks1, TCollision currentCollision, float posZ, int bin, float eventWeight) // function to fill the Output functions (sparse) and the delta eta and delta phi histograms (use buffer, only for mixevent)
737+
{
738+
float triggerWeight = 1.0f;
739+
float associatedWeight = 1.0f;
740+
// loop over all validCollisions in buffer
741+
for (const auto& collision : validCollisions[bin]) {
742+
double fSampleIndex = gRandom->Uniform(0, cfgSampleSize);
743+
744+
// Cache efficiency for particles (too many FindBin lookups)
745+
if (mEfficiency) {
746+
efficiencyAssociatedCache.clear();
747+
efficiencyAssociatedCache.reserve(collision.assocParticles.size());
748+
for (const auto& track2 : currentCollision.assocParticles) {
749+
float weff = 1.;
750+
getEfficiencyCorrection(weff, track2.eta, track2.pt, posZ);
751+
efficiencyAssociatedCache.push_back(weff);
752+
}
753+
}
754+
755+
// loop over all tracks
756+
for (const auto& track1 : tracks1) {
757+
if (!trackSelected(track1))
758+
continue;
759+
if (!getEfficiencyCorrection(triggerWeight, track1.eta(), track1.pt(), posZ))
760+
continue;
761+
762+
int index = 0;
763+
for (const auto& track2 : collision.assocParticles) {
764+
if (mEfficiency) {
765+
associatedWeight = efficiencyAssociatedCache[index];
766+
index++;
767+
}
768+
if (cfgUsePtOrder && cfgUsePtOrderInMixEvent && track1.pt() <= track2.pt)
769+
continue; // For pt-differential correlations in mixed events, skip if the trigger pt is less than the associate pt
770+
771+
float deltaPhi = RecoDecay::constrainAngle(track1.phi() - track2.phi, -PIHalf);
772+
float deltaEta = track1.eta() - track2.eta;
773+
774+
mixed->getPairHist()->Fill(step, fSampleIndex, posZ, track1.pt(), track2.pt, deltaPhi, deltaEta, eventWeight * triggerWeight * associatedWeight);
775+
registry.fill(HIST("deltaEta_deltaPhi_mixed"), deltaPhi, deltaEta, eventWeight * triggerWeight * associatedWeight);
776+
}
777+
}
778+
}
779+
}
780+
781+
template <CorrelationContainer::CFStep step, typename TTracks, typename TCollision>
782+
void fillCorrelationsCasc(TTracks tracks1, TCollision currentCollision, float posX, float posY, float posZ, int bin, float eventWeight) // function to fill the Output functions (sparse) and the delta eta and delta phi histograms (use buffer, only for mixevent)
783+
{
784+
float triggerWeight = 1.0f;
785+
float associatedWeight = 1.0f;
786+
// loop over all validCollisions in buffer
787+
for (const auto& collision : validCollisions[bin]) {
788+
double fSampleIndex = gRandom->Uniform(0, cfgSampleSize);
789+
790+
// Cache efficiency for particles (too many FindBin lookups)
791+
if (mEfficiency) {
792+
efficiencyAssociatedCache.clear();
793+
efficiencyAssociatedCache.reserve(collision.assocParticles.size());
794+
for (const auto& track2 : currentCollision.assocParticles) {
795+
float weff = 1.;
796+
getEfficiencyCorrection(weff, track2.eta, track2.pt, posZ);
797+
efficiencyAssociatedCache.push_back(weff);
798+
}
799+
}
800+
801+
// loop over all tracks
802+
for (const auto& track1 : tracks1) {
803+
if (!cascSelected(track1, posX, posY, posZ))
804+
continue;
805+
if (!getEfficiencyCorrection(triggerWeight, track1.eta(), track1.pt(), posZ))
806+
continue;
807+
808+
int index = 0;
809+
for (const auto& track2 : collision.assocParticles) {
810+
if (mEfficiency) {
811+
associatedWeight = efficiencyAssociatedCache[index];
812+
index++;
813+
}
814+
if (cfgUsePtOrder && cfgUsePtOrderInMixEvent && track1.pt() <= track2.pt)
815+
continue; // For pt-differential correlations in mixed events, skip if the trigger pt is less than the associate pt
816+
817+
float deltaPhi = RecoDecay::constrainAngle(track1.phi() - track2.phi, -PIHalf);
818+
float deltaEta = track1.eta() - track2.eta;
819+
820+
mixed->getPairHist()->Fill(step, fSampleIndex, posZ, track1.pt(), track2.pt, deltaPhi, deltaEta, eventWeight * triggerWeight * associatedWeight);
821+
registry.fill(HIST("deltaEta_deltaPhi_mixed"), deltaPhi, deltaEta, eventWeight * triggerWeight * associatedWeight);
822+
}
823+
}
824+
}
825+
}
826+
683827
template <CorrelationContainer::CFStep step, typename TTracks, typename TTracksAssoc>
684828
void fillCorrelations(TTracks tracks1, TTracksAssoc tracks2, float posZ, int system, int magneticField, float cent, float eventWeight) // function to fill the Output functions (sparse) and the delta eta and delta phi histograms
685829
{
@@ -698,7 +842,7 @@ struct CascDiHadronCorr {
698842
registry.fill(HIST("Centrality_used"), cent);
699843
registry.fill(HIST("Nch_used"), tracks1.size());
700844

701-
int fSampleIndex = gRandom->Uniform(0, cfgSampleSize);
845+
double fSampleIndex = gRandom->Uniform(0, cfgSampleSize);
702846

703847
float triggerWeight = 1.0f;
704848
float associatedWeight = 1.0f;
@@ -785,7 +929,7 @@ struct CascDiHadronCorr {
785929
registry.fill(HIST("Nch_used"), tracks1.size());
786930
}
787931

788-
int fSampleIndex = gRandom->Uniform(0, cfgSampleSize);
932+
double fSampleIndex = gRandom->Uniform(0, cfgSampleSize);
789933

790934
float triggerWeight = 1.0f;
791935
float associatedWeight = 1.0f;
@@ -798,6 +942,10 @@ struct CascDiHadronCorr {
798942
continue;
799943
if (system == SameEvent) {
800944
registry.fill(HIST("Trig_hist"), fSampleIndex, posZ, track1.pt(), eventWeight * triggerWeight);
945+
if (!cfgCentTableUnavailable && cfgOutputXi)
946+
registry.fill(HIST("Invmass"), track1.pt(), track1.mXi(), track1.eta(), cent, eventWeight * triggerWeight);
947+
if (!cfgCentTableUnavailable && cfgOutputOmega)
948+
registry.fill(HIST("Invmass"), track1.pt(), track1.mOmega(), track1.eta(), cent, eventWeight * triggerWeight);
801949
}
802950

803951
auto bachelor = track1.template bachelor_as<DaughterTracks>();
@@ -892,7 +1040,7 @@ struct CascDiHadronCorr {
8921040
registry.fill(HIST("Centrality_used"), cent);
8931041
registry.fill(HIST("Nch_used"), tracks1.size());
8941042

895-
int fSampleIndex = gRandom->Uniform(0, cfgSampleSize);
1043+
double fSampleIndex = gRandom->Uniform(0, cfgSampleSize);
8961044

8971045
float triggerWeight = 1.0f;
8981046
float associatedWeight = 1.0f;
@@ -965,7 +1113,7 @@ struct CascDiHadronCorr {
9651113
template <CorrelationContainer::CFStep step, typename TTracks, typename TTracksAssoc>
9661114
void fillMCCorrelations(TTracks tracks1, TTracksAssoc tracks2, float posZ, int system, float eventWeight) // function to fill the Output functions (sparse) and the delta eta and delta phi histograms
9671115
{
968-
int fSampleIndex = gRandom->Uniform(0, cfgSampleSize);
1116+
double fSampleIndex = gRandom->Uniform(0, cfgSampleSize);
9691117

9701118
float triggerWeight = 1.0f;
9711119
float associatedWeight = 1.0f;
@@ -1312,6 +1460,112 @@ struct CascDiHadronCorr {
13121460
}
13131461
PROCESS_SWITCH(CascDiHadronCorr, processMixedCasc, "Process mixed events", true);
13141462

1463+
// the process for filling the mixed events by buffer
1464+
void processMixedBuffer(FilteredCollisions::iterator const& collision, FilteredTracks const& tracks)
1465+
{
1466+
ValidCollision currentCollision;
1467+
int nBinsMult = 0;
1468+
int binMult = 0;
1469+
int nBinsVtxZ = 0;
1470+
int binVtxZ = 0;
1471+
currentCollision.pvz = collision.posZ();
1472+
currentCollision.mult = tracks.size();
1473+
1474+
nBinsMult = registry.get<TH1>(HIST("Nch"))->GetNbinsX();
1475+
binMult = registry.get<TH1>(HIST("Nch"))->GetXaxis()->FindBin(tracks.size()) - 1;
1476+
nBinsVtxZ = registry.get<TH1>(HIST("zVtx"))->GetNbinsX();
1477+
binVtxZ = registry.get<TH1>(HIST("zVtx"))->GetXaxis()->FindBin(collision.posZ()) - 1;
1478+
if (binVtxZ < 0 || binVtxZ > nBinsVtxZ - 1 || binMult < 0 || binMult > nBinsMult - 1)
1479+
return;
1480+
int bin = binMult * nBinsVtxZ + binVtxZ;
1481+
1482+
if (!collision.sel8())
1483+
return;
1484+
1485+
if (cfgSelCollByNch && tracks.size() < cfgCutMultMin)
1486+
return;
1487+
1488+
float cent = -1;
1489+
if (!cfgCentTableUnavailable) {
1490+
cent = getCentrality(collision);
1491+
}
1492+
if (cfgUseAdditionalEventCut && !eventSelected(collision, tracks.size(), cent, false))
1493+
return;
1494+
1495+
if (!cfgSelCollByNch && !cfgCentTableUnavailable && (cent < cfgCutCentMin || cent >= cfgCutCentMax))
1496+
return;
1497+
1498+
float weightCent = 1.0f;
1499+
if (!cfgCentTableUnavailable)
1500+
getCentralityWeight(weightCent, cent);
1501+
1502+
for (const auto& track : tracks) {
1503+
if (!trackSelected(track))
1504+
continue;
1505+
currentCollision.addValidParticle(track.eta(), track.phi(), track.pt(), 0, 1, 1, 1);
1506+
}
1507+
1508+
fillCorrelations<CorrelationContainer::kCFStepReconstructed>(tracks, currentCollision, collision.posZ(), bin, weightCent);
1509+
if (validCollisions[bin].size() >= static_cast<size_t>(cfgMixEventNumMin)) {
1510+
validCollisions[bin].erase(validCollisions[bin].begin());
1511+
}
1512+
validCollisions[bin].push_back(currentCollision);
1513+
}
1514+
1515+
PROCESS_SWITCH(CascDiHadronCorr, processMixedBuffer, "Process mixed events", false);
1516+
1517+
void processMixedBufferCasc(FilteredCollisions::iterator const& collision, FilteredTracks const& tracks, aod::CascDatas const& cascades, DaughterTracks const&)
1518+
{
1519+
ValidCollision currentCollision;
1520+
int nBinsMult = 0;
1521+
int binMult = 0;
1522+
int nBinsVtxZ = 0;
1523+
int binVtxZ = 0;
1524+
currentCollision.pvz = collision.posZ();
1525+
currentCollision.mult = tracks.size();
1526+
1527+
nBinsMult = registry.get<TH1>(HIST("Nch"))->GetNbinsX();
1528+
binMult = registry.get<TH1>(HIST("Nch"))->GetXaxis()->FindBin(tracks.size()) - 1;
1529+
nBinsVtxZ = registry.get<TH1>(HIST("zVtx"))->GetNbinsX();
1530+
binVtxZ = registry.get<TH1>(HIST("zVtx"))->GetXaxis()->FindBin(collision.posZ()) - 1;
1531+
if (binVtxZ < 0 || binVtxZ > nBinsVtxZ - 1 || binMult < 0 || binMult > nBinsMult - 1)
1532+
return;
1533+
int bin = binMult * nBinsVtxZ + binVtxZ;
1534+
1535+
if (!collision.sel8())
1536+
return;
1537+
1538+
if (cfgSelCollByNch && tracks.size() < cfgCutMultMin)
1539+
return;
1540+
1541+
float cent = -1;
1542+
if (!cfgCentTableUnavailable) {
1543+
cent = getCentrality(collision);
1544+
}
1545+
if (cfgUseAdditionalEventCut && !eventSelected(collision, tracks.size(), cent, false))
1546+
return;
1547+
1548+
if (!cfgSelCollByNch && !cfgCentTableUnavailable && (cent < cfgCutCentMin || cent >= cfgCutCentMax))
1549+
return;
1550+
1551+
float weightCent = 1.0f;
1552+
if (!cfgCentTableUnavailable)
1553+
getCentralityWeight(weightCent, cent);
1554+
1555+
for (const auto& track : tracks) {
1556+
if (!trackSelected(track))
1557+
continue;
1558+
currentCollision.addValidParticle(track.eta(), track.phi(), track.pt(), 0, 1, 1, 1);
1559+
}
1560+
1561+
fillCorrelationsCasc<CorrelationContainer::kCFStepReconstructed>(cascades, currentCollision, collision.posX(), collision.posY(), collision.posZ(), bin, weightCent);
1562+
if (validCollisions[bin].size() >= static_cast<size_t>(cfgMixEventNumMin)) {
1563+
validCollisions[bin].erase(validCollisions[bin].begin());
1564+
}
1565+
validCollisions[bin].push_back(currentCollision);
1566+
}
1567+
PROCESS_SWITCH(CascDiHadronCorr, processMixedBufferCasc, "Process mixed events", true);
1568+
13151569
int getSpecies(int pdgCode)
13161570
{
13171571
switch (std::abs(pdgCode)) {

0 commit comments

Comments
 (0)