From 3ddc3205fb620dedfbef51a2c9f83856ebac36a3 Mon Sep 17 00:00:00 2001
From: hrocho <hrocho@vodacionline.sk>
Date: Tue, 31 Jul 2018 18:52:53 +0200
Subject: [PATCH] reactions2

---
 reactions.cpp           |  8 ++++----
 tests/test_reaction.cpp | 11 ++++++++++-
 2 files changed, 14 insertions(+), 5 deletions(-)

diff --git a/reactions.cpp b/reactions.cpp
index 2499f8f..74efb1d 100644
--- a/reactions.cpp
+++ b/reactions.cpp
@@ -22,16 +22,17 @@ double nonreaction_rate(Projectile &projectile, const Material &target, const Co
 
     int ap = lround(projectile.A);
     int zp = lround(projectile.Z);
-    int zt = target.get_element(0).Z;
-    int at = abundance::get_isotope_a(zt,0); // most abundand natural isotope mass
 
     auto data = _storage.Get(projectile,target,c);
     Interpolator range_spline(energy_table.values,data.range.data(),energy_table.num);
     if(energy_out(projectile.T, target.thickness(), range_spline) < emin_reaction)return -1.0;
+    
     auto sigma_r = [&](double th){
         double stn_sum=0.0, sum=0.0;
         double e = energy_out(projectile.T, th, range_spline);
         for(unsigned int i = 0;i<target.ncomponents();i++){
+            int zt = target.get_element(i).Z;
+            int at = abundance::get_isotope_a(zt,0); // most abundand natural isotope mass
             stn_sum += target.molar_fraction(i);
             sum += target.molar_fraction(i)*SigmaR_Kox(ap,zp,e,at,zt); 
         }
@@ -50,9 +51,8 @@ double nonreaction_rate(Projectile &projectile, const Material &target, const Co
         cs = target.number_density_cm2()*(cs0 + cs1)/2.0;
     }
     else{
-        cs = catima::integrator.integrate(sigma_r,0,target.number_density_cm2());
+        cs = Avogadro*catima::integrator.integrate(sigma_r,0,target.thickness())/target.M();
     }
-    
     return exp(-cs*0.0001);
     }
 
diff --git a/tests/test_reaction.cpp b/tests/test_reaction.cpp
index 7ea291b..58c13dd 100644
--- a/tests/test_reaction.cpp
+++ b/tests/test_reaction.cpp
@@ -11,11 +11,12 @@ const lest::test specification[] =
 {
     CASE("reaction"){
         catima::Projectile proj{12,6,6,870};
+        catima::Projectile proj2{238,92,92,500};
         auto c = catima::get_material(6);
         auto h = catima::get_material(1);
         catima::Material water({{0,8,2},{0,1,1}},1.0);
         c.thickness(2.0);
-        double r;
+        double r,r2;
 
         r= catima::nonreaction_rate(proj, c);
         EXPECT(r == approx(0.92,0.02));
@@ -40,6 +41,14 @@ const lest::test specification[] =
         r= catima::nonreaction_rate(proj, c);
         EXPECT(r == -1.0);
 
+        proj.T=870;
+        water.thickness(1);
+        r= catima::nonreaction_rate(proj2, water);
+        r2= catima::nonreaction_rate(proj, water);
+        EXPECT( (r > 0 && r<1.0) );
+        EXPECT( (r2 > 0 && r2<1.0) );
+        EXPECT( r2>r );
+
     }
 };