source: sasmodels/sasmodels/models/cylinder.c @ 592343f

core_shell_microgelscostrafo411magnetic_modelticket-1257-vesicle-productticket_1156ticket_1265_superballticket_822_more_unit_tests
Last change on this file since 592343f was 592343f, checked in by wojciech, 8 years ago

sas_J1c translated to sas_2J1x_x

  • Property mode set to 100644
File size: 2.1 KB
Line 
1double form_volume(double radius, double length);
2double fq(double q, double sn, double cn,double radius, double length);
3double orient_avg_1D(double q, double radius, double length);
4double Iq(double q, double sld, double solvent_sld, double radius, double length);
5double Iqxy(double qx, double qy, double sld, double solvent_sld,
6    double radius, double length, double theta, double phi);
7
8#define INVALID(v) (v.radius<0 || v.length<0)
9
10double form_volume(double radius, double length)
11{
12    return M_PI*radius*radius*length;
13}
14
15double fq(double q, double sn, double cn, double radius, double length)
16{
17    // precompute qr and qh to save time in the loop
18    const double qr = q*radius;
19    const double qh = q*0.5*length; 
20    return sas_2J1x_x(qr*sn) * sas_sinx_x(qh*cn);
21}
22
23double orient_avg_1D(double q, double radius, double length)
24{
25    // translate a point in [-1,1] to a point in [0, pi/2]
26    const double zm = M_PI_4;
27    const double zb = M_PI_4; 
28
29    double total = 0.0;
30    for (int i=0; i<76 ;i++) {
31        const double alpha = Gauss76Z[i]*zm + zb;
32        double sn, cn; // slots to hold sincos function output
33        // alpha(theta,phi) the projection of the cylinder on the detector plane
34        SINCOS(alpha, sn, cn);
35        total += Gauss76Wt[i] * square( fq(q, sn, cn, radius, length) ) * sn;
36    }
37    // translate dx in [-1,1] to dx in [lower,upper]
38    return total*zm;
39}
40
41double Iq(double q,
42    double sld,
43    double solvent_sld,
44    double radius,
45    double length)
46{
47    const double s = (sld - solvent_sld) * form_volume(radius, length);
48    return 1.0e-4 * s * s * orient_avg_1D(q, radius, length);
49}
50
51
52double Iqxy(double qx, double qy,
53    double sld,
54    double solvent_sld,
55    double radius,
56    double length,
57    double theta,
58    double phi)
59{
60    double q, sin_alpha, cos_alpha;
61    ORIENT_SYMMETRIC(qx, qy, theta, phi, q, sin_alpha, cos_alpha);
62    //printf("sn: %g cn: %g\n", sin_alpha, cos_alpha);
63    const double s = (sld-solvent_sld) * form_volume(radius, length);
64    const double form = fq(q, sin_alpha, cos_alpha, radius, length);
65    return 1.0e-4 * square(s * form);
66}
Note: See TracBrowser for help on using the repository browser.