From 08dd75ef0fa7b3f97e91d663490e8d87aa3e97c9 Mon Sep 17 00:00:00 2001
From: Felix Riehn <friehn@lip.pt>
Date: Thu, 1 Aug 2024 12:35:32 +0200
Subject: [PATCH] cleanup tests, add test for cross section bug

---
 tests/modules/testPythia8Interaction.inl | 96 ++++++++----------------
 1 file changed, 33 insertions(+), 63 deletions(-)

diff --git a/tests/modules/testPythia8Interaction.inl b/tests/modules/testPythia8Interaction.inl
index 71e35186a..1ccc00fcc 100644
--- a/tests/modules/testPythia8Interaction.inl
+++ b/tests/modules/testPythia8Interaction.inl
@@ -35,57 +35,18 @@ SECTION("pythia interaction configurations") {
   REQUIRE(collision.canInteract(Code::PiMinus));
   REQUIRE(collision.canInteract(Code::PiPlus));
   REQUIRE_FALSE(collision.canInteract(Code::Electron));
+  REQUIRE_FALSE(collision.canInteract(Code::MuPlus));
 }
 
 corsika::units::si::HEPMomentumType P0 = 400_TeV;
 
-SECTION("pythia interaction: proton") {
+SECTION("pythia interaction") {
   // test some combinations of valid target and projectile particles
   // so far only hadron-hadron and hadron-Nucleus is allowed
-  Code const target = GENERATE(Code::Proton, Code::Nitrogen, Code::Oxygen, Code::Argon);
-  Code const projectile = Code::Proton;
-
-  CORSIKA_LOG_INFO("testing: {}-{}", projectile, target);
-  REQUIRE(
-      collision.getCrossSection(
-          projectile, target,
-          {calculate_total_energy(P0, get_mass(projectile)), {rootCS, {0_eV, 0_eV, P0}}},
-          {get_mass(target), {rootCS, {0_eV, 0_eV, 0_eV}}}) > 0_mb);
-
-  collision.doInteraction(view, projectile, target,
-                          {corsika::calculate_total_energy(P0, get_mass(projectile)),
-                           {rootCS, {0_eV, 0_eV, P0}}},
-                          {get_mass(target), {rootCS, {0_eV, 0_eV, 0_eV}}});
-  REQUIRE(view.getSize() >= 2);
-}
-
-SECTION("pythia interaction: PiPlus") {
-  // test some combinations of valid target and projectile particles
-  //
   Code const target = GENERATE(Code::Nitrogen, Code::Oxygen, Code::Argon);
-  Code const projectile = Code::PiPlus;
+  Code const projectile = GENERATE(Code::Proton, Code::PiPlus, Code::KPlus);
 
-  CORSIKA_LOG_INFO("testing: {}-{}", projectile, target);
-  REQUIRE(
-      collision.getCrossSection(
-          projectile, target,
-          {calculate_total_energy(P0, get_mass(projectile)), {rootCS, {0_eV, 0_eV, P0}}},
-          {get_mass(target), {rootCS, {0_eV, 0_eV, 0_eV}}}) > 0_mb);
-
-  collision.doInteraction(view, projectile, target,
-                          {corsika::calculate_total_energy(P0, get_mass(projectile)),
-                           {rootCS, {0_eV, 0_eV, P0}}},
-                          {get_mass(target), {rootCS, {0_eV, 0_eV, 0_eV}}});
-  REQUIRE(view.getSize() >= 2);
-}
-
-SECTION("pythia interaction: KPlus") {
-  // test some combinations of valid target and projectile particles
-  // so far only hadron-hadron and hadron-Nucleus is allowed
-  Code const target = GENERATE(Code::Nitrogen, Code::Oxygen, Code::Argon);
-  Code const projectile = Code::KPlus;
-
-  CORSIKA_LOG_INFO("testing: {}-{}", projectile, target);
+  CORSIKA_LOG_INFO("testing: {} - {}", projectile, target);
   REQUIRE(
       collision.getCrossSection(
           projectile, target,
@@ -100,17 +61,12 @@ SECTION("pythia interaction: KPlus") {
 }
 
 SECTION("pythia interaction: angantyr fail") {
-  // test some combinations of valid target and projectile particles
-  // so far only hadron-hadron and hadron-Nucleus is allowed
+  // when initializing pythia with angantyr, kaon-proton or pion-proton interactions do
+  // not seem to work
   Code const target = Code::Proton;
   Code const projectile = GENERATE(Code::KPlus, Code::PiPlus);
 
   CORSIKA_LOG_INFO("testing: {}-{}", projectile, target);
-  REQUIRE(
-      collision.getCrossSection(
-          projectile, target,
-          {calculate_total_energy(P0, get_mass(projectile)), {rootCS, {0_eV, 0_eV, P0}}},
-          {get_mass(target), {rootCS, {0_eV, 0_eV, 0_eV}}}) > 0_mb);
 
   REQUIRE_THROWS(
       collision.doInteraction(view, projectile, target,
@@ -119,6 +75,22 @@ SECTION("pythia interaction: angantyr fail") {
                               {get_mass(target), {rootCS, {0_eV, 0_eV, 0_eV}}}));
 }
 
+SECTION("pythia cross section bug") {
+  // there is a problem with the interpolation of the energies.
+  // 80GeV lab should work!
+  Code const target = Code::Oxygen;
+  Code const projectile = Code::Proton;
+
+  corsika::units::si::HEPMomentumType P1 = 80_GeV;
+
+  CORSIKA_LOG_INFO("testing: {} - {}", projectile, target);
+  REQUIRE(
+      collision.getCrossSection(
+          projectile, target,
+          {calculate_total_energy(P1, get_mass(projectile)), {rootCS, {0_eV, 0_eV, P1}}},
+          {get_mass(target), {rootCS, {0_eV, 0_eV, 0_eV}}}) > 0_mb);
+}
+
 SECTION("pythia too low energy") {
 
   // this is a projectile neutron with very little energy
@@ -189,25 +161,23 @@ SECTION("pythia wrong projectile") {
   corsika::pythia8::Interaction collision;
   REQUIRE(collision.getCrossSectionInelEla(
               Code::Electron, Code::Electron,
-              {sqrt(static_pow<2>(Electron::mass) + static_pow<2>(100_GeV)),
-               {rootCS, {0_eV, 0_eV, 100_GeV}}},
+              {sqrt(static_pow<2>(Electron::mass) + static_pow<2>(P0)),
+               {rootCS, {0_eV, 0_eV, P0}}},
               {Proton::mass, {rootCS, {0_eV, 0_eV, 0_eV}}}) == std::tuple{0_mb, 0_mb});
 
-  REQUIRE_THROWS(
-      collision.doInteraction(*secViewPtr, Code::Helium, Code::Nitrogen,
-                              {sqrt(static_pow<2>(Helium::mass) + static_pow<2>(100_GeV)),
-                               {rootCS, {0_eV, 0_eV, 100_GeV}}},
-                              {Nitrogen::mass, {rootCS, {0_eV, 0_eV, 0_eV}}}));
+  REQUIRE_THROWS(collision.doInteraction(
+      *secViewPtr, Code::Helium, Code::Nitrogen,
+      {sqrt(static_pow<2>(Helium::mass) + static_pow<2>(P0)), {rootCS, {0_eV, 0_eV, P0}}},
+      {Nitrogen::mass, {rootCS, {0_eV, 0_eV, 0_eV}}}));
 
   // gamma+p not possible
-  REQUIRE(
-      collision.getCrossSection(
-          Code::H0, Code::Proton,
-          {calculate_total_energy(100_GeV, H0::mass), {rootCS, {0_eV, 0_eV, 100_GeV}}},
-          {Proton::mass, {rootCS, {0_eV, 0_eV, 0_eV}}}) == CrossSectionType::zero());
+  REQUIRE(collision.getCrossSection(
+              Code::H0, Code::Proton,
+              {calculate_total_energy(P0, H0::mass), {rootCS, {0_eV, 0_eV, P0}}},
+              {Proton::mass, {rootCS, {0_eV, 0_eV, 0_eV}}}) == CrossSectionType::zero());
 
   REQUIRE(collision.getCrossSectionInelEla(
-              Code::Photon, Code::Proton, {100_GeV, {rootCS, {0_eV, 0_eV, 100_GeV}}},
+              Code::Photon, Code::Proton, {P0, {rootCS, {0_eV, 0_eV, P0}}},
               {Proton::mass, {rootCS, {0_eV, 0_eV, 0_eV}}}) ==
           std::make_tuple(CrossSectionType::zero(), CrossSectionType::zero()));
 }
\ No newline at end of file
-- 
GitLab