Skip to content

Commit c9097c4

Browse files
dsekihatdsekihat
authored
[PWGEM/Dilepton] update HF->ll for ULS-LS in MC (AliceO2Group#14354)
1 parent 0b29a89 commit c9097c4

File tree

1 file changed

+76
-34
lines changed

1 file changed

+76
-34
lines changed

‎PWGEM/Dilepton/Utils/MCUtilities.h‎

Lines changed: 76 additions & 34 deletions
Original file line numberDiff line numberDiff line change
@@ -191,6 +191,50 @@ bool isWeakDecayFromBeautyHadron(T const& mcParticle, U const& mcParticles)
191191
}
192192
}
193193
//_______________________________________________________________________
194+
template <typename TMCParticle, typename TMCParticles>
195+
bool isDiquark(TMCParticle const& p1, TMCParticle const& p2, TMCParticles const& mcParticles)
196+
{
197+
bool isDiquark = false;
198+
199+
int motherid1 = p1.mothersIds()[0]; // first mother index
200+
while (motherid1 > -1) {
201+
if (motherid1 < mcParticles.size()) { // protect against bad mother indices. why is this needed?
202+
auto mp1 = mcParticles.iteratorAt(motherid1);
203+
if (mp1.has_mothers()) {
204+
motherid1 = mp1.mothersIds()[0];
205+
} else {
206+
motherid1 = -999;
207+
}
208+
} else {
209+
LOGF(info, "Mother label(%d) exceeds the McParticles size(%d)", motherid1, mcParticles.size());
210+
break;
211+
}
212+
213+
int motherid2 = p2.mothersIds()[0]; // first mother index
214+
while (motherid2 > -1) {
215+
if (motherid2 < mcParticles.size()) { // protect against bad mother indices. why is this needed?
216+
auto mp2 = mcParticles.iteratorAt(motherid2);
217+
if (mp2.has_mothers()) {
218+
motherid2 = mp2.mothersIds()[0];
219+
} else {
220+
motherid2 = -999;
221+
}
222+
} else {
223+
LOGF(info, "Mother label(%d) exceeds the McParticles size(%d)", motherid2, mcParticles.size());
224+
break;
225+
}
226+
227+
if (motherid1 == motherid2) { // common mother is found.
228+
auto common_mp = mcParticles.iteratorAt(motherid1);
229+
int mp_pdg = common_mp.pdgCode();
230+
isDiquark = (1100 < std::abs(mp_pdg) && std::abs(mp_pdg) < 5600) && std::to_string(mp_pdg)[std::to_string(mp_pdg).length() - 2] == '0';
231+
}
232+
233+
} // end of motherid2
234+
} // end of motherid1
235+
236+
return isDiquark;
237+
}
194238
//_______________________________________________________________________
195239
template <typename TMCParticle1, typename TMCParticle2>
196240
int FindCommonMotherFrom2ProngsWithoutPDG(TMCParticle1 const& p1, TMCParticle2 const& p2)
@@ -518,6 +562,10 @@ int IsHF(TMCParticle1 const& p1, TMCParticle2 const& p2, TMCParticles const& mcp
518562
return static_cast<int>(EM_HFeeType::kUndef); // this never happens in correlated HF->ee decays
519563
}
520564

565+
if (isDiquark(p1, p2, mcparticles)) {
566+
return static_cast<int>(EM_HFeeType::kUndef); // remove diquark
567+
}
568+
521569
// store all mother1 relation
522570
std::vector<int> mothers_id1;
523571
std::vector<int> mothers_pdg1;
@@ -560,12 +608,10 @@ int IsHF(TMCParticle1 const& p1, TMCParticle2 const& p2, TMCParticles const& mcp
560608
}
561609
}
562610

563-
// require correlation between q-qbar. (not q-q) // need statusCode
564-
565-
auto mpfh1 = mcparticles.iteratorAt(find1stHadron(p1, mcparticles));
566-
auto mpfh2 = mcparticles.iteratorAt(find1stHadron(p2, mcparticles));
567-
bool isFOat1stDecay1 = isFlavorOscillationB(mpfh1); // oscillation occured at 1st hb decay.
568-
bool isFOat1stDecay2 = isFlavorOscillationB(mpfh2); // oscillation occured at 1st hb decay.
611+
// auto mpfh1 = mcparticles.iteratorAt(find1stHadron(p1, mcparticles));
612+
// auto mpfh2 = mcparticles.iteratorAt(find1stHadron(p2, mcparticles));
613+
// bool isFOat1stDecay1 = isFlavorOscillationB(mpfh1); // oscillation occured at 1st hb decay.
614+
// bool isFOat1stDecay2 = isFlavorOscillationB(mpfh2); // oscillation occured at 1st hb decay.
569615

570616
bool is_direct_from_b1 = isWeakDecayFromBeautyHadron(p1, mcparticles);
571617
bool is_direct_from_b2 = isWeakDecayFromBeautyHadron(p2, mcparticles);
@@ -574,7 +620,7 @@ int IsHF(TMCParticle1 const& p1, TMCParticle2 const& p2, TMCParticles const& mcp
574620
bool is_c_from_b1 = isWeakDecayFromCharmHadron(p1, mcparticles) && IsFromBeauty(p1, mcparticles) > 0;
575621
bool is_c_from_b2 = isWeakDecayFromCharmHadron(p2, mcparticles) && IsFromBeauty(p2, mcparticles) > 0;
576622

577-
if (is_prompt_c1 && is_prompt_c2 && mpfh1.pdgCode() * mpfh2.pdgCode() < 0) { // charmed mesons never oscillate. Be careful with D(2007)*0 -> D0 e+ e-
623+
if (is_prompt_c1 && is_prompt_c2) { // don't check sign of first hadrons. // charmed mesons never oscillate. Be careful with D(2007)*0 -> D0 e+ e-
578624
mothers_id1.clear();
579625
mothers_pdg1.clear();
580626
mothers_id2.clear();
@@ -586,11 +632,11 @@ int IsHF(TMCParticle1 const& p1, TMCParticle2 const& p2, TMCParticles const& mcp
586632
return static_cast<int>(EM_HFeeType::kCe_Ce); // cc->ee, decay type = 0
587633
}
588634

589-
bool b2l_b2l_case0 = is_direct_from_b1 && is_direct_from_b2 && mpfh1.pdgCode() * mpfh2.pdgCode() < 0 && !isFOat1stDecay1 && !isFOat1stDecay2; // bbbar -> ll ULS
590-
bool b2l_b2l_case1 = is_direct_from_b1 && is_direct_from_b2 && mpfh1.pdgCode() * mpfh2.pdgCode() > 0 && static_cast<bool>(isFOat1stDecay1 ^ isFOat1stDecay2); // bbbar -> ll LS
591-
bool b2l_b2l_case2 = is_direct_from_b1 && is_direct_from_b2 && mpfh1.pdgCode() * mpfh2.pdgCode() < 0 && isFOat1stDecay1 && isFOat1stDecay2; // bbbar -> ll ULS
635+
// bool b2l_b2l_case0 = is_direct_from_b1 && is_direct_from_b2 && mpfh1.pdgCode() * mpfh2.pdgCode() < 0 && !isFOat1stDecay1 && !isFOat1stDecay2; // bbbar -> ll ULS
636+
// bool b2l_b2l_case1 = is_direct_from_b1 && is_direct_from_b2 && mpfh1.pdgCode() * mpfh2.pdgCode() > 0 && static_cast<bool>(isFOat1stDecay1 ^ isFOat1stDecay2); // bbbar -> ll LS
637+
// bool b2l_b2l_case2 = is_direct_from_b1 && is_direct_from_b2 && mpfh1.pdgCode() * mpfh2.pdgCode() < 0 && isFOat1stDecay1 && isFOat1stDecay2; // bbbar -> ll ULS
592638

593-
if (b2l_b2l_case0 || b2l_b2l_case1 || b2l_b2l_case2) {
639+
if (is_direct_from_b1 && is_direct_from_b2) { // analyzer should do ULS - LS to take into flavor oscillation account.
594640
mothers_id1.clear();
595641
mothers_pdg1.clear();
596642
mothers_id2.clear();
@@ -602,11 +648,11 @@ int IsHF(TMCParticle1 const& p1, TMCParticle2 const& p2, TMCParticles const& mcp
602648
return static_cast<int>(EM_HFeeType::kBe_Be); // bb->ee, decay type = 2
603649
}
604650

605-
bool b2c2l_b2c2l_case0 = is_c_from_b1 && is_c_from_b2 && mpfh1.pdgCode() * mpfh2.pdgCode() < 0 && !isFOat1stDecay1 && !isFOat1stDecay2; // bbbar -> ll ULS
606-
bool b2c2l_b2c2l_case1 = is_c_from_b1 && is_c_from_b2 && mpfh1.pdgCode() * mpfh2.pdgCode() > 0 && static_cast<bool>(isFOat1stDecay1 ^ isFOat1stDecay2); // bbbar -> ll LS
607-
bool b2c2l_b2c2l_case2 = is_c_from_b1 && is_c_from_b2 && mpfh1.pdgCode() * mpfh2.pdgCode() < 0 && isFOat1stDecay1 && isFOat1stDecay2; // bbbar -> ll ULS
651+
// bool b2c2l_b2c2l_case0 = is_c_from_b1 && is_c_from_b2 && mpfh1.pdgCode() * mpfh2.pdgCode() < 0 && !isFOat1stDecay1 && !isFOat1stDecay2; // bbbar -> ll ULS
652+
// bool b2c2l_b2c2l_case1 = is_c_from_b1 && is_c_from_b2 && mpfh1.pdgCode() * mpfh2.pdgCode() > 0 && static_cast<bool>(isFOat1stDecay1 ^ isFOat1stDecay2); // bbbar -> ll LS
653+
// bool b2c2l_b2c2l_case2 = is_c_from_b1 && is_c_from_b2 && mpfh1.pdgCode() * mpfh2.pdgCode() < 0 && isFOat1stDecay1 && isFOat1stDecay2; // bbbar -> ll ULS
608654

609-
if (b2c2l_b2c2l_case0 || b2c2l_b2c2l_case1 || b2c2l_b2c2l_case2) {
655+
if (is_c_from_b1 && is_c_from_b2) {
610656
mothers_id1.clear();
611657
mothers_pdg1.clear();
612658
mothers_id2.clear();
@@ -619,14 +665,12 @@ int IsHF(TMCParticle1 const& p1, TMCParticle2 const& p2, TMCParticles const& mcp
619665
}
620666

621667
if ((is_direct_from_b1 && is_c_from_b2) || (is_direct_from_b2 && is_c_from_b1)) {
622-
// No pair sign oscillation due to B0(s) oscillation for the same mother.
623668
for (const auto& mid1 : mothers_id1) {
624669
for (const auto& mid2 : mothers_id2) {
625670
if (mid1 == mid2) {
626671
auto common_mp = mcparticles.iteratorAt(mid1);
627672
int mp_pdg = common_mp.pdgCode();
628-
bool is_mp_diquark = (1100 < std::abs(mp_pdg) && std::abs(mp_pdg) < 5600) && std::to_string(mp_pdg)[std::to_string(mp_pdg).length() - 2] == '0';
629-
if (!is_mp_diquark && std::abs(mp_pdg) < 1e+9 && (std::to_string(std::abs(mp_pdg))[std::to_string(std::abs(mp_pdg)).length() - 3] == '5' || std::to_string(std::abs(mp_pdg))[std::to_string(std::abs(mp_pdg)).length() - 4] == '5')) {
673+
if (std::abs(mp_pdg) < 1e+9 && (std::to_string(std::abs(mp_pdg))[std::to_string(std::abs(mp_pdg)).length() - 3] == '5' || std::to_string(std::abs(mp_pdg))[std::to_string(std::abs(mp_pdg)).length() - 4] == '5')) {
630674
mothers_id1.clear();
631675
mothers_pdg1.clear();
632676
mothers_id2.clear();
@@ -635,27 +679,25 @@ int IsHF(TMCParticle1 const& p1, TMCParticle2 const& p2, TMCParticles const& mcp
635679
mothers_pdg1.shrink_to_fit();
636680
mothers_id2.shrink_to_fit();
637681
mothers_pdg2.shrink_to_fit();
638-
return static_cast<int>(EM_HFeeType::kBCe_Be_SameB); // b->c->e and b->e, decay type = 3
682+
return static_cast<int>(EM_HFeeType::kBCe_Be_SameB); // b->c->e and b->e, decay type = 3 // No pair sign oscillation due to B0(s) oscillation for the same mother.
639683
}
640684
}
641685
} // end of motherid2
642686
} // end of motherid1
643687

644-
bool b2c2l_b2l_diffb_case0 = mpfh1.pdgCode() * mpfh2.pdgCode() < 0 && !isFOat1stDecay1 && !isFOat1stDecay2; // bbbar -> ll LS
645-
bool b2c2l_b2l_diffb_case1 = mpfh1.pdgCode() * mpfh2.pdgCode() > 0 && static_cast<bool>(isFOat1stDecay1 ^ isFOat1stDecay2); // bbbar -> ll ULS
646-
bool b2c2l_b2l_diffb_case2 = mpfh1.pdgCode() * mpfh2.pdgCode() < 0 && isFOat1stDecay1 && isFOat1stDecay2; // bbbar -> ll LS
647-
648-
if (b2c2l_b2l_diffb_case0 || b2c2l_b2l_diffb_case1 || b2c2l_b2l_diffb_case2) {
649-
mothers_id1.clear();
650-
mothers_pdg1.clear();
651-
mothers_id2.clear();
652-
mothers_pdg2.clear();
653-
mothers_id1.shrink_to_fit();
654-
mothers_pdg1.shrink_to_fit();
655-
mothers_id2.shrink_to_fit();
656-
mothers_pdg2.shrink_to_fit();
657-
return static_cast<int>(EM_HFeeType::kBCe_Be_DiffB); // b->c->e and b->e, decay type = 4
658-
}
688+
// bool b2c2l_b2l_diffb_case0 = mpfh1.pdgCode() * mpfh2.pdgCode() < 0 && !isFOat1stDecay1 && !isFOat1stDecay2; // bbbar -> ll LS
689+
// bool b2c2l_b2l_diffb_case1 = mpfh1.pdgCode() * mpfh2.pdgCode() > 0 && static_cast<bool>(isFOat1stDecay1 ^ isFOat1stDecay2); // bbbar -> ll ULS
690+
// bool b2c2l_b2l_diffb_case2 = mpfh1.pdgCode() * mpfh2.pdgCode() < 0 && isFOat1stDecay1 && isFOat1stDecay2; // bbbar -> ll LS
691+
692+
mothers_id1.clear();
693+
mothers_pdg1.clear();
694+
mothers_id2.clear();
695+
mothers_pdg2.clear();
696+
mothers_id1.shrink_to_fit();
697+
mothers_pdg1.shrink_to_fit();
698+
mothers_id2.shrink_to_fit();
699+
mothers_pdg2.shrink_to_fit();
700+
return static_cast<int>(EM_HFeeType::kBCe_Be_DiffB); // b->c->e and b->e, decay type = 4
659701
}
660702

661703
mothers_id1.clear();

0 commit comments

Comments
 (0)