gwm17/catima
1
0
Fork 0
mirror of https://github.com/gwm17/catima.git synced 2024-11-22 18:28:51 -05:00
Code Issues Projects Releases Wiki Activity

Merge pull request #2 from hrosiak/ls2

Browse Source 
interpolation limits fixed for LS
This commit is contained in:
Andrej Prochazka 2017-07-26 19:34:22 +02:00 committed by GitHub
parent 3afadfc15b fd334455bb
commit d53b18e31c
3 changed files with 18 additions and 4 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

4
calculations.cpp
View File

@ -588,15 +588,15 @@ double radiation_length(int z, int m){
return lr; return lr;
} }
double precalculated_lindhard(const Projectile &p){ double precalculated_lindhard(const Projectile &p){
double T = p.T; double T = p.T;
int z = (int)p.Z ; int z = (int)p.Z ;
if(z>LS_MAX_Z)z=LS_MAX_Z; if(z>LS_MAX_Z)z=LS_MAX_Z;
//if(p.T<ls_coefficients::ls_energy_points[0])T=ls_coefficients::ls_energy_points[0];
if(p.T<ls_coefficients::ls_energy_table(0))T=ls_coefficients::ls_energy_table(0); 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); double da = (p.A - element_atomic_weight(z))/element_atomic_weight(z);
z = z-1; 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); //catima::Interpolator ls_a(ls_coefficients::ls_energy_points,ls_coefficients::ls_coefficients_a[z],LS_NUM_ENERGY_POINTS,interpolation_t::linear);
//catima::Interpolator ls_ahi(ls_coefficients::ls_energy_points,ls_coefficients::ls_coefficients_ahi[z],LS_NUM_ENERGY_POINTS,interpolation_t::linear); //catima::Interpolator ls_ahi(ls_coefficients::ls_energy_points,ls_coefficients::ls_coefficients_ahi[z],LS_NUM_ENERGY_POINTS,interpolation_t::linear);
//catima::Interpolator ls_a(ls_coefficients::ls_energy_table.values,ls_coefficients::ls_coefficients_a[z],LS_NUM_ENERGY_POINTS,interpolation_t::cspline); //catima::Interpolator ls_a(ls_coefficients::ls_energy_table.values,ls_coefficients::ls_coefficients_a[z],LS_NUM_ENERGY_POINTS,interpolation_t::cspline);

14
calculations.h
View File

@ -65,17 +65,31 @@ namespace catima{
/** /**
* returns linhard correction (L) calulated from tabulated precalculated data * returns linhard correction (L) calulated from tabulated precalculated data
* if energy is less than minimal calculated energy the LS coefficient of at minimal
* calculated energy is returned and similar for highest caclulated energy limit
*/ */
double precalculated_lindhard(const Projectile &p); double precalculated_lindhard(const Projectile &p);
/** /**
* returns linhard energy loss straggling correction (X) calulated from tabulated precalculated data * returns linhard energy loss straggling correction (X) calulated from tabulated precalculated data
* if energy is less than minimal calculated energy the X coefficient of at minimal
* calculated energy is returned and similar for highest caclulated energy limit
*/ */
double precalculated_lindhard_X(const Projectile &p); double precalculated_lindhard_X(const Projectile &p);
/**
* this function is not used and is not tested
*/
double energy_straggling_firsov(double z1,double energy, double z2, double m2); double energy_straggling_firsov(double z1,double energy, double z2, double m2);
/**
* electronic energy loss for low energy, should be like SRIM
*/
double sezi_dedx_e(const Projectile &p, const Target &t); double sezi_dedx_e(const Projectile &p, const Target &t);
/**
* electronic energy loss of protons for low energy, should be like SRIM
*/
double sezi_p_se(double energy,const Target &t); double sezi_p_se(double energy,const Target &t);

4
storage.h
View File

@ -60,8 +60,8 @@ namespace catima{
double EnergyTable_interpolate(const EnergyTable<N> &table, double xval, double *y){ double EnergyTable_interpolate(const EnergyTable<N> &table, double xval, double *y){
double r; double r;
double lxval = log(xval)/M_LN10; double lxval = log(xval)/M_LN10;
if(lxval<table.values[0] || lxval>table.values[table.num-1])return 0.0; if(xval<table.values[0] || xval>table.values[table.num-1])return 0.0;
if(lxval==table.values[table.num-1])return y[table.num-1]; if(xval==table.values[table.num-1])return y[table.num-1];
int i = (int)(lxval/table.step); int i = (int)(lxval/table.step);
double linstep = table.values[i+1] - table.values[i]; double linstep = table.values[i+1] - table.values[i];
double x = 1.0 - ((xval - table.values[i])/linstep); double x = 1.0 - ((xval - table.values[i])/linstep);