From e8d2973a79a4bde7512a360fd6f0a11dcee61bbe Mon Sep 17 00:00:00 2001
From: hrocho <hrocho@vodacionline.sk>
Date: Wed, 9 Oct 2019 08:19:07 +0200
Subject: [PATCH] config change

---
 calculations.cpp | 122 ++++++++++----------
 config.h         |  27 +----
 constants.h      |  13 +--
 data_atima.h     | 288 +----------------------------------------------
 4 files changed, 72 insertions(+), 378 deletions(-)

diff --git a/calculations.cpp b/calculations.cpp
index f0875e1..cbd3cd5 100644
--- a/calculations.cpp
+++ b/calculations.cpp
@@ -9,7 +9,7 @@
 #include "catima/generated_LS_coeff.h"
 #include "catima/nucdata.h"
 #include "catima/storage.h"
-#include "catima/srim.h" 
+#include "catima/srim.h"
 #ifdef GLOBAL
 extern "C"
 {
@@ -35,11 +35,11 @@ double dedx_n(const Projectile &p, const Material &mat){
         w = mat.weight_fraction(i);
         sum += w*dedx_n(p,t);
     }
-    return sum;   
+    return sum;
 }
 
 double dedx_n(const Projectile &p, const Target &t){
-    
+
     double zpowers = pow(p.Z,0.23)+pow(t.Z,0.23);
     double asum = p.A + t.A;
     double epsilon = 32.53*t.A*1000*p.T*p.A/(p.Z*t.Z*asum*zpowers); //projectile energy is converted from MeV/u to keV
@@ -69,7 +69,7 @@ double bethek_dedx_e(Projectile &p,const Material &mat, const Config &c){
         w = mat.weight_fraction(i);
         sum += w*bethek_dedx_e(p,t,c,mat.I());
     }
-    return sum;   
+    return sum;
 }
 
 double bethek_dedx_e(Projectile &p, const Target &t, const Config &c, double I){
@@ -80,7 +80,7 @@ double bethek_dedx_e(Projectile &p, const Target &t, const Config &c, double I){
     double gamma=1.0 + p.T/atomic_mass_unit;
     double beta2=1.0-1.0/(gamma*gamma);
     double beta = sqrt(beta2);
-    double zp_eff = z_effective(p,t,c);    
+    double zp_eff = z_effective(p,t,c);
     assert(zp_eff>=0);
     double Ipot = (I>0.0)?I:ipot(t.Z);
     assert(Ipot>0);
@@ -91,21 +91,21 @@ double bethek_dedx_e(Projectile &p, const Target &t, const Config &c, double I){
     if(!(c.corrections&corrections::no_shell_correction) &&  eta>=0.13){ //shell corrections
         double cor = (+0.422377*pow(eta,-2)
                     +0.0304043*pow(eta,-4)
-                    -0.00038106*pow(eta,-6))*1e-6*pow(Ipot,2) 
-                  +(+3.858019*pow(eta,-2) 
+                    -0.00038106*pow(eta,-6))*1e-6*pow(Ipot,2)
+                  +(+3.858019*pow(eta,-2)
                     -0.1667989*(pow(eta,-4))
-                    +0.00157955*(pow(eta,-6)))*1.0e-9*pow(Ipot,3); 
+                    +0.00157955*(pow(eta,-6)))*1.0e-9*pow(Ipot,3);
         f2 = f2 -cor/t.Z;
     }
     f2+=2*log(gamma) -beta2;
-    
+
     double barkas=1.0;
     if(!(c.corrections&corrections::no_barkas)){
         barkas = bethek_barkas(zp_eff,eta,t.Z);
         }
-    
+
     double delta = bethek_density_effect(beta, t.Z);
-    
+
     double LS = 0.0;
     if(!(c.corrections&corrections::no_lindhard)){
         //double LS = bethek_lindhard(p);
@@ -118,7 +118,7 @@ double bethek_dedx_e(Projectile &p, const Target &t, const Config &c, double I){
         result += pair_production(p,t);
         result += bremsstrahlung(p,t);
     }
-    
+
     return result;
 }
 
@@ -149,12 +149,12 @@ double bethek_density_effect(double beta, int zt){
     double x = log(beta * gamma) / 2.3025851;
     int i;
     double del = 0;
-    
+
     if(zt>97){  // check if data are available, if not take highest z data
         zt=97;
         }
     i = zt-1;
-    
+
     if (x < density_effect::x0[i] ){
         if(density_effect::del_0[i] > 0.)del = density_effect::del_0[i] * pow(10.0,(2.*(x-density_effect::x0[i])));
         }
@@ -176,7 +176,7 @@ double bethek_lindhard(const Projectile &p){
     double beta_gamma_R = beta*gamma*rho;
     double sum = 0;
     int n=1;
-    
+
     if(gamma < 10.0/rho){
         double dk[3];
         double dmk = 0;
@@ -199,7 +199,7 @@ double bethek_lindhard(const Projectile &p){
                 std::complex<double> cpiske(0.0,(M_PI*(l-sk)/2.0) - lngamma(csketa).imag());
                 std::complex<double> cedr = cexir*std::exp(cpiske);
                 double H=0;
-                
+
                 //std::complex<double> ceds(0.0,0.0);
                 // finite struct part
                 std::complex<double> cmsketa (-sk + 1.0, eta);
@@ -215,10 +215,10 @@ double bethek_lindhard(const Projectile &p){
                 std::complex<double> clambda_s = cexis*std::exp(cmbeta_gamma_R)*hyperg(cmsketa,cm2sk_1,c2beta_gamma_R);
                 std::complex<double> cGrGs = lngamma(cm2sk_1);
                 double GrGs = clambda_r.imag()/clambda_s.imag();
-                GrGs *= exp( lngamma(csketa).real() 
-                             - lngamma(cmsketa).real() 
+                GrGs *= exp( lngamma(csketa).real()
+                             - lngamma(cmsketa).real()
                              - lngamma(c2sk_1).real()
-                             + cGrGs.real() 
+                             + cGrGs.real()
                              + 2.0*sk*log(2.0*beta_gamma_R));
                 if(cos(cGrGs.imag()) < 1.0)GrGs*=-1;
                 if(fabs(GrGs)>1.0e-9){
@@ -247,16 +247,16 @@ double bethek_lindhard(const Projectile &p){
                     }
                     double figi= (k>0) ? asum/bsum : bsum/asum;
                     H = (FrGr - figi)/(figi-FsGs)* GrGs;
-                    
+
                 }
-                else 
+                else
                     H = 0;
-                
-                
+
+
                 dk[i] = std::arg(cedr + H*ceds);
             }
             if(n>1)dk[2] = dmk;
-        
+
             double sdm2 = sin(dk[2]-dk[1]);
             double term1 = k0*(k0+1.0)*sdm2*sdm2/(eta*eta*(2.0*k0 + 1.0));
             if(n>1){
@@ -270,13 +270,13 @@ double bethek_lindhard(const Projectile &p){
             n += 1;
             dmk = dk[1];
             dkm1 = dk[0];
-            
+
         }// end of while n<100
-        
+
     }
     else{ // ultrarelativistic limit
         sum = -log(beta_gamma_R) - 0.2;
-    }   
+    }
     return sum + (0.5*beta2);
 }
 
@@ -290,14 +290,14 @@ double bethek_lindhard_X(const Projectile &p){
     double beta_gamma_R = beta*gamma*rho;
     double sum = 0;
     int n=1;
-    
+
     if(1){
         double dk[4];
         double dmk = 0;
         double dmkp1 = 0;
         double dkm1 = 0;
         double dkm2 = 0;
-        
+
         while(n<1000){
             double k0 = n;
             //int max = (n==1)?4:2;
@@ -318,7 +318,7 @@ double bethek_lindhard_X(const Projectile &p){
                 std::complex<double> cpiske(0.0,(M_PI*(l-sk)/2.0) - lngamma(csketa).imag());
                 std::complex<double> cedr = cexir*std::exp(cpiske);
                 double H=0;
-                
+
                 //std::complex<double> ceds(0.0,0.0);
                 // finite struct part
                 std::complex<double> cmsketa (-sk + 1.0, eta);
@@ -334,10 +334,10 @@ double bethek_lindhard_X(const Projectile &p){
                 std::complex<double> clambda_s = cexis*std::exp(cmbeta_gamma_R)*hyperg(cmsketa,cm2sk_1,c2beta_gamma_R);
                 std::complex<double> cGrGs = lngamma(cm2sk_1);
                 double GrGs = clambda_r.imag()/clambda_s.imag();
-                GrGs *= exp( lngamma(csketa).real() 
-                             - lngamma(cmsketa).real() 
+                GrGs *= exp( lngamma(csketa).real()
+                             - lngamma(cmsketa).real()
                              - lngamma(c2sk_1).real()
-                             + cGrGs.real() 
+                             + cGrGs.real()
                              + 2.0*sk*log(2.0*beta_gamma_R));
                 if(cos(cGrGs.imag()) < 1.0)GrGs*=-1;
                 if(fabs(GrGs)>1.0e-9){
@@ -366,38 +366,38 @@ double bethek_lindhard_X(const Projectile &p){
                     }
                     double figi= (k>0) ? asum/bsum : bsum/asum;
                     H = (FrGr - figi)/(figi-FsGs)* GrGs;
-                    
+
                 }
-                else 
+                else
                     H = 0;
-                
-                
+
+
                 dk[i] = std::arg(cedr + H*ceds);
             }
             if(n>1)dk[2] = dmk;
-        
+
             double strterm1p = 0;
             double strterm1n = 0;
             double strterm2 = 0;
             double strterm3 = 0;
             double eta2 = eta*eta;
-            
+
             double sdm2 = sin(dk[0]-dkm2);
             if(n>2){
                 strterm1p = sdm2*sdm2*(k0-1)*(k0-2)/((2.0*k0 - 1.0)*(2.0*k0-3.0));
             }
-            
+
             sdm2 = sin(dk[2]-dk[3]);
             strterm1n = sdm2*sdm2*(-k0-1)*(-k0-2)/((-2.0*k0 - 1.0)*(-2.0*k0-3.0));
-            
+
             if(n>1){
                 double sd2 = sin(dk[0]-dmkp1);
                 strterm2 += (k0-1.0)*sd2*sd2/((2.0*k0 - 3.0)*(4.0*k0*k0 - 1.0));
             }
-            
+
             double sdd = sin(dk[0]-dk[1]);
             strterm3 = sdd*sdd*(k0+1.0)*((1/(4.0*k0*k0 -1.0))+(1/(4*(k0+1.0)*(k0+1.0) - 1.0)))/(2.0*k0 + 1.0);
-            
+
             //sum += k0*(strterm1p + strterm1n + (strterm2*2) + strterm3)/eta2;
             sum += k0*(strterm1p + strterm1n + (strterm2*2) + strterm3)/eta2;
             sum += - (2.0/k0);
@@ -407,12 +407,12 @@ double bethek_lindhard_X(const Projectile &p){
             dkm2 = dkm1;
             dkm1 = dk[0];
             dmkp1 = dk[2];
-            
+
         }// end of while n<100
-        
+
     }
     else{ // ultrarelativistic limit
-        
+
     }
     double res = 2*bethek_lindhard(p) - sum - beta2;
     return (res>=0)?res:0.0;
@@ -447,7 +447,7 @@ double bremsstrahlung(const Projectile &p, const Target &t){
 double sezi_dedx_e(const Projectile &p, const Material &mat, const Config &c){
     double w;
     double sum=0.0;
-    bool use95 = config_lowenergy(c) == low_energy_types::srim_95;
+    bool use95 = c.low_energy == low_energy_types::srim_95;
     for(int i=0;i<mat.ncomponents();i++){
         auto t = mat.get_element(i);
         w = mat.weight_fraction(i);
@@ -507,13 +507,13 @@ double radiation_length(int z, double m){
     if(z==29){return 12.86;}
     if(z==32){return 12.25;}
     if(z==50){return 8.82;}
-    if(z==54){return 8.48;} 
+    if(z==54){return 8.48;}
     if(z==74){return 6.76;}
     if(z==78){return 6.54;}
     if(z==79){return 6.46;}
     if(z==82){return 6.37;}
     if(z==92){return 6.00;}
-    
+
     double z2 = z*z;
     double z_13 = 1.0/pow(z,1./3.);
     double z_23 = z_13*z_13;
@@ -530,7 +530,7 @@ double precalculated_lindhard(const Projectile &p){
     int z = (int)p.Z ;
     if(z>LS_MAX_Z)z=LS_MAX_Z;
     if(p.T<ls_coefficients::ls_energy_table(0))T=ls_coefficients::ls_energy_table(0);
-    
+
     double da = (p.A - element_atomic_weight(z))/element_atomic_weight(z);
     z = z-1;
     //catima::Interpolator ls_a(ls_coefficients::ls_energy_points,ls_coefficients::ls_coefficients_a[z],LS_NUM_ENERGY_POINTS,interpolation_t::linear);
@@ -544,7 +544,7 @@ double precalculated_lindhard(const Projectile &p){
     //return ls_a(T)+(dif*da/ls_coefficients::a_rel_increase);
     double dif = v2 - v1;
     return v1+(dif*da/ls_coefficients::a_rel_increase);
-    
+
 }
 
 double precalculated_lindhard_X(const Projectile &p){
@@ -556,12 +556,12 @@ double precalculated_lindhard_X(const Projectile &p){
 		return 1.0;
     double da = (p.A - element_atomic_weight(z))/element_atomic_weight(z);
     z = z-1;
-    
+
     //catima::Interpolator ls_X_a(ls_coefficients::ls_energy_table.values,ls_coefficients::ls_X_coefficients_a[z],LS_NUM_ENERGY_POINTS,interpolation_t::linear);
     //catima::Interpolator ls_X_ahi(ls_coefficients::ls_energy_table.values,ls_coefficients::ls_X_coefficients_ahi[z],LS_NUM_ENERGY_POINTS,interpolation_t::linear);
     double v1 = EnergyTable_interpolate(ls_coefficients::ls_energy_table,T,ls_coefficients::ls_X_coefficients_a[z]);
     double v2 = EnergyTable_interpolate(ls_coefficients::ls_energy_table,T,ls_coefficients::ls_X_coefficients_ahi[z]);
-    
+
     //double dif = ls_X_ahi(T) - ls_X_a(T);
     //return ls_X_a(T)+(dif*da/ls_coefficients::a_rel_increase);
     double dif = v2 - v1;
@@ -576,7 +576,7 @@ double dedx_variance(Projectile &p, Target &t, const Config &c){
     double zp_eff = z_effective(p,t,c);
     double f = domega2dx_constant*pow(zp_eff,2)*t.Z/t.A;
 
-    if(config_omega(c) == omega_types::atima){
+    if( (c.calculation == omega_types::atima) ){
         cor = 24.89 * pow(t.Z,1.2324)/(electron_mass*1e6 * beta2)*
 			log( 2.0*electron_mass*1e6*beta2/(33.05*pow(t.Z,1.6364)));
 	    cor = std::max(cor, 0.0 );
@@ -603,7 +603,7 @@ double z_effective(const Projectile &p,const Target &t, const Config &c){
     else if(c.z_effective == z_eff_type::anthony_landorf){
         return z_eff_Anthony_Landford(p.Z, beta, t.Z);
     }
-    
+
     else if(c.z_effective == z_eff_type::hubert){
         return z_eff_Hubert(p.Z, p.T, t.Z);
     }
@@ -625,13 +625,13 @@ double z_effective(const Projectile &p,const Target &t, const Config &c){
 }
 
 double z_eff_Pierce_Blann(double z, double beta){
-    return z*(1.0-exp(-0.95*fine_structure_inverted*beta/pow(z,2.0/3.0))); 
+    return z*(1.0-exp(-0.95*fine_structure_inverted*beta/pow(z,2.0/3.0)));
 }
 
 double z_eff_Anthony_Landford(double pz, double beta, double tz){
     double B = 1.18-tz*(7.5e-03 - 4.53e-05*tz);
     double A = 1.16-tz*(1.91e-03 - 1.26e-05*tz);
-    return pz*(1.0-(A*exp(-fine_structure_inverted*B*beta/pow(pz,2.0/3.0)))); 
+    return pz*(1.0-(A*exp(-fine_structure_inverted*B*beta/pow(pz,2.0/3.0))));
 }
 
 double z_eff_Hubert(double pz, double E, double tz){
@@ -658,13 +658,13 @@ double z_eff_Hubert(double pz, double E, double tz){
         x3 = 0.314 + 0.01072*lntz;
         x4 = 0.5218 + 0.02521*lntz;
     }
-     
+
     return pz*(1-x1*exp(-x2*catima::power(E,x3)*catima::power(pz,-x4)));
 }
 
 double z_eff_Winger(double pz, double beta, double tz){
 	double c0,c1,c2,c3,c4,c5,c6,c7,c8,c9,c10,c11,c12,c13;
-	double x, lnz, lnzt, a0,a1,a2,a3,a4; 
+	double x, lnz, lnzt, a0,a1,a2,a3,a4;
 
 	c0 = 0.4662;
 	c1 = 0.5491;
@@ -697,7 +697,7 @@ double z_eff_Winger(double pz, double beta, double tz){
 	a3 = -c10 * exp( c11*lnz - c12*lnz*lnz*lnz);
 	a4 = -c13;
 
-    return pz * (1. - exp(a0 +a1*x +a2*x*x +a3*x*x*x +a4*x*x*x*x));	
+    return pz * (1. - exp(a0 +a1*x +a2*x*x +a3*x*x*x +a4*x*x*x*x));
 	}
 
 double z_eff_global(double pz, double E, double tz){
@@ -742,7 +742,7 @@ double z_eff_atima14(double pz, double T, double tz){
         qglobal = z_eff_global(pz,T,tz);
         qglobal = (qglobal - qpb)*c1/catima::power(tz+1,c2)*(1-exp(-c3*T)) + qpb;
     }
-    
+
     emax = 1500.;
     emin = 1000.;
     if(T>emax){
diff --git a/config.h b/config.h
index 5e0256a..6e30fa6 100644
--- a/config.h
+++ b/config.h
@@ -3,9 +3,9 @@
 #define CONFIG
 #include <cstring>
 namespace catima{
-    
-    /** 
-      * \enum z_eff_type 
+
+    /**
+      * \enum z_eff_type
       * enum to select formulat to calculate effective charge of the Projectile
       */
     enum z_eff_type:unsigned char {
@@ -30,7 +30,7 @@ namespace catima{
     };
 
     /**
-      * enum to select which dEdx straggling options
+      * enum to select dEdx straggling options
       */
     enum omega_types:unsigned char{
         atima = 0,
@@ -38,7 +38,7 @@ namespace catima{
     };
 
     /**
-      * enum to select which how low energy part is calculated
+      * enum to select how low energy part is calculated
       */
     enum low_energy_types:unsigned char{
         srim_85 = 0,
@@ -57,25 +57,10 @@ namespace catima{
 
         unsigned char corrections = 0;
         unsigned char calculation = 1;
+        unsigned char low_energy = 0;
     };
 
-    inline void set_config_lowenergy(Config c, low_energy_types lt){
-      c.calculation = c.calculation & (lt<<2);
-    }
 
-    inline unsigned char config_lowenergy(const Config c){
-      return (c.calculation>>2) & 0x7;
-    }
-
-    inline void set_config_omega(Config c, omega_types ot){
-      c.calculation = c.calculation & ot;
-    }
-
-    inline unsigned char config_omega(const Config c){
-      return c.calculation & 0x3;
-    }
-
-    
     extern Config default_config;
 }
 
diff --git a/constants.h b/constants.h
index 96ce7c8..6032f92 100644
--- a/constants.h
+++ b/constants.h
@@ -4,22 +4,16 @@
 #include "catima/build_config.h"
 namespace catima {
 
+// config
 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 = 60; // number of datapoints which can be stored in cache
-constexpr double numeric_epsilon = std::numeric_limits<double>::epsilon();
-
-/// required integration precision (relative units)
-/*
-constexpr double int_eps_range = 0.001;
-constexpr double int_eps_range_str = 0.001;
-constexpr double int_eps_ang_str = 0.01;
-constexpr double int_eps_tof = 0.01;
-*/
+constexpr double numeric_epsilon = 10*std::numeric_limits<double>::epsilon();
 
 constexpr double thin_target_limit = 1 - 1e-3;
+
 #ifdef REACTIONS
 constexpr double emin_reaction = 30.0;
 constexpr bool reactions = true;
@@ -27,6 +21,7 @@ constexpr bool reactions = true;
 constexpr bool reactions = false;
 #endif
 
+// constants
 constexpr double PI = 3.1415926535897932384626433832795;
 constexpr double Avogadro = 6.022140857; // 10^23
 constexpr double electron_mass = 0.510998928;   // MeV/c^2
diff --git a/data_atima.h b/data_atima.h
index 47c5e7c..f302e25 100644
--- a/data_atima.h
+++ b/data_atima.h
@@ -3,294 +3,8 @@
 
 // this is taken from atima
 
-namespace catima{ 
-/*
-constexpr double proton_stopping_coef[92][8] = { // proton in material stopping coefficient
-{  .0091827, .0053496, .69741, .48493, 316.07,1.0143, 9329.3, .0539890}, //H
-{  .11393, .0051984, 1.0822, .39252,  1081.0, 1.0645, 4068.5, .0176990}, //He
-{  .85837, .0050147, 1.6044, .38844,  1337.3, 1.0470, 2659.2, .01898},
-{  .8781,  .0051049, 5.4232, .2032,   1200.6, 1.0211, 1401.8, .0385290},
-{ 1.4608,  .0048836, 2.338,  .44249, 1801.3, 1.0352, 1784.1, .02024},
-{ 3.2579,  .0049148, 2.7156, .36473, 2092.2, 1.0291, 2643.6, .0182370}, //C
-{  .59674, .0050837, 4.2073, .30612,  2394.2, 1.0255, 4892.1, .0160060},
-{  .75253, .0050314, 4.0824, .30067, 2455.8, 1.0181, 5069.7, .0174260}, //O
-{ 1.226,   .0051385, 3.2246, .32703, 2525.0, 1.0142, 7563.6, .0194690},
-{ 1.0332,  .0051645, 3.004,  .33889, 2338.6, .99997, 6991.2, .0217990}, //Ne
-{ 6.0972,  .0044292, 3.1929, .45763, 1363.3, .95182, 2380.6, .0818350},
-{14.013,   .0043646, 2.2641, .36326, 2187.4, .99098, 6264.8, .0462},
-{ .039001, .0045415, 5.5463, .39562, 1589.2, .95316, 816.16, .0474840},
-{ 2.072,   .0044516, 3.5585, .53933, 1515.2, .93161, 1790.3, .0351980},
-{17.575,   .0038346, .078694, 1.2388,2806.0, .97284, 1037.6, .0128790},
-{16.126,   .0038315, .054164, 1.3104,2813.3, .96587, 1251.4, .0118470},
-{ 3.217,   .0044579, 3.6696, .5091,  2734.6, .96253, 2187.5, .0169070},
-{ 2.0379,  .0044775, 3.0743, .54773, 3505.0, .97575, 1714.00, .0117010},
-{ .74171, .0043051, 1.1515, .95083,   917.21, .8782, 389.93, .18926},
-{ 9.1316,  .0043809, 5.4611, .31327, 3891.8, .97933, 6267.9, .0151960}, //Ca
-{ 7.2247,  .0043718, 6.1017, .37511, 2829.2, .95218, 6376.1, .0203980},
-{  .147,   .0048456, 6.3485, .41057, 2164.1, .94028, 5292.6, .0502630},
-{ 5.0611,  .0039867, 2.6174, .57957, 2218.9, .92361, 6323.00, .0256690},
-{  .53267, .0042968, .39005, 1.2725, 1872.7, .90776, 64.166, .0301070},
-{  .47697, .0043038, .31452, 1.3289, 1920.5, .90649, 45.576, .0274690},
-{  .027426, .0035443, .031563, 2.1755,  1919.5, .90099, 23.902, .0253630},
-{  .16383, .0043042, .073454, 1.8592, 1918.4, .89678, 27.61, .0231840},
-{ 4.2562,  .0043737, 1.5606, .72067, 1546.8, .87958, 302.02, .0409440}, //Ni
-{ 2.3508,  .0043237, 2.882,  .50113, 1837.7, .89992, 2377.00, .04965},
-{ 3.1095,  .0038455, .11477, 1.5037, 2184.7, .89309, 67.306, .0165880},
-{15.322,   .0040306, .65391, .67668, 3001.7, .92484, 3344.2, .0163660},
-{ 3.6932,  .0044813, 8.608,  .27638,  2982.7, .9276, 3166.6, .0308740},
-{ 7.1373,  .0043134, 9.4247, .27937,  2725.8, .91597, 3166.1, .0250080},
-{ 4.8979,  .0042937, 3.7793, .50004,  2824.5, .91028, 1282.4, .0170610},
-{ 1.3683,  .0043024, 2.5679, .60822,  6907.8, .9817, 628.01, .0068055},
-{ 1.8301,  .0042983, 2.9057, .6038,   4744.6, .94722, 936.64, .0092242},
-{  .42056, .0041169, .01695, 2.3616, 2252.7, .89192, 39.752, .0277570},
-{30.78,    .0037736, .55813, .76816, 7113.2, .97697, 1604.4, .0065268},
-{11.576,   .0042119, 7.0244, .37764, 4713.5, .94264, 2493.2, .01127},
-{ 6.2406,  .0041916, 5.2701, .49453, 4234.6, .93232, 2063.9, .0118440},
-{  .33073, .0041243, 1.7246, 1.1062, 1930.2, .86907, 27.416, .0382080},
-{  .017747, .0041715, .14586, 1.7305,1803.6, .86315, 29.669, .0321230},
-{ 3.7229,  .0041768, 4.6286, .56769, 1678.0, .86202, 3094.00, .06244},
-{  .13998, .0041329, .25573, 1.4241, 1919.3, .86326, 72.797, .0322350},
-{  .2859,  .0041386, .31301, 1.3424, 1954.8, .86175, 115.18, .0293420},
-{ .76,    .0042179, 3.386,  .76285,  1867.4, .85805, 69.994, .0364480},
-{ 6.3957,  .0041935, 5.4689, .41378, 1712.6, .85397,18493.00, .0564710},
-{ 3.4717,  .0041344, 3.2337, .63788, 1116.4, .81959, 4766.0, .1179},
-{ 2.5265,  .0042282, 4.532,  .53562, 1030.8, .81652,16252.0, .19722},
-{ 7.3683,  .0041007, 4.6791, .51428, 1160.0, .82454,17965.0, .13316},
-{ 7.7197,  .004388, 3.242,   .68434, 1428.1, .83398, 1786.7, .0665120},
-{16.78,    .0041918, 9.3198, .29568, 3370.9, .90289, 7431.7, .02616},
-{ 4.2132,  .0042098, 4.6753, .57945, 3503.9, .89261, 1468.9, .0143590},
-{ 4.0818,  .004214, 4.4425,  .58393, 3945.3, .90281, 1340.5, .0134140}, //Xe
-{.18517, .0036215, .00058788,3.5315, 2931.3, .88936, 26.18, .0263930},
-{ 4.8248,  .0041458, 6.0934, .57026, 2300.1, .86359, 2980.7, .0386790},
-{  .49857, .0041054, 1.9775, .95877, 786.55, .78509, 806.6, .40882},
-{ 3.2754,  .0042177, 5.768,  .54054, 6631.3, .94282, 744.07, .0083026},
-{ 2.9978,  .0040901, 4.5299, .62025, 2161.2, .85669, 1268.6, .0430310},
-{ 2.8701,  .004096, 4.2568,  .6138,  2130.4, .85235, 1704.1, .0393850},
-{10.853,   .0041149, 5.8907, .46834, 2857.2, .8755, 3654.2, .0299550},
-{ 3.6407,  .0041782, 4.8742, .57861, 1267.7, .82211, 3508.2, .24174},
-{17.645,   .0040992, 6.5855, .32734, 3931.3, .90754, 5156.7, .0362780},
-{ 7.5309,  .0040814, 4.9389, .50679, 2519.7, .85819, 3314.6, .0305140},
-{ 5.4742,  .0040829, 4.897,  .51113, 2340.1, .85296, 2342.7, .0356620},
-{ 4.2661,  .0040667, 4.5032, .55257, 2076.4, .84151, 1666.6, .0408010},
-{ 6.8313,  .0040486, 4.3987, .51675, 2003.0, .83437, 1410.4, .03478},
-{ 1.2707,  .0040553, 4.6295, .57428, 1626.3, .81858, 995.68, .0553190},
-{ 5.7561,  .0040491, 4.357,  .52496, 2207.3, .83796, 1579.5, .0271650},
-{14.127,   .0040596, 5.8304, .37755, 3645.9, .87823, 3411.8, .0163920},
-{ 6.6948,  .0040603, 4.9361, .47961, 2719.0, .85249, 1885.8, .0197130},
-{ 3.0619,  .0040511, 3.5803, .59082, 2346.1, .83713, 1222.0, .0200720},
-{10.811,   .0033008, 1.3767, .76512, 2003.7, .82269, 1110.6, .0249580},
-{ 2.7101,  .0040961, 1.2289, .98598, 1232.4, .79066, 155.42, .0472940},
-{  .52345, .0040244, 1.4038, .8551,  1461.4, .79677, 503.34, .0367890},
-{  .4616,  .0040203, 1.3014, .87043, 1473.5, .79687, 443.09, .0363010},
-{  .97814, .0040374, 2.0127, .7225,  1890.8, .81747, 930.7, .02769},
-{ 3.2086,  .004051, 3.6658,  .53618, 3091.2, .85602, 1508.1, .0154010},
-{ 2.0035,  .0040431, 7.4882, .3561,  4464.3, .88836, 3966.5, .0128390},
-{15.43,    .0039432, 1.1237, .70703, 4595.7, .88437, 1576.5, .0088534},
-{ 3.1512,  .0040524, 4.0996, .5425,  3246.3, .85772, 1691.8, .0150580},
-{ 7.1896,  .0040588, 8.6927, .35842, 4760.6, .88833, 2888.3, .0110290}, //Pb
-{ 9.3209,  .004054, 11.543,  .32027, 4866.2, .89124, 3213.4, .0119350},
-{29.242,   .0036195, .16864, 1.1226, 5688.0, .89812, 1033.3, .0071303},
-{ 1.8522,  .0039973, 3.1556, .65096, 3755.0, .86383, 1602.0, .0120420},
-{ 3.222,   .0040041, 5.9024, .52678, 4040.2, .86804, 1658.4, .0117470},
-{ 9.3412,  .0039661, 7.921,  .42977, 5180.9, .88773, 2173.2, .0092007},
-{36.183,   .0036003, .58341, .86747, 6990.2, .91082, 1417.1, .0062187},
-{ 5.9284,  .0039695, 6.4082, .52122, 4619.5, .88083, 2323.5, .0116270},
-{ 5.2454,  .0039744, 6.7969, .48542, 4586.3, .87794, 2481.5, .0112820},
-{33.702,   .0036901, .47257, .89235, 5295.7, .8893, 2053.3, .0091908},
-{2.7589,  .0039806, 3.2092, .66122,  2505.4, .82863, 2065.1, .0228160} //U
-};
+namespace catima{
 
-const double  atima_vfermi[92] = {
-1.0309, 
-0.15976,
-0.59782,
-1.0781, 
-1.0486, 
-1.00,   
-1.058,  
-0.93942,
-0.74562,
-0.3424, 
-0.45259,
-0.71074,
-0.90519,
-0.97411,
-0.97184,
-0.89852,
-0.70827,
-0.39816,
-0.36552,
-0.62712,
-0.81707,
-0.9943, 
-1.1423, 
-1.2381, 
-1.1222, 
-0.92705,
-1.0047, 
-1.2,    
-1.0661, 
-0.97411,
-0.84912,
-0.95,   
-1.0903, 
-1.0429, 
-0.49715,
-0.37755,
-0.35211,
-0.57801,
-0.77773,
-1.0207, 
-1.029,  
-1.2542, 
-1.122,  
-1.1241, 
-1.0882, 
-1.2709, 
-1.2542, 
-0.90094,
-0.74093,
-0.86054,
-0.93155,
-1.0047, 
-0.55379,
-0.43289,
-0.32636,
-0.5131, 
-0.6950, 
-0.72591,
-0.71202,
-0.67413,
-0.71418,
-0.71453,
-0.5911, 
-0.70263,
-0.68049,
-0.68203,
-0.68121,
-0.68532,
-0.68715,
-0.61884,
-0.71801,
-0.83048,
-1.1222, 
-1.2381, 
-1.045,  
-1.0733, 
-1.0953, 
-1.2381, 
-1.2879, 
-0.78654,
-0.66401,
-0.84912,
-0.88433,
-0.80746,
-0.43357,
-0.41923,
-0.43638,
-0.51464,
-0.73087,
-0.81065,
-1.9578, 
-1.0257};
-
-
-
-const double atima_lambda_screening[92]= {
- 1.00,
- 1.00,
- 1.10,
- 1.06,
- 1.01,
- 1.03,
- 1.04,
- 0.99,
- 0.95,
- 0.90,
- 0.82,
- 0.81,
- 0.83,
- 0.88,
- 1.00,
- 0.95,
- 0.97,
- 0.99,
- 0.98,
- 0.97,
- 0.98,
- 0.97,
- 0.96,
- 0.93,
- 0.91,
- 0.9, 
- 0.88,
- 0.9, 
- 0.9, 
- 0.9, 
- 0.9, 
- 0.85,
- 0.9, 
- 0.9, 
- 0.91,
- 0.92,
- 0.9, 
- 0.9, 
- 0.9, 
- 0.9, 
- 0.9, 
- 0.88,
- 0.9, 
- 0.88,
- 0.88,
- 0.9, 
- 0.9, 
- 0.88,
- 0.9, 
- 0.9, 
- 0.9, 
- 0.9, 
- 0.96,
- 1.2, 
- 0.9, 
- 0.88,
- 0.88,
- 0.85,
- 0.90,
- 0.90,
- 0.92,
- 0.95,
- 0.99,
- 1.03,
- 1.05,
- 1.07,
- 1.08,
- 1.10,
- 1.08,
- 1.08,
- 1.08,
- 1.08,
- 1.09,
- 1.09,
- 1.10,
- 1.11,
- 1.12,
- 1.13,
- 1.14,
- 1.15,
- 1.17,
- 1.20,
- 1.18,
- 1.17,
- 1.17,
- 1.16,
- 1.16,
- 1.16,
- 1.16,
- 1.16,
- 1.16,
- 1.16};
-
-*/
 }
 namespace density_effect{
 const double x0[97]= {