[PWGEM/PhotonMeson] Add clusterizer selection to EM tasks

PWGEM/PhotonMeson/Core/EMCPhotonCut.cxx

@@ -15,11 +15,17 @@
 
 #include "Framework/Logger.h"
 #include "PWGEM/PhotonMeson/Core/EMCPhotonCut.h"
+#include "PWGJE/DataModel/EMCALClusters.h"
 
 ClassImp(EMCPhotonCut);
 
-const char* EMCPhotonCut::mCutNames[static_cast<int>(EMCPhotonCut::EMCPhotonCuts::kNCuts)] = {"Energy", "NCell", "M02", "Timing", "TrackMatching", "Exotic"};
+const char* EMCPhotonCut::mCutNames[static_cast<int>(EMCPhotonCut::EMCPhotonCuts::kNCuts)] = {"Definition", "Energy", "NCell", "M02", "Timing", "TrackMatching", "Exotic"};
 
+void EMCPhotonCut::SetClusterizer(std::string clusterDefinitionString)
+{
+  mDefinition = static_cast<int>(o2::aod::emcalcluster::getClusterDefinitionFromString(clusterDefinitionString));
+  LOG(info) << "EMCal Photon Cut, set cluster definition to: " << mDefinition << " (" << clusterDefinitionString << ")";
+}
 void EMCPhotonCut::SetMinE(float min)
 {
   mMinE = min;
@@ -72,6 +78,9 @@ void EMCPhotonCut::print() const
   LOG(info) << "EMCal Photon Cut:";
   for (int i = 0; i < static_cast<int>(EMCPhotonCuts::kNCuts); i++) {
     switch (static_cast<EMCPhotonCuts>(i)) {
+      case EMCPhotonCuts::kDefinition:
+        LOG(info) << mCutNames[i] << " > " << mDefinition;
+        break;
       case EMCPhotonCuts::kEnergy:
         LOG(info) << mCutNames[i] << " > " << mMinE;
         break;

PWGEM/PhotonMeson/Core/EMCPhotonCut.h

@@ -34,7 +34,8 @@ class EMCPhotonCut : public TNamed
 
   enum class EMCPhotonCuts : int {
     // cluster cut
-    kEnergy = 0,
+    kDefinition = 0,
+    kEnergy,
     kNCell,
     kM02,
     kTiming,
@@ -49,6 +50,9 @@ class EMCPhotonCut : public TNamed
   template <typename T, typename Cluster>
   bool IsSelected(Cluster const& cluster) const
   {
+    if (!IsSelectedEMCal(EMCPhotonCuts::kDefinition, cluster)) {
+      return false;
+    }
     if (!IsSelectedEMCal(EMCPhotonCuts::kEnergy, cluster)) {
       return false;
     }
@@ -75,6 +79,9 @@ class EMCPhotonCut : public TNamed
   bool IsSelectedEMCal(const EMCPhotonCuts& cut, Cluster const& cluster) const
   {
     switch (cut) {
+      case EMCPhotonCuts::kDefinition:
+        return cluster.definition() == mDefinition;
+
       case EMCPhotonCuts::kEnergy:
         return cluster.e() > mMinE;
 
@@ -113,6 +120,7 @@ class EMCPhotonCut : public TNamed
   }
 
   // Setters
+  void SetClusterizer(std::string clusterDefinitionString = "kV3Default");
   void SetMinE(float min = 0.7f);
   void SetMinNCell(int min = 1);
   void SetM02Range(float min = 0.1f, float max = 0.7f);
@@ -128,6 +136,7 @@ class EMCPhotonCut : public TNamed
 
  private:
   // EMCal cluster cuts
+  int mDefinition{10};      ///< clusterizer definition
   float mMinE{0.7f};        ///< minimum energy
   int mMinNCell{1};         ///< minimum number of cells per cluster
   float mMinM02{0.1f};      ///< minimum M02 for a cluster

PWGEM/PhotonMeson/Core/Pi0EtaToGammaGamma.h

@@ -180,6 +180,7 @@ struct Pi0EtaToGammaGamma {
   EMCPhotonCut fEMCCut;
   struct : ConfigurableGroup {
     std::string prefix = "emccut_group";
+    Configurable<std::string> clusterDefinition{"clusterDefinition", "kV3Default", "Clusterizer to be selected, e.g. V3Default"};
     Configurable<float> minOpenAngle{"minOpenAngle", 0.0202, "apply min opening angle"};
     Configurable<float> EMC_minTime{"EMC_minTime", -20., "Minimum cluster time for EMCal time cut"};
     Configurable<float> EMC_maxTime{"EMC_maxTime", +25., "Maximum cluster time for EMCal time cut"};
@@ -422,6 +423,7 @@ struct Pi0EtaToGammaGamma {
 
     fEMCCut = EMCPhotonCut("fEMCCut", "fEMCCut");
 
+    fEMCCut.SetClusterizer(emccuts.clusterDefinition);
     fEMCCut.SetMinE(emccuts.EMC_minE);
     fEMCCut.SetMinNCell(emccuts.EMC_minNCell);
     fEMCCut.SetM02Range(emccuts.EMC_minM02, emccuts.EMC_maxM02);

PWGEM/PhotonMeson/Core/Pi0EtaToGammaGammaMC.h

@@ -168,6 +168,7 @@ struct Pi0EtaToGammaGammaMC {
   EMCPhotonCut fEMCCut;
   struct : ConfigurableGroup {
     std::string prefix = "emccut_group";
+    Configurable<std::string> clusterDefinition{"clusterDefinition", "kV3Default", "Clusterizer to be selected, e.g. V3Default"};
     Configurable<float> minOpenAngle{"minOpenAngle", 0.0202, "apply min opening angle"};
     Configurable<float> EMC_minTime{"EMC_minTime", -20., "Minimum cluster time for EMCal time cut"};
     Configurable<float> EMC_maxTime{"EMC_maxTime", +25., "Maximum cluster time for EMCal time cut"};
@@ -385,6 +386,7 @@ struct Pi0EtaToGammaGammaMC {
 
     fEMCCut = EMCPhotonCut("fEMCCut", "fEMCCut");
 
+    fEMCCut.SetClusterizer(emccuts.clusterDefinition);
     fEMCCut.SetMinE(emccuts.EMC_minE);
     fEMCCut.SetMinNCell(emccuts.EMC_minNCell);
     fEMCCut.SetM02Range(emccuts.EMC_minM02, emccuts.EMC_maxM02);

PWGEM/PhotonMeson/DataModel/gammaTables.h

@@ -459,7 +459,7 @@ DECLARE_SOA_COLUMN(TrackPt, trackpt, std::vector<float>);
 DECLARE_SOA_DYNAMIC_COLUMN(Pt, pt, [](float e, float eta, float m = 0) -> float { return sqrt(e * e - m * m) / cosh(eta); }); //! cluster pt, mass to be given as argument when getter is called!
 } // namespace emccluster
 DECLARE_SOA_TABLE(SkimEMCClusters, "AOD", "SKIMEMCCLUSTERS", //! table of skimmed EMCal clusters
-                  o2::soa::Index<>, skimmedcluster::CollisionId, skimmedcluster::E, skimmedcluster::Eta, skimmedcluster::Phi,
+                  o2::soa::Index<>, skimmedcluster::CollisionId, emccluster::Definition, skimmedcluster::E, skimmedcluster::Eta, skimmedcluster::Phi,
                   skimmedcluster::M02, skimmedcluster::NCells, skimmedcluster::Time, emccluster::IsExotic, emccluster::TrackEta,
                   emccluster::TrackPhi, emccluster::TrackP, emccluster::TrackPt, emccluster::Pt<skimmedcluster::E, skimmedcluster::Eta>);
 using SkimEMCCluster = SkimEMCClusters::iterator;

PWGEM/PhotonMeson/TableProducer/skimmerGammaCalo.cxx

@@ -137,7 +137,7 @@ struct skimmerGammaCalo {
       historeg.fill(HIST("hCaloClusterEOut"), emccluster.energy());
       historeg.fill(HIST("hCaloClusterFilter"), 4);
 
-      tableGammaEMCReco(emccluster.collisionId(), emccluster.energy(), emccluster.eta(), emccluster.phi(), emccluster.m02(),
+      tableGammaEMCReco(emccluster.collisionId(), emccluster.definition(), emccluster.energy(), emccluster.eta(), emccluster.phi(), emccluster.m02(),
                         emccluster.nCells(), emccluster.time(), emccluster.isExotic(), vEta, vPhi, vP, vPt);
     }
   }

PWGEM/PhotonMeson/Tasks/emcalQC.cxx

@@ -76,6 +76,7 @@ struct emcalQC {
   EMCPhotonCut fEMCCut;
   struct : ConfigurableGroup {
     std::string prefix = "emccut_group";
+    Configurable<std::string> clusterDefinition{"clusterDefinition", "kV3Default", "Clusterizer to be selected, e.g. V3Default"};
     Configurable<float> minOpenAngle{"minOpenAngle", 0.0202, "apply min opening angle"};
     Configurable<float> EMC_minTime{"EMC_minTime", -20., "Minimum cluster time for EMCal time cut"};
     Configurable<float> EMC_maxTime{"EMC_maxTime", +25., "Maximum cluster time for EMCal time cut"};
@@ -115,6 +116,7 @@ struct emcalQC {
     const float f = emccuts.EMC_TM_Phi->at(2);
     LOGF(info, "EMCal track matching parameters : a = %f, b = %f, c = %f, d = %f, e = %f, f = %f", a, b, c, d, e, f);
 
+    fEMCCut.SetClusterizer(emccuts.clusterDefinition);
     fEMCCut.SetMinE(emccuts.EMC_minE);
     fEMCCut.SetMinNCell(emccuts.EMC_minNCell);
     fEMCCut.SetM02Range(emccuts.EMC_minM02, emccuts.EMC_maxM02);
@@ -193,9 +195,13 @@ struct emcalQC {
 
       auto clusters_per_coll = clusters.sliceBy(perCollision, collision.collisionId());
       fRegistry.fill(HIST("Cluster/before/hNgamma"), clusters_per_coll.size(), collision.weight());
-      int ng = 0;
+      int ngBefore = 0;
+      int ngAfter = 0;
       for (auto& cluster : clusters_per_coll) {
         // Fill the cluster properties before applying any cuts
+        if (!fEMCCut.IsSelectedEMCal(EMCPhotonCut::EMCPhotonCuts::kDefinition, cluster))
+          continue;
+        ngBefore++;
         o2::aod::pwgem::photonmeson::utils::clusterhistogram::fillClusterHistograms<0>(&fRegistry, cluster, cfgDo2DQA, collision.weight());
 
         // Apply cuts one by one and fill in hClusterQualityCuts histogram
@@ -205,10 +211,10 @@ struct emcalQC {
         bool survivesIsSelectedEMCalCuts = true;                        // Survives "manual" cuts listed in this task
         bool survivesIsSelectedCuts = fEMCCut.IsSelected<int>(cluster); // Survives the cutlist defines in EMCPhotonCut.h, which is also used in the Pi0Eta task
 
-        for (int icut = 0; icut < static_cast<int>(EMCPhotonCut::EMCPhotonCuts::kNCuts); icut++) { // Loop through different cut observables
+        for (int icut = 1; icut < static_cast<int>(EMCPhotonCut::EMCPhotonCuts::kNCuts); icut++) { // Loop through different cut observables, start at 1 to ignore ClusterDefinition
           EMCPhotonCut::EMCPhotonCuts specificcut = static_cast<EMCPhotonCut::EMCPhotonCuts>(icut);
           if (!fEMCCut.IsSelectedEMCal(specificcut, cluster)) { // Check whether cluster passes this cluster requirement, if not, fill why in the next row
-            fRegistry.fill(HIST("Cluster/hClusterQualityCuts"), icut + 1, cluster.e(), collision.weight());
+            fRegistry.fill(HIST("Cluster/hClusterQualityCuts"), icut, cluster.e(), collision.weight());
             survivesIsSelectedEMCalCuts = false;
           }
         }
@@ -220,10 +226,11 @@ struct emcalQC {
         if (survivesIsSelectedCuts) {
           o2::aod::pwgem::photonmeson::utils::clusterhistogram::fillClusterHistograms<1>(&fRegistry, cluster, cfgDo2DQA, collision.weight());
           fRegistry.fill(HIST("Cluster/hClusterQualityCuts"), 7., cluster.e(), collision.weight());
-          ng++;
+          ngAfter++;
         }
       }
-      fRegistry.fill(HIST("Cluster/after/hNgamma"), ng, collision.weight());
+      fRegistry.fill(HIST("Cluster/before/hNgamma"), ngBefore, collision.weight());
+      fRegistry.fill(HIST("Cluster/after/hNgamma"), ngAfter, collision.weight());
     } // end of collision loop
   } // end of process