From 251f1b8c1ef961142fd7ba61ab12d0fd37264ee8 Mon Sep 17 00:00:00 2001
From: hrocho <hrocho@vodacionline.sk>
Date: Thu, 1 Nov 2018 01:49:48 +0100
Subject: [PATCH] h3_u3

---
 calculations.cpp       |  6 +++---
 catima.pyx             |  1 +
 constants.h            |  2 +-
 storage.h              | 10 +++++++---
 tests/test_storage.cpp |  4 +++-
 5 files changed, 15 insertions(+), 8 deletions(-)

diff --git a/calculations.cpp b/calculations.cpp
index 7dd632d..20faa16 100644
--- a/calculations.cpp
+++ b/calculations.cpp
@@ -113,7 +113,7 @@ double bethek_dedx_e(Projectile &p, const Target &t, const Config &c, double I){
     double result  = (f2)*barkas + LS - delta/2.;
     result *=f1;
 
-    if( (p.T>50000.0) && (!c.dedx&corrections::no_highenergy)){
+    if( (p.T>50000.0) && !(c.dedx&corrections::no_highenergy)){
         result += pair_production(p,t);
         result += bremsstrahlung(p,t);
     }
@@ -413,8 +413,8 @@ double bethek_lindhard_X(const Projectile &p){
     else{ // ultrarelativistic limit
         
     }
-    return 2*bethek_lindhard(p) - sum - beta2;
-    //return sum;
+    double res = 2*bethek_lindhard(p) - sum - beta2;
+    return (res>=0)?res:0.0;
 }
 
 double pair_production(const Projectile &p, const Target &t){
diff --git a/catima.pyx b/catima.pyx
index 2117054..962828c 100644
--- a/catima.pyx
+++ b/catima.pyx
@@ -297,6 +297,7 @@ class corrections(IntEnum):
         no_barkas = 1
         no_lindhard = 2
         no_shell_correction = 4
+        no_highenergy = 8
 
 cdef class Config:
     cdef catimac.Config cbase
diff --git a/constants.h b/constants.h
index 0b3ebfb..24c365d 100644
--- a/constants.h
+++ b/constants.h
@@ -8,7 +8,7 @@ constexpr double Ezero = 1E-3; // lowest E to calculate, below taken as 0
 constexpr double logEmin = -3; // log of minimum energy
 constexpr double logEmax = 7.0;  // log of max energy
 constexpr int max_datapoints = 500; // how many datapoints between logEmin and logEmax
-constexpr int max_storage_data = 100; // number of datapoints which can be stored in cache
+constexpr int max_storage_data = 60; // number of datapoints which can be stored in cache
 constexpr double numeric_epsilon = std::numeric_limits<double>::epsilon();
 
 /// required integration precision (relative units)
diff --git a/storage.h b/storage.h
index c34cda6..223a6ea 100644
--- a/storage.h
+++ b/storage.h
@@ -53,8 +53,10 @@ namespace catima{
     int index(double v)const noexcept{
         double lxval = (log(v/values[0])/M_LN10);
         if(v<values[0] || step==0.0)return -1;
-        if(v>=values[N-1])return N-1;
-        return static_cast<int> (std::floor(lxval/step));
+        if(v>=values[N-1]-numeric_epsilon)return N-1;
+        int i = static_cast<int> (std::floor(lxval/step));
+        if(v >= values[i+1]-numeric_epsilon)i++; // this is correction for floating point precision
+        return i;
     };
 	std::size_t num;
     };
@@ -65,7 +67,9 @@ namespace catima{
 	int EnergyTable_index(const EnergyTable<N> &table, double val){
 		if(val<table.values[0] || val>table.values[table.num-1])return -1;
 		double lxval = (log(val/table.values[0])/M_LN10);
-        return static_cast<int>( std::floor(lxval/table.step));
+        int i = static_cast<int>( std::floor(lxval/table.step));
+        if(val >= table.values[i+1]-numeric_epsilon)i++; // this is correction for floating point precision
+        return i;
 	}
 
 	template<int N>
diff --git a/tests/test_storage.cpp b/tests/test_storage.cpp
index 6415121..12e0d97 100644
--- a/tests/test_storage.cpp
+++ b/tests/test_storage.cpp
@@ -131,11 +131,13 @@ const lest::test specification[] =
       
       EXPECT(EnergyTable_index(catima::energy_table, 0.0)==-1);
 
-      for(int i:{5,10,100,498}){
+      for(int i=0;i<catima::max_datapoints-1;i++){
           val = catima::energy_table.values[i];
           dif = catima::energy_table.values[i+1] - val;
           EXPECT(EnergyTable_index(catima::energy_table, val)==i);
           EXPECT(EnergyTable_index(catima::energy_table, val+0.5*dif)==i);
+          EXPECT(catima::energy_table.index(val)==i);
+          EXPECT(catima::energy_table.index(val+0.5*dif)==i);
       }
     }
 };