gwm17/catima
1
0
Fork 0
mirror of https://github.com/gwm17/catima.git synced 2025-04-11 10:28:50 -04:00
Code Issues Projects Releases Wiki Activity

dx2 v1

Browse Source
This commit is contained in:
hrocho 2020-07-30 01:37:20 +02:00
parent b600864bcc
commit 17c77bde5f
2 changed files with 19 additions and 10 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

23
catima.cpp
View File

@ -215,7 +215,7 @@ Result calculate(Projectile &p, const Material &t, const Config &c){
//Interpolator range_straggling_spline(energy_table.values,data.range_straggling.data(),energy_table.num); //Interpolator range_straggling_spline(energy_table.values,data.range_straggling.data(),energy_table.num);
spline_type range_straggling_spline = get_range_straggling_spline(data); spline_type range_straggling_spline = get_range_straggling_spline(data);
spline_type angular_variance_spline = get_angular_variance_spline(data);
if(res.Eout<Ezero){ if(res.Eout<Ezero){
res.dEdxo = 0.0; res.dEdxo = 0.0;
res.sigma_a = 0.0; res.sigma_a = 0.0;
@ -223,25 +223,19 @@ Result calculate(Projectile &p, const Material &t, const Config &c){
res.sigma_E = 0.0; res.sigma_E = 0.0;
} }
else{ else{
res.dEdxo = p.A/range_spline.derivative(res.Eout); res.dEdxo = p.A/range_spline.derivative(res.Eout);
#ifdef THIN_TARGET_APPROXIMATION #ifdef THIN_TARGET_APPROXIMATION
if(thin_target_limit*res.Ein<res.Eout){ if(thin_target_limit*res.Ein<res.Eout){
double edif = (res.Ein-res.Eout); double edif = (res.Ein-res.Eout);
double s1 = range_straggling_spline.derivative(T); double s1 = range_straggling_spline.derivative(T);
double s2 = range_straggling_spline.derivative(res.Eout); double s2 = range_straggling_spline.derivative(res.Eout);
res.sigma_E = res.dEdxo*sqrt(edif*0.5*(s1+s2))/p.A; res.sigma_E = res.dEdxo*sqrt(edif*0.5*(s1+s2))/p.A;
//Interpolator angular_variance_spline(energy_table.values,data.angular_variance.data(),energy_table.num);
spline_type angular_variance_spline = get_angular_variance_spline(data);
s1 = angular_variance_spline.derivative(T); s1 = angular_variance_spline.derivative(T);
s2 = angular_variance_spline.derivative(res.Eout); s2 = angular_variance_spline.derivative(res.Eout);
res.sigma_a = sqrt(0.5*(s1+s2)*edif); res.sigma_a = sqrt(0.5*(s1+s2)*edif);
} }
else{ else{
res.sigma_E = res.dEdxo*sqrt(range_straggling_spline(T) - range_straggling_spline(res.Eout))/p.A; res.sigma_E = res.dEdxo*sqrt(range_straggling_spline(T) - range_straggling_spline(res.Eout))/p.A;
//Interpolator angular_variance_spline(energy_table.values,data.angular_variance.data(),energy_table.num);
spline_type angular_variance_spline = get_angular_variance_spline(data);
res.sigma_a = sqrt(angular_variance_spline(T) - angular_variance_spline(res.Eout)); res.sigma_a = sqrt(angular_variance_spline(T) - angular_variance_spline(res.Eout));
} }
@ -258,8 +252,14 @@ Result calculate(Projectile &p, const Material &t, const Config &c){
res.tof = calculate_tof_from_E(p,res.Eout,t); res.tof = calculate_tof_from_E(p,res.Eout,t);
} }
} }
res.sigma_r = sqrt(range_straggling_spline(T)); res.sigma_r = sqrt(range_straggling_spline(T));
res.Eloss = (res.Ein - res.Eout)*p.A; res.Eloss = (res.Ein - res.Eout)*p.A;
auto fx2 = [&](double x)->double{
double tt = range_spline(T)-range_spline(x);
return (tt-t.thickness())*(tt-t.thickness())*angular_variance_spline.derivative(x);
};
res.sigma_x = sqrt(integrator.integrate(fx2,res.Eout, res.Ein));
#ifdef REACTIONS #ifdef REACTIONS
res.sp = nonreaction_rate(p,t,c); res.sp = nonreaction_rate(p,t,c);
#endif #endif
@ -280,6 +280,7 @@ MultiResult calculate(Projectile &p, const Layers &layers, const Config &c){
res.total_result.sigma_E += r.sigma_E*r.sigma_E; res.total_result.sigma_E += r.sigma_E*r.sigma_E;
res.total_result.tof += r.tof; res.total_result.tof += r.tof;
res.total_result.Eout = r.Eout; res.total_result.Eout = r.Eout;
res.total_result.sigma_x += r.sigma_x*r.sigma_x;
#ifdef REACTIONS #ifdef REACTIONS
res.total_result.sp = (r.sp>=0.0)?res.total_result.sp*r.sp:-1; res.total_result.sp = (r.sp>=0.0)?res.total_result.sp*r.sp:-1;
#endif #endif
@ -288,6 +289,8 @@ MultiResult calculate(Projectile &p, const Layers &layers, const Config &c){
if(e>Ezero){ if(e>Ezero){
res.total_result.sigma_a = sqrt(res.total_result.sigma_a); res.total_result.sigma_a = sqrt(res.total_result.sigma_a);
res.total_result.sigma_E = sqrt(res.total_result.sigma_E); res.total_result.sigma_E = sqrt(res.total_result.sigma_E);
res.total_result.sigma_x = sqrt(res.total_result.sigma_x);
} }
else{ else{
res.total_result.sigma_a = 0.0; res.total_result.sigma_a = 0.0;

6
structures.h
View File

@ -161,6 +161,11 @@ namespace catima{
*/ */
double thickness() const {return th;}; double thickness() const {return th;};
/**
* @return returns thickness in cm
*/
double thickness_cm() const {return th/rho;};
/** /**
* sets thickness in g/cm^2 * sets thickness in g/cm^2
*/ */
@ -233,6 +238,7 @@ namespace catima{
double sigma_E=0.0; double sigma_E=0.0;
double sigma_a=0.0; double sigma_a=0.0;
double sigma_r=0.0; double sigma_r=0.0;
double sigma_x=0.0;
double tof=0.0; double tof=0.0;
#ifdef REACTIONS #ifdef REACTIONS
double sp = 1.0; double sp = 1.0;