source: sasmodels/sasmodels/models/cylinder.c @ e9b1663d

core_shell_microgelscostrafo411magnetic_modelrelease_v0.94release_v0.95ticket-1257-vesicle-productticket_1156ticket_1265_superballticket_822_more_unit_tests
Last change on this file since e9b1663d was e9b1663d, checked in by Paul Kienzle <pkienzle@…>, 8 years ago

pull from master

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