[a807206] | 1 | double form_volume(double radius_polar, double radius_equatorial); |
---|
| 2 | double Iq(double q, double sld, double sld_solvent, double radius_polar, double radius_equatorial); |
---|
| 3 | double Iqxy(double qx, double qy, double sld, double sld_solvent, |
---|
| 4 | double radius_polar, double radius_equatorial, double theta, double phi); |
---|
[ce27e21] | 5 | |
---|
[a807206] | 6 | double _ellipsoid_kernel(double q, double radius_polar, double radius_equatorial, double sin_alpha); |
---|
| 7 | double _ellipsoid_kernel(double q, double radius_polar, double radius_equatorial, double sin_alpha) |
---|
[ce27e21] | 8 | { |
---|
[a807206] | 9 | double ratio = radius_polar/radius_equatorial; |
---|
| 10 | const double u = q*radius_equatorial*sqrt(1.0 |
---|
[3f832f9] | 11 | + sin_alpha*sin_alpha*(ratio*ratio - 1.0)); |
---|
[9c461c7] | 12 | const double f = sph_j1c(u); |
---|
[513efc5] | 13 | |
---|
[5d4777d] | 14 | return f*f; |
---|
[ce27e21] | 15 | } |
---|
| 16 | |
---|
[a807206] | 17 | double form_volume(double radius_polar, double radius_equatorial) |
---|
[ce27e21] | 18 | { |
---|
[a807206] | 19 | return M_4PI_3*radius_polar*radius_equatorial*radius_equatorial; |
---|
[ce27e21] | 20 | } |
---|
| 21 | |
---|
[994d77f] | 22 | double Iq(double q, |
---|
| 23 | double sld, |
---|
[a807206] | 24 | double sld_solvent, |
---|
| 25 | double radius_polar, |
---|
| 26 | double radius_equatorial) |
---|
[ce27e21] | 27 | { |
---|
[50e1e40] | 28 | // translate a point in [-1,1] to a point in [0, 1] |
---|
| 29 | const double zm = 0.5; |
---|
| 30 | const double zb = 0.5; |
---|
[994d77f] | 31 | double total = 0.0; |
---|
[ce27e21] | 32 | for (int i=0;i<76;i++) { |
---|
[3f832f9] | 33 | //const double sin_alpha = (Gauss76Z[i]*(upper-lower) + upper + lower)/2; |
---|
[50e1e40] | 34 | const double sin_alpha = Gauss76Z[i]*zm + zb; |
---|
[a807206] | 35 | total += Gauss76Wt[i] * _ellipsoid_kernel(q, radius_polar, radius_equatorial, sin_alpha); |
---|
[ce27e21] | 36 | } |
---|
[50e1e40] | 37 | // translate dx in [-1,1] to dx in [lower,upper] |
---|
| 38 | const double form = total*zm; |
---|
[a807206] | 39 | const double s = (sld - sld_solvent) * form_volume(radius_polar, radius_equatorial); |
---|
[50e1e40] | 40 | return 1.0e-4 * s * s * form; |
---|
[ce27e21] | 41 | } |
---|
| 42 | |
---|
[994d77f] | 43 | double Iqxy(double qx, double qy, |
---|
| 44 | double sld, |
---|
[a807206] | 45 | double sld_solvent, |
---|
| 46 | double radius_polar, |
---|
| 47 | double radius_equatorial, |
---|
[994d77f] | 48 | double theta, |
---|
| 49 | double phi) |
---|
[ce27e21] | 50 | { |
---|
[994d77f] | 51 | double sn, cn; |
---|
[ce27e21] | 52 | |
---|
[994d77f] | 53 | const double q = sqrt(qx*qx + qy*qy); |
---|
[ce27e21] | 54 | SINCOS(theta*M_PI_180, sn, cn); |
---|
[3f832f9] | 55 | // TODO: check if this is actually sin(alpha), not cos(alpha) |
---|
[994d77f] | 56 | const double cos_alpha = cn*cos(phi*M_PI_180)*(qx/q) + sn*(qy/q); |
---|
[a807206] | 57 | const double form = _ellipsoid_kernel(q, radius_polar, radius_equatorial, cos_alpha); |
---|
| 58 | const double s = (sld - sld_solvent) * form_volume(radius_polar, radius_equatorial); |
---|
[ce27e21] | 59 | |
---|
[994d77f] | 60 | return 1.0e-4 * form * s * s; |
---|
[ce27e21] | 61 | } |
---|
| 62 | |
---|