diff --git a/Documentation/Examples/cascade_example.cc b/Documentation/Examples/cascade_example.cc
index 232e7d002945f41a613266c0738015a2dd34717e..db52594fc2a13c3df511711847264aa318b2a532 100644
--- a/Documentation/Examples/cascade_example.cc
+++ b/Documentation/Examples/cascade_example.cc
@@ -226,7 +226,7 @@ int main() {
   // setup particle stack, and add primary particle
   setup::Stack stack;
   stack.Clear();
-  const hep::EnergyType E0 = 1_TeV;
+  const hep::EnergyType E0 = 100_GeV;
   {
     auto particle = stack.NewParticle();
     particle.SetPID(Code::Proton);
diff --git a/Processes/Sibyll/Interaction.h b/Processes/Sibyll/Interaction.h
index eb45a570dd45dd679ff99437948a484aa7663c75..8fd2795867f1afcb0c04c2333a7ace0cebe8cdca 100644
--- a/Processes/Sibyll/Interaction.h
+++ b/Processes/Sibyll/Interaction.h
@@ -64,8 +64,8 @@ namespace corsika::process::sibyll {
       // FOR NOW: assume target is oxygen
       const int kTarget = corsika::particles::Oxygen::GetNucleusA();
 
-      hep::EnergyType Etot =
-          p.GetEnergy() + kTarget * corsika::particles::Proton::GetMass();
+      const hep::MassType nucleon_mass =  0.5*(corsika::particles::Proton::GetMass()+corsika::particles::Neutron::GetMass());
+      hep::EnergyType Etot = p.GetEnergy() + nucleon_mass;
       MomentumVector Ptot(rootCS, {0.0_GeV, 0.0_GeV, 0.0_GeV});
       // FOR NOW: assume target is at rest
       MomentumVector pTarget(rootCS, {0.0_GeV, 0.0_GeV, 0.0_GeV});
@@ -146,7 +146,8 @@ namespace corsika::process::sibyll {
             GetNucleusA(); // env.GetTargetParticle().GetPID();
 
         // FOR NOW: target is always at rest
-        const EnergyType Etarget = 0_GeV + corsika::particles::Proton::GetMass();
+	const hep::MassType nucleon_mass =  0.5*(corsika::particles::Proton::GetMass()+corsika::particles::Neutron::GetMass());
+        const EnergyType Etarget = 0_GeV + nucleon_mass;
         const auto pTarget = MomentumVector(rootCS, 0_GeV, 0_GeV, 0_GeV);
         cout << "target momentum (GeV/c): " << pTarget.GetComponents() / 1_GeV << endl;
         cout << "beam momentum (GeV/c): " << p.GetMomentum().GetComponents() / 1_GeV