clean-up sibyll

Processes/Sibyll/Decay.h

@@ -21,11 +21,12 @@ namespace corsika::process {
       // i.e. corsika::particles::ListOfParticles()
       std::cout << "Sibyll: setting hadrons unstable.." << std::endl;
       // make ALL particles unstable, then set EM stable
-      for (auto& p : corsika2sibyll) {
+      for (auto const& p : corsika2sibyll) {
         // std::cout << (int)p << std::endl;
-        const int sibCode = (int)p;
+        const int sibCode = static_cast<int>(p);
         // skip unknown and antiparticles
-        if (sibCode < 1) continue;
+        if (sibCode < 1)
+          continue;
         // std::cout << "Sibyll: Decay: setting " << ConvertFromSibyll(
         // static_cast<SibyllCode> ( sibCode ) ) << " unstable" << std::endl;
         s_csydec_.idb[sibCode - 1] = abs(s_csydec_.idb[sibCode - 1]);
@@ -84,6 +85,7 @@ namespace corsika::process {
     class Decay : public corsika::process::DecayProcess<Decay> {
     public:
       Decay() {}
+      
       void Init() {
         setHadronsUnstable();
         setTrackedParticlesStable();
@@ -97,9 +99,11 @@ namespace corsika::process {
         // i.e. corsika::particles::ListOfParticles()
         for (auto& p : corsika2sibyll) {
           // std::cout << (int)p << std::endl;
-          const int sibCode = (int)p;
+          const int sibCode = static_cast<int>(p);
           // skip unknown and antiparticles
-          if (sibCode < 1) continue;
+          if (sibCode < 1)
+            continue;
+            
           std::cout << "Sibyll: Decay: setting "
                     << ConvertFromSibyll(static_cast<SibyllCode>(sibCode)) << " stable"
                     << std::endl;
@@ -115,23 +119,15 @@ namespace corsika::process {
         corsika::units::hep::EnergyType E = p.GetEnergy();
         corsika::units::hep::MassType m = corsika::particles::GetMass(p.GetPID());
 
-        // const MassDensityType density = 1.25e-3 * kilogram / (1_cm * 1_cm * 1_cm);
-
-        const double gamma = E / m;
-
         const TimeType t0 = particles::GetLifetime(p.GetPID());
+        const double gamma = E / m;
+        corsika::units::si::TimeType const lifetime = gamma * t0;
+        
         cout << "Decay: code: " << (p.GetPID()) << endl;
         cout << "Decay: MinStep: t0: " << t0 << endl;
         cout << "Decay: MinStep: gamma: " << gamma << endl;
-        // cout << "Decay: MinStep: density: " << density << endl;
-        // return as column density
-        // const double x0 = density * t0 * gamma * constants::c / kilogram * 1_cm * 1_cm;
-        // cout << "Decay: MinStep: x0: " << x0 << endl;
-        corsika::units::si::TimeType const lifetime = gamma * t0;
         cout << "Decay: MinStep: tau: " << lifetime << endl;
-        // int a = 1;
-        // const double x = -x0 * log(s_rndm_(a));
-        // cout << "Decay: next decay: " << x << endl;
+        
         return lifetime;
       }
 

Processes/Sibyll/ParticleConversion.h

@@ -34,9 +34,8 @@ namespace corsika::process::sibyll {
   }
 
   SibyllCode constexpr ConvertToSibyll(corsika::particles::Code pCode) {
-    //~ assert(handledBySibyll(pCode));
     return static_cast<SibyllCode>(
-        corsika2sibyll[static_cast<corsika::particles::CodeIntType>(pCode)]);
+      corsika2sibyll[static_cast<corsika::particles::CodeIntType>(pCode)]);
   }
 
   corsika::particles::Code constexpr ConvertFromSibyll(SibyllCode pCode) {
@@ -44,8 +43,7 @@ namespace corsika::process::sibyll {
   }
 
   int constexpr ConvertToSibyllRaw(corsika::particles::Code pCode) {
-    return (int)static_cast<corsika::process::sibyll::SibyllCodeIntType>(
-        corsika::process::sibyll::ConvertToSibyll(pCode));
+    return static_cast<int>(ConvertToSibyll(pCode));
   }
 
   int constexpr GetSibyllXSCode(corsika::particles::Code pCode) {

Processes/Sibyll/sibyll2.3c.cc

 #include <corsika/process/sibyll/sibyll2.3c.h>
 
 #include <corsika/random/RNGManager.h>
+#include <random>
 
 double s_rndm_(int&) {
-  static corsika::random::RNG& rmng =
+  static corsika::random::RNG& rng =
       corsika::random::RNGManager::GetInstance().GetRandomStream("s_rndm");
-  ;
-  return rmng() / (double)rmng.max();
+  
+  std::uniform_real_distribution<double> dist;
+  return dist(rng);
 }

Stack/SuperStupidStack/SuperStupidStack.h

@@ -70,12 +70,6 @@ namespace corsika::stack {
       }
       Point GetPosition() const { return GetStackData().GetPosition(GetIndex()); }
       TimeType GetTime() const { return GetStackData().GetTime(GetIndex()); }
-
-#warning this does not really work, nor makes sense:
-      Vector<SpeedType::dimension_type> GetDirection() const {
-        auto P = GetMomentum();
-        return P / P.norm() * 1e10 * (units::si::meter / units::si::second);
-      }
     };
 
     /**