diff --git a/calculations.cpp b/calculations.cpp index 955ac4d..c32a9a3 100644 --- a/calculations.cpp +++ b/calculations.cpp @@ -486,6 +486,15 @@ double angular_scattering_variance(const Projectile &p, const Target &t){ return 198.81 * pow(p.Z,2)/(lr*pow(_p*beta,2)); } +double angular_scattering_variance(const Projectile &p, const Material &mat){ + if(p.T<=0)return 0.0; + double e=p.T; + double _p = p_from_T(e,p.A); + double beta = _p/((e+atomic_mass_unit)*p.A); + double X0 = radiation_length(mat); + return 198.81 * pow(p.Z,2)/(X0*pow(_p*beta,2)); +} + /// radiation lengths are taken from Particle Data Group 2014 double radiation_length(int z, double m){ double lr = 0; @@ -524,6 +533,16 @@ double radiation_length(int z, double m){ return lr; } +double radiation_length(const Material &material){ + double sum = 0; + for(int i=0;idouble{ - double range = range_spline(T); double e =energy_out(T,x*t.density(),range_spline); return (rrange-x)*(rrange-x)*da2dx(p(e), t, c)*t.density(); }; diff --git a/structures.h b/structures.h index 057fb7c..4a6a811 100644 --- a/structures.h +++ b/structures.h @@ -73,6 +73,7 @@ namespace catima{ double th=0; double molar_mass=0; double i_potential=0; + double _X0=0; std::vectoratoms; public: @@ -187,6 +188,17 @@ namespace catima{ */ double I() const {return i_potential;}; + /** + * set radiation length in g/cm2 unit + * @param value - radiation length in g/cm2 unit, if 0 default radiation length will be calculated + */ + Material& X0(double value){_X0 = value;return *this;} + + /** + * return radiation length in g/cm2 unit + */ + double X0(){return _X0;} + /** * return number density of atoms/molecules per cm3 in 10^23 units */ @@ -215,7 +227,7 @@ namespace catima{ double number_density_cm2(int i)const{ if(i>=atoms.size())return 0.0; return number_density_cm2()*molar_fraction(i); - } + } friend bool operator==(const Material &a, const Material&b); }; diff --git a/tests/test_calculations.cpp b/tests/test_calculations.cpp index 6d461ba..9a00099 100644 --- a/tests/test_calculations.cpp +++ b/tests/test_calculations.cpp @@ -404,6 +404,7 @@ using namespace std; auto water = catima::get_material(catima::material::Water); auto res2 = catima::calculate(p(600),water,600); CHECK(res2.dEdxi == approx(92.5).epsilon(2)); + CHECK(catima::radiation_length(water)==approx(36.1,0.2)); } TEST_CASE("z_eff"){ using namespace catima;