fixed small problems and examples

corsika/detail/modules/sibyll/Decay.inl

@@ -35,11 +35,12 @@ namespace corsika::sibyll {
 
   inline bool Decay::canHandleDecay(const Code vParticleCode) {
     // if known to sibyll and not proton or neutrino it can decay
-    if (vParticleCode == Code::Proton || vParticleCode == Code::AntiProton ||
-        vParticleCode == Code::NuE || vParticleCode == Code::NuMu ||
-        vParticleCode == Code::NuTau || vParticleCode == Code::NuEBar ||
-        vParticleCode == Code::NuMuBar || vParticleCode == Code::NuTauBar ||
-        vParticleCode == Code::Electron || vParticleCode == Code::Positron)
+    if (is_nucleus(vParticleCode) || vParticleCode == Code::Proton ||
+        vParticleCode == Code::AntiProton || vParticleCode == Code::NuE ||
+        vParticleCode == Code::NuMu || vParticleCode == Code::NuTau ||
+        vParticleCode == Code::NuEBar || vParticleCode == Code::NuMuBar ||
+        vParticleCode == Code::NuTauBar || vParticleCode == Code::Electron ||
+        vParticleCode == Code::Positron)
       return false;
     else if (corsika::sibyll::convertToSibyllRaw(
                  vParticleCode)) // non-zero for particles known to sibyll

corsika/modules/sibyll/ParticleConversion.hpp

@@ -48,13 +48,13 @@ namespace corsika::sibyll {
     return corsikaCode;
   }
 
-  int constexpr convertToSibyllRaw(corsika::Code pCode) {
-    return static_cast<int>(convertToSibyll(pCode));
+  int constexpr convertToSibyllRaw(Code const code) {
+    return static_cast<int>(convertToSibyll(code));
   }
 
-  int constexpr getSibyllXSCode(corsika::Code pCode) {
+  int constexpr getSibyllXSCode(Code const code) {
     return static_cast<SibyllXSClassIntType>(
-        corsika2sibyllXStype[static_cast<corsika::CodeIntType>(pCode)]);
+        corsika2sibyllXStype[static_cast<corsika::CodeIntType>(code)]);
   }
 
   bool constexpr canInteract(corsika::Code pCode) { return getSibyllXSCode(pCode) > 0; }

examples/cascade_example.cpp

@@ -8,6 +8,7 @@
 
 #include <corsika/framework/core/Cascade.hpp>
 #include <corsika/framework/process/ProcessSequence.hpp>
+#include <corsika/framework/process/SwitchProcessSequence.hpp>
 #include <corsika/framework/core/PhysicalUnits.hpp>
 #include <corsika/framework/random/RNGManager.hpp>
 #include <corsika/framework/geometry/Sphere.hpp>
@@ -146,8 +147,11 @@ int main() {
   BetheBlochPDG eLoss{showerAxis};
 
   // assemble all processes into an ordered process list
-  auto sequence =
-      make_sequence(stackInspect, sibyll, sibyllNuc, decay, eLoss, cut, trackWriter);
+  auto sequence = make_sequence(
+      stackInspect,
+      make_select([](auto const& particle) { return is_nucleus(particle.getPID()); },
+                  sibyllNuc, sibyll),
+      decay, eLoss, cut, trackWriter);
 
   // define air shower object, run simulation
   Cascade EAS(env, tracking, sequence, output, stack);

examples/corsika.cpp

@@ -292,7 +292,8 @@ int main(int argc, char** argv) {
   InteractionCounter sibyllCounted(sibyll);
   corsika::sibyll::NuclearInteraction sibyllNuc(sibyll, env);
   InteractionCounter sibyllNucCounted(sibyllNuc);
-  auto heModelCounted = make_sequence(sibyllNucCounted, sibyllCounted);
+  auto heModelCounted = make_select([](auto const& p) { return is_nucleus(p.getPID()); },
+                                    sibyllNucCounted, sibyllCounted);
 
   corsika::pythia8::Decay decayPythia;
 

examples/hybrid_MC.cpp

@@ -158,10 +158,10 @@ int main(int argc, char** argv) {
   setup::Stack stack;
   stack.clear();
   unsigned short const A = std::stoi(std::string(argv[1]));
-  unsigned short Z = std::stoi(std::string(argv[2]));
-  const Code beamCode = get_nucleus_code(A, Z);
+  unsigned short const Z = std::stoi(std::string(argv[2]));
+  Code const beamCode = get_nucleus_code(A, Z);
   auto const mass = get_mass(beamCode);
-  const HEPEnergyType E0 = 1_GeV * std::stof(std::string(argv[3]));
+  HEPEnergyType const E0 = 1_GeV * std::stof(std::string(argv[3]));
   double theta = 0.;
   auto const thetaRad = theta / 180. * M_PI;
 

examples/mars.cpp

@@ -313,7 +313,8 @@ int main(int argc, char** argv) {
   InteractionCounter sibyllCounted(sibyll);
   corsika::sibyll::NuclearInteraction sibyllNuc(sibyll, env);
   InteractionCounter sibyllNucCounted(sibyllNuc);
-  auto heModelCounted = make_sequence(sibyllNucCounted, sibyllCounted);
+  auto heModelCounted = make_select([](auto const& p) { return is_nucleus(p.getPID()); },
+                                    sibyllNucCounted, sibyllCounted);
 
   corsika::pythia8::Decay decayPythia;
 

examples/vertical_EAS.cpp

@@ -251,8 +251,10 @@ int main(int argc, char** argv) {
         : cutE_(cutE) {}
     bool operator()(const Particle& p) { return (p.getEnergy() < cutE_); }
   };
-  auto hadronSequence = make_select(EnergySwitch(55_GeV), urqmdCounted,
-                                    make_sequence(sibyllNucCounted, sibyllCounted));
+  auto hadronSequence =
+      make_select(EnergySwitch(55_GeV), urqmdCounted,
+                  make_select([](auto const& p) { return is_nucleus(p.getPID()); },
+                              sibyllNucCounted, sibyllCounted));
   auto decaySequence = make_sequence(decayPythia, decaySibyll);
 
   // directory for outputs

tests/modules/testSibyll.cpp

@@ -60,7 +60,8 @@ TEST_CASE("Sibyll", "modules") {
   }
 
   SECTION("cross-section type") {
-
+    CHECK(corsika::sibyll::getSibyllXSCode(Code::Helium) == 0);
+    CHECK(corsika::sibyll::getSibyllXSCode(Code::Proton) == 1);
     CHECK(corsika::sibyll::getSibyllXSCode(Code::Electron) == 0);
     CHECK(corsika::sibyll::getSibyllXSCode(Code::K0Long) == 3);
     CHECK(corsika::sibyll::getSibyllXSCode(Code::SigmaPlus) == 1);