Resolve "baryon and hadron number as particle property"

Framework/Particles/ParticleProperties.h

@@ -51,6 +51,10 @@ namespace corsika::particles {
   corsika::units::si::TimeType constexpr GetLifetime(Code const);
 
   bool constexpr IsNucleus(Code const);
+  bool constexpr IsHadron(Code const);
+  bool constexpr IsEM(Code const);
+  bool constexpr IsMuon(Code const);
+  bool constexpr IsNeutrino(Code const);
   int constexpr GetNucleusA(Code const);
   int constexpr GetNucleusZ(Code const);
 
@@ -100,6 +104,21 @@ namespace corsika::particles {
     return detail::lifetime[static_cast<CodeIntType>(p)] * corsika::units::si::second;
   }
 
+  bool constexpr IsHadron(Code const p) {
+    return detail::isHadron[static_cast<CodeIntType>(p)];
+  }
+
+  bool constexpr IsEM(Code c) {
+    return c == Code::Electron || c == Code::Positron || c == Code::Gamma;
+  }
+
+  bool constexpr IsMuon(Code c) { return c == Code::MuPlus || c == Code::MuMinus; }
+
+  bool constexpr IsNeutrino(Code c) {
+    return c == Code::NuE || c == Code::NuMu || c == Code::NuTau || c == Code::NuEBar ||
+           c == Code::NuMuBar || c == Code::NuTauBar;
+  }
+
   bool constexpr IsNucleus(Code const p) {
     return detail::isNucleus[static_cast<CodeIntType>(p)];
   }

Framework/Particles/pdxml_reader.py

@@ -171,7 +171,11 @@ def read_pythia_db(filename, particle_db, classnames):
             c_id = classnames[pdg]
         else:
             c_id = c_identifier_camel(name) # the camel case names
-        
+
+        hadron = False
+        if abs(pdg) > 100:
+            hadron = True
+            
         particle_db[c_id] = {
             "name" : name,
             "antiName" : antiName,
@@ -181,6 +185,7 @@ def read_pythia_db(filename, particle_db, classnames):
             "lifetime" : lifetime,
             "ngc_code" : next(counter),
             "isNucleus" : False,
+            "isHadron" : hadron,
         }
     
     return particle_db
@@ -214,6 +219,7 @@ def read_nuclei_db(filename, particle_db, classnames):
             "A" : A,
             "Z" : Z,
             "isNucleus" : True,
+            "isHadron" : True,
         }
     
     return particle_db
@@ -350,6 +356,15 @@ def gen_properties(particle_db):
             string += "  {tau:e}, \n".format(tau = p['lifetime'])
             #string += "  {tau:e} * corsika::units::si::second, \n".format(tau = p['lifetime'])
     string += "};\n"
+    
+    # is Hadron flag
+    string += "static constexpr std::array<bool, size> isHadron = {\n"
+    for p in particle_db.values():
+        value = 'false'
+        if p['isHadron']:
+            value = 'true'
+        string += "  {val},\n".format(val = value)
+    string += "};\n"
 
     
     ### nuclear data ###

Framework/Particles/testParticles.cc

@@ -84,6 +84,47 @@ TEST_CASE("ParticleProperties", "[Particles]") {
             (Approx(2.1970332555864364e-06).epsilon(1e-5)));
   }
 
+  SECTION("Particle groups: electromagnetic") {
+    REQUIRE(IsEM(Code::Gamma));
+    REQUIRE(IsEM(Code::Electron));
+    REQUIRE_FALSE(IsEM(Code::MuPlus));
+    REQUIRE_FALSE(IsEM(Code::NuE));
+    REQUIRE_FALSE(IsEM(Code::Proton));
+    REQUIRE_FALSE(IsEM(Code::PiPlus));
+    REQUIRE_FALSE(IsEM(Code::Oxygen));
+  }
+
+  SECTION("Particle groups: hadrons") {
+    REQUIRE_FALSE(IsHadron(Code::Gamma));
+    REQUIRE_FALSE(IsHadron(Code::Electron));
+    REQUIRE_FALSE(IsHadron(Code::MuPlus));
+    REQUIRE_FALSE(IsHadron(Code::NuE));
+    REQUIRE(IsHadron(Code::Proton));
+    REQUIRE(IsHadron(Code::PiPlus));
+    REQUIRE(IsHadron(Code::Oxygen));
+  }
+
+  SECTION("Particle groups: muons") {
+    REQUIRE_FALSE(IsMuon(Code::Gamma));
+    REQUIRE_FALSE(IsMuon(Code::Electron));
+    REQUIRE(IsMuon(Code::MuPlus));
+    REQUIRE_FALSE(IsMuon(Code::NuE));
+    REQUIRE_FALSE(IsMuon(Code::Proton));
+    REQUIRE_FALSE(IsMuon(Code::PiPlus));
+    REQUIRE_FALSE(IsMuon(Code::Oxygen));
+  }
+
+  SECTION("Particle groups: neutrinos") {
+    REQUIRE_FALSE(IsNeutrino(Code::Gamma));
+    REQUIRE_FALSE(IsNeutrino(Code::Electron));
+    REQUIRE_FALSE(IsNeutrino(Code::MuPlus));
+    REQUIRE(IsNeutrino(Code::NuE));
+    REQUIRE_FALSE(IsNeutrino(Code::Proton));
+    REQUIRE_FALSE(IsNeutrino(Code::PiPlus));
+    REQUIRE_FALSE(IsNeutrino(Code::Oxygen));
+  }
+
+
   SECTION("Nuclei") {
     REQUIRE_FALSE(IsNucleus(Code::Gamma));
     REQUIRE(IsNucleus(Code::Argon));