Skip to content
Merged
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
32 changes: 32 additions & 0 deletions PWGEM/PhotonMeson/TableProducer/associateMCinfoPhoton.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -65,6 +65,7 @@ struct AssociateMCInfoPhoton {
Configurable<float> margin_z_gen{"margin_z_gen", 15.f, "margin for Z of true photon conversion point to store generated information"};
Configurable<float> max_rxy_gen{"max_rxy_gen", 100, "max rxy to store generated information"};
Configurable<bool> requireGammaGammaDecay{"requireGammaGammaDecay", false, "require gamma gamma decay for generated pi0 and eta meson"};
Configurable<float> cfg_max_eta_photon{"cfg_max_eta_gamma", 0.8, "max eta for photons at PV"};

HistogramRegistry registry{"EMMCEvent"};

Expand Down Expand Up @@ -98,6 +99,7 @@ struct AssociateMCInfoPhoton {

for (uint i = 0; i < NParticleNames; i++) {
registry.add<TH2>(Form("Generated/h2PtY_%s", ParticleNames[i].data()), Form("Generated %s", ParticleNames[i].data()), kTH2F, {axisPt, axisRapidity}, true);
registry.add<TH2>(Form("Generated/h2PtY_Accepted_%s", ParticleNames[i].data()), Form("Generated %s", ParticleNames[i].data()), kTH2F, {axisPt, axisRapidity}, true);
}

// reserve space for generated vectors if that process enabled
Expand Down Expand Up @@ -170,6 +172,15 @@ struct AssociateMCInfoPhoton {
for (const auto& mcParticle : groupedMcParticles) { // store necessary information for denominator of efficiency
if ((mcParticle.isPhysicalPrimary() || mcParticle.producedByGenerator()) && std::fabs(mcParticle.y()) < 0.9f && mcParticle.pt() < 20.f) {
auto binNumber = hBinFinder->FindBin(mcParticle.pt(), std::fabs(mcParticle.y())); // caution: pack

bool isMesonAccepted = false;
if (mcParticle.has_daughters()) {
auto lDaughters = mcParticle.daughters_as<aod::McParticles>();

if (areTwoPhotonDaughtersAccepted(lDaughters, cfg_max_eta_photon)) {
isMesonAccepted = true;
}
}
switch (std::abs(mcParticle.pdgCode())) {
case PDG_t::kGamma:
registry.fill(HIST("Generated/h2PtY_Gamma"), mcParticle.pt(), std::fabs(mcParticle.y()));
Expand All @@ -180,12 +191,18 @@ struct AssociateMCInfoPhoton {
continue;
registry.fill(HIST("Generated/h2PtY_Pi0"), mcParticle.pt(), std::fabs(mcParticle.y()));
genPi0[binNumber]++;
if (isMesonAccepted) {
registry.fill(HIST("Generated/h2PtY_Accepted_Pi0"), mcParticle.pt(), std::fabs(mcParticle.y()));
}
break;
case Pdg::kEta:
if (requireGammaGammaDecay && !isGammaGammaDecay(mcParticle, mcParticles))
continue;
registry.fill(HIST("Generated/h2PtY_Eta"), mcParticle.pt(), std::fabs(mcParticle.y()));
genEta[binNumber]++;
if (isMesonAccepted) {
registry.fill(HIST("Generated/h2PtY_Accepted_Eta"), mcParticle.pt(), std::fabs(mcParticle.y()));
}
break;
default:
break;
Expand Down Expand Up @@ -492,6 +509,21 @@ struct AssociateMCInfoPhoton {
fCounters[1] = 0;
} // end of skimmingMC

template <typename Daughters>
inline bool areTwoPhotonDaughtersAccepted(const Daughters& lDaughters, float maxEta)
{
if (lDaughters.size() != 2) {
return false;
}

auto it0 = lDaughters.begin();
auto it1 = it0;
++it1;

return (std::fabs(it0.eta()) < maxEta &&
std::fabs(it1.eta()) < maxEta);
}

void processMC_PCM(MyCollisionsMC const& collisions, aod::BCs const& bcs, aod::McCollisions const& mccollisions, aod::McParticles const& mcParticles, TracksMC const& o2tracks, aod::V0PhotonsKF const& v0photons, aod::V0Legs const& v0legs)
{
skimmingMC<kPCM>(collisions, bcs, mccollisions, mcParticles, o2tracks, nullptr, v0photons, v0legs, nullptr, nullptr, nullptr);
Expand Down
Loading