Skip to content

PWGEM/PhotonMeson: reduce skimmed MC data size #5219

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
dsekihat merged 1 commit into from
Mar 18, 2024
Merged

PWGEM/PhotonMeson: reduce skimmed MC data size #5219

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
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
8 changes: 4 additions & 4 deletions PWGEM/PhotonMeson/Core/HistogramsLibrary.cxx
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -490,9 +490,9 @@ void o2::aod::pwgem::photon::histogram::DefineHistograms(THashList* list, const
list->Add(new TH1F("hPt_ConvertedPhoton", "converted photon pT;p_{T} (GeV/c)", 2000, 0.0f, 20));
list->Add(new TH1F("hY_ConvertedPhoton", "converted photon y;rapidity y", 40, -2.0f, 2.0f));
list->Add(new TH1F("hPhi_ConvertedPhoton", "converted photon #varphi;#varphi (rad.)", 180, 0, 2 * M_PI));
list->Add(new TH2F("hPhotonRxy", "conversion point in XY MC;V_{x} (cm);V_{y} (cm)", 2000, -100.0f, 100.0f, 2000, -100.0f, 100.0f));
list->Add(new TH2F("hPhotonRZ", "conversion point in RZ MC;V_{z} (cm);R_{xy} (cm)", 2000, -100.0f, 100.0f, 1000, 0.f, 100.0f));
list->Add(new TH2F("hPhotonPhivsRxy", "conversion point of #varphi vs. R_{xy} MC;#varphi (rad.);R_{xy} (cm);N_{e}", 360, 0.0f, 2 * M_PI, 100, 0, 100));
list->Add(new TH2F("hPhotonRxy", "conversion point in XY MC;V_{x} (cm);V_{y} (cm)", 800, -100.0f, 100.0f, 800, -100.0f, 100.0f));
list->Add(new TH2F("hPhotonRZ", "conversion point in RZ MC;V_{z} (cm);R_{xy} (cm)", 1000, -100.0f, 100.0f, 400, 0.f, 100.0f));
list->Add(new TH2F("hPhotonPhivsRxy", "conversion point of #varphi vs. R_{xy} MC;#varphi (rad.);R_{xy} (cm);N_{e}", 360, 0.0f, 2 * M_PI, 400, 0, 100));
}

if (TString(subGroup) == "Pi0Eta") {
Expand All @@ -507,7 +507,7 @@ void o2::aod::pwgem::photon::histogram::DefineHistograms(THashList* list, const
}
}
if (TString(subGroup) == "dielectron") {
TH2F* hMvsPt = new TH2F("hMvsPt", "m_{ee} vs. p_{T,ee};m_{ee} (GeV/c^{2});p_{T,ee} (GeV/c)", 110, 0, 1.1f, 1000, 0, 10.f);
TH2F* hMvsPt = new TH2F("hMvsPt", "m_{ee} vs. p_{T,ee};m_{ee} (GeV/c^{2});p_{T,ee} (GeV/c)", 400, 0, 4.0f, 1000, 0, 10.f);
hMvsPt->Sumw2();
list->Add(hMvsPt);
} else if (TString(subGroup) == "dimuon") {
Expand Down
40 changes: 26 additions & 14 deletions PWGEM/PhotonMeson/TableProducer/associateMCinfo.cxx
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -48,7 +48,7 @@ struct AssociateMCInfo {
Produces<o2::aod::EMPrimaryMuonMCLabels> emprimarymuonmclabels;

Configurable<float> max_rxy_gen{"max_rxy_gen", 100, "max rxy to store generated information"};
Configurable<float> max_Y_gen{"max_Y_gen", 1.5, "max rapidity Y to store generated information"};
Configurable<float> max_Y_gen{"max_Y_gen", 0.9, "max rapidity Y to store generated information"};
Configurable<float> margin_z_gen{"margin_z_gen", 15.f, "margin for Z of true photon conversion point to store generated information"};

HistogramRegistry registry{"EMMCEvent"};
Expand All @@ -58,6 +58,8 @@ struct AssociateMCInfo {
auto hEventCounter = registry.add<TH1>("hEventCounter", "hEventCounter", kTH1I, {{6, 0.5f, 6.5f}});
hEventCounter->GetXaxis()->SetBinLabel(1, "all");
hEventCounter->GetXaxis()->SetBinLabel(2, "has mc collision");
registry.add<TH2>("PCM/hXY", "hRxy;X (cm);Y (cm)", kTH2F, {{400, -100, +100}, {400, -100, +100}});
registry.add<TH2>("PCM/hRZ", "hRxy;R (cm);Z (cm)", kTH2F, {{400, -100, +100}, {200, 0, +100}});
}

Preslice<aod::McParticles> perMcCollision = aod::mcparticle::mcCollisionId;
Expand Down Expand Up @@ -92,7 +94,6 @@ struct AssociateMCInfo {
// make an entry for this MC event only if it was not already added to the table
if (!(fEventLabels.find(mcCollision.globalIndex()) != fEventLabels.end())) {
mcevents(mcCollision.globalIndex(), mcCollision.generatorsID(), mcCollision.posX(), mcCollision.posY(), mcCollision.posZ(), mcCollision.t(), mcCollision.impactParameter());
// mcevents(mcCollision.generatorsID(), mcCollision.posX(), mcCollision.posY(), mcCollision.posZ(), mcCollision.t(), mcCollision.impactParameter());
fEventLabels[mcCollision.globalIndex()] = fCounters[1];
fCounters[1]++;
}
Expand All @@ -103,7 +104,7 @@ struct AssociateMCInfo {
auto groupedMcTracks = mcTracks.sliceBy(perMcCollision, mcCollision.globalIndex());

for (auto& mctrack : groupedMcTracks) { // store necessary information for denominator of efficiency
if (mctrack.pt() < 1e-3 || abs(mctrack.y()) > max_Y_gen || abs(mctrack.vz()) > 250 || sqrt(pow(mctrack.vx(), 2) + pow(mctrack.vy(), 2)) > max_rxy_gen) {
if (mctrack.pt() < 1e-3 || abs(mctrack.vz()) > 250 || sqrt(pow(mctrack.vx(), 2) + pow(mctrack.vy(), 2)) > max_rxy_gen) {
continue;
}
int pdg = mctrack.pdgCode();
Expand Down Expand Up @@ -138,14 +139,30 @@ struct AssociateMCInfo {
}
// LOGF(info,"index = %d , mc track pdg = %d , producedByGenerator = %d , isPhysicalPrimary = %d", mctrack.index(), mctrack.pdgCode(), mctrack.producedByGenerator(), mctrack.isPhysicalPrimary());

if (abs(pdg) == 11) { // electrons // extra check of production vertex for secondary electrons to reduce data size.
if (!mctrack.isPhysicalPrimary() && !mctrack.producedByGenerator() && sqrt(pow(mctrack.vx(), 2) + pow(mctrack.vy(), 2)) < abs(mctrack.vz()) * std::tan(2 * std::atan(std::exp(-max_Y_gen))) - margin_z_gen) {
continue;
}
} else {
if (!mctrack.isPhysicalPrimary() && !mctrack.producedByGenerator()) {
if (!(mctrack.isPhysicalPrimary() || mctrack.producedByGenerator())) { // neither physical primary nor producedByGenerator
if (abs(pdg) == 11) { // one more check for secondary electrons. i.e. gamma->ee

if (sqrt(pow(mctrack.vx(), 2) + pow(mctrack.vy(), 2)) < abs(mctrack.vz()) * std::tan(2 * std::atan(std::exp(-max_Y_gen))) - margin_z_gen) {
continue;
}

if (mctrack.has_mothers()) {
auto mp = mctrack.template mothers_first_as<aod::McParticles>(); // mother particle of electron
int pdg_mother = mp.pdgCode();
if (pdg_mother != 22 || !(mp.isPhysicalPrimary() || mp.producedByGenerator())) { // mother of electron is not photon, or not physical primary, or not producedByGenerator
continue;
}
}
} else { // not physical primary, not producedByGenerator, not electrons
continue;
}
} else if (abs(mctrack.y()) > max_Y_gen) { // physical primary or producedByGenerator, but outside of acceptance.
continue;
}

if (abs(pdg) == 11 && !(mctrack.isPhysicalPrimary() || mctrack.producedByGenerator())) { // only for quick check, only secondary electrons should appear.
registry.fill(HIST("PCM/hXY"), mctrack.vx(), mctrack.vy());
registry.fill(HIST("PCM/hRZ"), mctrack.vz(), sqrt(pow(mctrack.vx(), 2) + pow(mctrack.vy(), 2)));
}

// these are used as denominator for efficiency. (i.e. generated information)
Expand Down Expand Up @@ -326,12 +343,7 @@ struct AssociateMCInfo {

std::vector<int> mothers;
if (mctrack.has_mothers()) {
// LOGF(info, "mother ids size = %d", mctrack.mothersIds().size());
// int counter = 0;
for (auto& m : mctrack.mothersIds()) {
// LOGF(info, "counter = %d , mother id = %d", counter, m);
// counter++;

if (m < mcTracks.size()) { // protect against bad mother indices
if (fNewLabels.find(m) != fNewLabels.end()) {
mothers.push_back(fNewLabels.find(m)->second);
Expand Down