PWGEM/PhotonMeson: reduce skimmed data size (#5219)

Author: dsekihat

PWGEM/PhotonMeson/Core/HistogramsLibrary.cxx

@@ -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") {
@@ -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") {

PWGEM/PhotonMeson/TableProducer/associateMCinfo.cxx

@@ -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"};
@@ -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;
@@ -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]++;
       }
@@ -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();
@@ -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)
@@ -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);