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

backup2

Browse Source
This commit is contained in:
hrocho 2020-08-05 01:10:03 +02:00
parent be74892e11
commit b38fe06d1f
3 changed files with 26 additions and 11 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
calculations.cpp
View File

@ -486,7 +486,7 @@ double angular_scattering_variance(Projectile &p, Target &t){
return 198.81 * pow(p.Z,2)/(lr*pow(_p*beta,2)); return 198.81 * pow(p.Z,2)/(lr*pow(_p*beta,2));
} }
/// radioation lengths are taken frm Particle Data Group 2014 /// radiation lengths are taken from Particle Data Group 2014
double radiation_length(int z, double m){ double radiation_length(int z, double m){
double lr = 0; double lr = 0;
if(z==1){return 63.04;} if(z==1){return 63.04;}

30
catima.cpp
View File

@ -259,23 +259,35 @@ Result calculate(Projectile &p, const Material &t, const Config &c){
double range = range_spline(T); double range = range_spline(T);
double tt = range - range_spline(x); double tt = range - range_spline(x);
double t0 = std::min(range, t.thickness()); double t0 = std::min(range, t.thickness());
return (tt-t0)*(tt-t0)*angular_variance_spline.derivative(x); return (tt-t0)*(tt-t0)*da2dx(p,x, t, c);
}; };
auto fx2p = [&](double x)->double{
double range = range_spline(T);
double e =energy_out(T,x*t.density(),range_spline);
double t0 = std::min(range/t.density(), t.thickness_cm());
return (t0-x)*(t0-x)*da2dx(p,e, t, c)*t.density();
};
res.sigma_x = integrator_adaptive.integrate(fx2,res.Eout, res.Ein,1e-3, 1e-3,4)/t.density()/t.density(); //res.sigma_x = integrator_adaptive.integrate(fx2,res.Eout, res.Ein,1e-3, 1e-3,4)/t.density();
//res.sigma_x = integrator_adaptive.integrate(fx2p,0, t.thickness_cm(),1e-3, 1e-3,4)/t.density()/t.density(); res.sigma_x = integrator_adaptive.integrate(fx2p,0, t.thickness_cm(),1e-3, 1e-3,4);
res.sigma_x = sqrt(res.sigma_x); res.sigma_x = sqrt(res.sigma_x);
auto fx1 = [&](double x)->double{ auto fx1 = [&](double x)->double{
double range = range_spline(T); double range = range_spline(T);
double tt = range - range_spline(x); double tt = range - range_spline(x);
double t0 = std::min(range, t.thickness()); double t0 = std::min(range, t.thickness());
return (t0-tt)*angular_variance_spline.derivative(x); return (t0-tt)*da2dx(p,x, t, c);
};
auto fx1p = [&](double x)->double{
double range = range_spline(T);
double e =energy_out(T,x*t.density(),range_spline);
double t0 = std::min(range/t.density(), t.thickness_cm());
return (t0-x)*da2dx(p,e, t, c)*t.density();
}; };
res.cov = integrator_adaptive.integrate(fx1,res.Eout, res.Eout, 1e-6, 1e-3,4)/t.density(); //res.cov = integrator_adaptive.integrate(fx1,res.Eout, res.Eout, 1e-6, 1e-3,4);
res.cov = integrator_adaptive.integrate(fx1p,0, t.thickness_cm(), 1e-6, 1e-3,4);
p.T = T;
#ifdef REACTIONS #ifdef REACTIONS
res.sp = nonreaction_rate(p,t,c); res.sp = nonreaction_rate(p,t,c);
#endif #endif
@ -298,7 +310,9 @@ MultiResult calculate(Projectile &p, const Layers &layers, const Config &c){
res.total_result.tof += r.tof; res.total_result.tof += r.tof;
res.total_result.Eout = r.Eout; res.total_result.Eout = r.Eout;
double a2 = res.total_result.sigma_a; double a2 = res.total_result.sigma_a;
res.total_result.sigma_x += (2*m.thickness_cm()*res.total_result.cov) + (a2*m.thickness_cm()*m.thickness_cm()) + r.sigma_x*r.sigma_x; res.total_result.sigma_x += (2*m.thickness_cm()*res.total_result.cov)
+ (a2*m.thickness_cm()*m.thickness_cm())
+ r.sigma_x*r.sigma_x;
//res.total_result.sigma_x += (a2*m.thickness_cm()*m.thickness_cm()) + r.sigma_x*r.sigma_x; //res.total_result.sigma_x += (a2*m.thickness_cm()*m.thickness_cm()) + r.sigma_x*r.sigma_x;
res.total_result.cov += a2*m.thickness_cm() + r.cov; res.total_result.cov += a2*m.thickness_cm() + r.cov;
res.total_result.sigma_a += r.sigma_a*r.sigma_a; res.total_result.sigma_a += r.sigma_a*r.sigma_a;

5
tests/test_calculations.cpp
View File

@ -217,6 +217,7 @@ using namespace std;
{1.00794,1,2}, {1.00794,1,2},
{15.9994,8,1} {15.9994,8,1}
}); });
water2.thickness_cm(9.6);
water2.density(2.0); water2.density(2.0);
water.thickness_cm(9.6); water.thickness_cm(9.6);
@ -250,7 +251,7 @@ using namespace std;
l.add(water); l.add(water);
l.add(water); l.add(water);
l.add(water); l.add(water);
res2 = catima::calculate(p,l); res2 = catima::calculate(p(215),l);
CHECK(res2.total_result.sigma_x == approx(res.sigma_x).R(0.01)); CHECK(res2.total_result.sigma_x == approx(res.sigma_x).R(0.01));
catima::Layers lll; catima::Layers lll;
@ -259,7 +260,7 @@ using namespace std;
lll.add(water); lll.add(water);
lll.add(water); lll.add(water);
lll.add(water); lll.add(water);
res2 = catima::calculate(p,lll); res2 = catima::calculate(p(215),lll);
CHECK(res2.total_result.sigma_x == approx(res29.sigma_x).R(0.01)); CHECK(res2.total_result.sigma_x == approx(res29.sigma_x).R(0.01));
} }