Skip to content

Commit bdef18b

Browse files
JimunLeejimun_lee
andauthored
[PWGLF] Fixed reconstruction part of KstarInOO.cxx (#13283)
Co-authored-by: jimun_lee <jimun.lee@cern.ch>
1 parent aa4012a commit bdef18b

File tree

1 file changed

+141
-45
lines changed

1 file changed

+141
-45
lines changed

PWGLF/Tasks/Resonances/kstarInOO.cxx

Lines changed: 141 additions & 45 deletions
Original file line numberDiff line numberDiff line change
@@ -179,21 +179,25 @@ struct kstarInOO {
179179
}
180180

181181
if (cfgMcHistos) {
182+
histos.add("hPion_PID_Purity", "hPion_PID_Purity", kTH1F, {{3, -1.5, 1.5}});
183+
histos.add("hKaon_PID_Purity", "hKaon_PID_Purity", kTH1F, {{3, -1.5, 1.5}});
184+
histos.add("hSimplePion_PID_Purity", "hSimplePion_PID_Purity", kTH1F, {{3, -1.5, 1.5}});
185+
histos.add("hSimpleKaon_PID_Purity", "hSimpleKaon_PID_Purity", kTH1F, {{3, -1.5, 1.5}});
186+
182187
histos.add("nEvents_MC", "nEvents_MC", kTH1F, {{4, 0.0, 4.0}});
183188
histos.add("nEvents_MC_True", "nEvents_MC_True", kTH1F, {{4, 0.0, 4.0}});
184189

185-
histos.add("hMC_USS", "hMC_USS", kTHnSparseF, {cfgCentAxis, ptAxis, minvAxis});
186-
histos.add("hMC_LSS", "hMC_LSS", kTHnSparseF, {cfgCentAxis, ptAxis, minvAxis});
187-
histos.add("hMC_USS_Mix", "hMC_USS_Mix", kTHnSparseF, {cfgCentAxis, ptAxis, minvAxis});
188-
histos.add("hMC_LSS_Mix", "hMC_LSS_Mix", kTHnSparseF, {cfgCentAxis, ptAxis, minvAxis});
189-
histos.add("hMC_USS_True", "hMC_USS_True", kTHnSparseF, {cfgCentAxis, ptAxis, minvAxis});
190190
histos.add("hMC_kstar_True", "hMC_kstar_True", kTHnSparseF, {cfgCentAxis, ptAxis});
191191

192-
histos.add("hMC_USS_pion", "hMC_USS_pion", kTHnSparseF, {cfgCentAxis, ptAxis, minvAxis});
193-
histos.add("hMC_LSS_pion", "hMC_LSS_pion", kTHnSparseF, {cfgCentAxis, ptAxis, minvAxis});
194-
histos.add("hMC_USS_Mix_pion", "hMC_USS_Mix_pion", kTHnSparseF, {cfgCentAxis, ptAxis, minvAxis});
195-
histos.add("hMC_LSS_Mix_pion", "hMC_LSS_Mix_pion", kTHnSparseF, {cfgCentAxis, ptAxis, minvAxis});
196-
histos.add("hMC_USS_pion_True", "hMC_USS_pion_True", kTHnSparseF, {cfgCentAxis, ptAxis, minvAxis});
192+
histos.add("hMC_USS_True", "hMC_USS_True", kTHnSparseF, {cfgCentAxis, ptAxis, minvAxis});
193+
histos.add("hMC_USS_KPi", "hMC_USS_KPi", kTHnSparseF, {cfgCentAxis, ptAxis, minvAxis});
194+
histos.add("hMC_USS_PiK", "hMC_USS_PiK", kTHnSparseF, {cfgCentAxis, ptAxis, minvAxis});
195+
histos.add("hMC_LSS_KPi", "hMC_LSS_KPi", kTHnSparseF, {cfgCentAxis, ptAxis, minvAxis});
196+
histos.add("hMC_LSS_PiK", "hMC_LSS_PiK", kTHnSparseF, {cfgCentAxis, ptAxis, minvAxis});
197+
histos.add("hMC_USS_KPi_Mix", "hMC_USS_KPi_Mix", kTHnSparseF, {cfgCentAxis, ptAxis, minvAxis});
198+
histos.add("hMC_USS_PiK_Mix", "hMC_USS_PiK_Mix", kTHnSparseF, {cfgCentAxis, ptAxis, minvAxis});
199+
histos.add("hMC_USS_KPi_True", "hMC_USS_KPi_True", kTHnSparseF, {cfgCentAxis, ptAxis, minvAxis});
200+
histos.add("hMC_USS_PiK_True", "hMC_USS_PiK_True", kTHnSparseF, {cfgCentAxis, ptAxis, minvAxis});
197201
}
198202
} // end of init
199203

@@ -226,7 +230,6 @@ struct kstarInOO {
226230
histos.fill(HIST("hPosZ_BC"), event.posZ());
227231
histos.fill(HIST("hcentFT0C_BC"), event.centFT0C());
228232
}
229-
230233
if (!event.sel8())
231234
return false;
232235
if (std::abs(event.posZ()) > cfgEventVtxCut)
@@ -264,6 +267,7 @@ struct kstarInOO {
264267
histos.fill(HIST("hTPCChi2_BC"), track.tpcChi2NCl());
265268
histos.fill(HIST("QA_track_pT_BC"), track.pt());
266269
}
270+
267271
if (cfgTrackGlobalSel && !track.isGlobalTrack())
268272
return false;
269273
if (track.pt() < cfgTrackMinPt)
@@ -278,11 +282,11 @@ struct kstarInOO {
278282
return false;
279283
if (cfgTrackGlobalWoDCATrack && !track.isGlobalTrackWoDCA())
280284
return false;
281-
if (track.tpcNClsFindable() < cfgTracknFindableTPCClusters)
285+
if (cfgTracknFindableTPCClusters > 0 && track.tpcNClsFindable() < cfgTracknFindableTPCClusters)
282286
return false;
283287
if (track.tpcNClsCrossedRows() < cfgTracknTPCCrossedRows)
284288
return false;
285-
if (track.tpcCrossedRowsOverFindableCls() > cfgTracknRowsOverFindable)
289+
if (cfgTracknRowsOverFindable > 0 && track.tpcCrossedRowsOverFindableCls() > cfgTracknRowsOverFindable)
286290
return false;
287291
if (track.tpcChi2NCl() > cfgTracknTPCChi2)
288292
return false;
@@ -358,7 +362,7 @@ struct kstarInOO {
358362
} else {
359363
tofPIDPassed = true;
360364
}
361-
// TPC & TOF
365+
// TPC & TOF
362366
if (tpcPIDPassed && tofPIDPassed) {
363367
if (cfgTrackCutQA && QA) {
364368
histos.fill(HIST("QA_nSigma_pion_TPC_AC"), candidate.pt(), candidate.tpcNSigmaPi());
@@ -378,10 +382,8 @@ struct kstarInOO {
378382
auto centrality = collision1.centFT0C();
379383

380384
for (const auto& [trk1, trk2] : combinations(o2::soa::CombinationsFullIndexPolicy(tracks1, tracks2))) {
381-
if (std::fabs(trk1.signed1Pt()) <= 0.f || std::fabs(trk2.signed1Pt()) <= 0.f)
382-
continue;
383385

384-
auto [KstarPt, Minv] = minvReconstruction(trk1, trk2, QA);
386+
auto [KstarPt, Minv] = minvReconstruction(trk1, trk2, QA, false);
385387

386388
double conjugate = trk1.sign() * trk2.sign();
387389
if (cfgDataHistos) {
@@ -411,43 +413,49 @@ struct kstarInOO {
411413
auto tracks2 = pionMC->sliceByCached(aod::track::collisionId, collision2.globalIndex(), cache);
412414
auto centrality = collision1.centFT0C();
413415

416+
std::vector<int> mcMemory;
417+
std::vector<int> PIDPurityKey_Kaon;
418+
std::vector<int> PIDPurityKey_Pion;
419+
420+
double KstarPt_Kpi, Minv_Kpi;
421+
414422
for (const auto& [trk1, trk2] : combinations(o2::soa::CombinationsFullIndexPolicy(tracks1, tracks2))) {
415423
if (!trk1.has_mcParticle() || !trk2.has_mcParticle())
416424
continue;
417-
if (std::fabs(trk1.signed1Pt()) <= 0.f || std::fabs(trk2.signed1Pt()) <= 0.f)
425+
426+
// auto [KstarPt_Kpi, Minv_Kpi] = minvReconstruction(trk1, trk2, QA, false);
427+
// auto [KstarPt_piK, Minv_piK] = minvReconstruction(trk1, trk2, QA, true);
428+
429+
std::tie(KstarPt_Kpi, Minv_Kpi) = minvReconstruction(trk1, trk2, QA, false);
430+
std::tie(KstarPt_Kpi, Minv_Kpi) = minvReconstruction(trk1, trk2, QA, true);
431+
432+
if (Minv_Kpi < 0)
418433
continue;
419-
auto [KstarPt_Kpi, Minv_Kpi] = minvReconstruction(trk1, trk2, QA);
420434

421435
double conjugate = trk1.sign() * trk2.sign();
422436
if (cfgMcHistos) {
423-
if (Minv_Kpi > 0) {
424-
if (!IsMix) {
425-
if (conjugate < 0) {
426-
histos.fill(HIST("hMC_USS"), centrality, KstarPt_Kpi, Minv_Kpi);
427-
} else if (conjugate > 0) {
428-
histos.fill(HIST("hMC_LSS"), centrality, KstarPt_Kpi, Minv_Kpi);
429-
}
430-
} else {
431-
if (conjugate < 0) {
432-
histos.fill(HIST("hMC_USS_Mix"), centrality, KstarPt_Kpi, Minv_Kpi);
433-
} else if (conjugate > 0) {
434-
histos.fill(HIST("hMC_LSS_Mix"), centrality, KstarPt_Kpi, Minv_Kpi);
435-
}
437+
if (!IsMix) {
438+
if (conjugate < 0) {
439+
histos.fill(HIST("hMC_USS_KPi"), centrality, KstarPt_Kpi, Minv_Kpi);
440+
} else if (conjugate > 0) {
441+
histos.fill(HIST("hMC_LSS_KPi"), centrality, KstarPt_Kpi, Minv_Kpi);
442+
}
443+
} else {
444+
if (conjugate < 0) {
445+
histos.fill(HIST("hMC_USS_KPi_Mix"), centrality, KstarPt_Kpi, Minv_Kpi);
436446
}
437447
}
438448
}
439449
//======================
440450
// Gen MC
441451
auto particle1 = trk1.mcParticle();
442452
auto particle2 = trk2.mcParticle();
443-
if (std::fabs(particle1.pdgCode()) != 321)
444-
continue; // Not Kaon
445-
if (std::fabs(particle2.pdgCode()) != 211)
446-
continue; // Not Pion
447453

448454
if (!particle1.has_mothers() || !particle2.has_mothers()) {
449455
continue;
450456
}
457+
int mcindex1 = trk1.globalIndex();
458+
int mcindex2 = trk2.globalIndex();
451459

452460
std::vector<int> mothers1{};
453461
std::vector<int> mothers1PDG{};
@@ -471,27 +479,85 @@ struct kstarInOO {
471479
if (mothers1[0] != mothers2[0])
472480
continue; // Kaon and pion not from the same K*0
473481

482+
if (std::fabs(particle1.pdgCode()) != 211 && std::fabs(particle1.pdgCode()) != 321)
483+
continue;
484+
if (std::fabs(particle2.pdgCode()) != 211 && std::fabs(particle2.pdgCode()) != 321)
485+
continue;
486+
487+
double track1_mass, track2_mass;
488+
bool track1f{false}; // true means pion
489+
490+
if (std::fabs(particle1.pdgCode()) == 211) {
491+
track1f = true;
492+
track1_mass = massPi;
493+
} else {
494+
track1_mass = massKa;
495+
}
496+
497+
if (std::fabs(particle2.pdgCode()) == 211) {
498+
track2_mass = massPi;
499+
} else {
500+
track2_mass = massKa;
501+
}
502+
503+
if (track1_mass == track2_mass) {
504+
return;
505+
}
506+
507+
bool exists1 = std::find(mcMemory.begin(), mcMemory.end(), mcindex1) != mcMemory.end();
508+
bool exists2 = std::find(mcMemory.begin(), mcMemory.end(), mcindex2) != mcMemory.end();
509+
if (exists1 || exists2) {
510+
continue;
511+
} else {
512+
mcMemory.push_back(trk1.globalIndex());
513+
mcMemory.push_back(trk2.globalIndex());
514+
}
515+
516+
TLorentzVector lDecayDaughter1, lDecayDaughter2, lResonance;
517+
lDecayDaughter1.SetXYZM(trk1.px(), trk1.py(), trk1.pz(), track1_mass);
518+
lDecayDaughter2.SetXYZM(trk2.px(), trk2.py(), trk2.pz(), track2_mass);
519+
lResonance = lDecayDaughter1 + lDecayDaughter2;
520+
474521
if (cfgMcHistos) {
475-
histos.fill(HIST("hMC_USS_True"), centrality, KstarPt_Kpi, Minv_Kpi);
522+
histos.fill(HIST("hMC_USS_True"), centrality, lResonance.Pt(), lResonance.M());
523+
if (track1f) {
524+
histos.fill(HIST("hMC_USS_PiK_True"), centrality, lResonance.Pt(), lResonance.M());
525+
} else {
526+
histos.fill(HIST("hMC_USS_KPi_True"), centrality, lResonance.Pt(), lResonance.M());
527+
}
476528
}
477529
//======================
478530
} // for
479531
} // TrackSlicingMC
480532

481533
template <typename TracksType>
482-
std::pair<double, double> minvReconstruction(const TracksType& trk1, const TracksType& trk2, const bool QA)
534+
std::pair<double, double> minvReconstruction(const TracksType& trk1, const TracksType& trk2, const bool QA, const bool flip)
483535
{
484536
if (!trackSelection(trk1, false) || !trackSelection(trk2, false))
485537
return {-1.0, -1.0};
486-
if (!trackPIDKaon(trk1, QA) || !trackPIDPion(trk2, QA))
487-
return {-1.0, -1.0};
488-
if (trk1.globalIndex() == trk2.globalIndex())
538+
539+
if (!flip) {
540+
if (!trackPIDKaon(trk1, QA) || !trackPIDPion(trk2, QA)) {
541+
return {-1.0, -1.0};
542+
}
543+
} else {
544+
if (!trackPIDPion(trk1, false) || !trackPIDKaon(trk2, false))
545+
return {-1.0, -1.0};
546+
}
547+
548+
// if (trk1.globalIndex() == trk2.globalIndex())
549+
// return {-1.0, -1.0};
550+
if (trk1.index() >= trk2.index())
489551
return {-1.0, -1.0};
490552

491553
TLorentzVector lDecayDaughter1, lDecayDaughter2, lResonance;
492-
493-
lDecayDaughter1.SetXYZM(trk1.px(), trk1.py(), trk1.pz(), massKa);
494-
lDecayDaughter2.SetXYZM(trk2.px(), trk2.py(), trk2.pz(), massPi);
554+
if (!flip) {
555+
lDecayDaughter1.SetXYZM(trk1.px(), trk1.py(), trk1.pz(), massKa);
556+
lDecayDaughter2.SetXYZM(trk2.px(), trk2.py(), trk2.pz(), massPi);
557+
} else {
558+
lDecayDaughter1.SetXYZM(trk1.px(), trk1.py(), trk1.pz(), massPi);
559+
lDecayDaughter2.SetXYZM(trk2.px(), trk2.py(), trk2.pz(), massKa);
560+
}
495561
lResonance = lDecayDaughter1 + lDecayDaughter2;
496562

497563
if (std::abs(lResonance.Eta()) > cfgTrackMaxEta)
@@ -573,7 +639,7 @@ struct kstarInOO {
573639
//|
574640
//| MC STUFF (SE)
575641
//|
576-
//=======================================================
642+
//=========================================================
577643
int nEventsMC = 0;
578644
void processSameEventMC(EventCandidates::iterator const& collision, TrackCandidatesMC const& tracks, aod::McParticles const&)
579645
{
@@ -603,6 +669,36 @@ struct kstarInOO {
603669
if (!INELgt0)
604670
return;
605671

672+
auto tracks1 = kaonMC->sliceByCached(aod::track::collisionId, collision.globalIndex(), cache);
673+
for (const auto& kaon : tracks1) {
674+
if (!trackSelection(kaon, false))
675+
continue;
676+
if (!trackPIDKaon(kaon, false))
677+
continue;
678+
auto particle1 = kaon.mcParticle();
679+
if (std::fabs(particle1.pdgCode()) == 321)
680+
histos.fill(HIST("hSimpleKaon_PID_Purity"), 1); // histogram with two bins, -1.5, 1.5 fill 1 or -1
681+
else if (std::fabs(particle1.pdgCode()) == 211)
682+
histos.fill(HIST("hSimpleKaon_PID_Purity"), -1); // histogram with two bins, -1.5, 1.5 fill 1 or -1
683+
else
684+
histos.fill(HIST("hSimpleKaon_PID_Purity"), 0); // histogram with two bins, -1.5, 1.5 fill 1 or -1
685+
}
686+
687+
auto tracks2 = pionMC->sliceByCached(aod::track::collisionId, collision.globalIndex(), cache);
688+
for (const auto& pion : tracks2) {
689+
if (!trackSelection(pion, false))
690+
continue;
691+
if (!trackPIDPion(pion, false))
692+
continue;
693+
auto particle2 = pion.mcParticle();
694+
if (std::fabs(particle2.pdgCode()) == 211)
695+
histos.fill(HIST("hSimplePion_PID_Purity"), 1); // histogram with two bins, -1.5, 1.5 fill 1 or -1
696+
else if (std::fabs(particle2.pdgCode()) == 321)
697+
histos.fill(HIST("hSimplePion_PID_Purity"), -1); // histogram with two bins, -1.5, 1.5 fill 1 or -1
698+
else
699+
histos.fill(HIST("hSimplePion_PID_Purity"), 0); // histogram with two bins, -1.5, 1.5 fill 1 or -1
700+
}
701+
606702
if (cfgMcHistos) {
607703
histos.fill(HIST("nEvents_MC"), 1.5);
608704
}

0 commit comments

Comments
 (0)