Changeset 11ca2ab in sasmodels for sasmodels/models/barbell.c


Ignore:
Timestamp:
Oct 14, 2016 3:13:55 AM (8 years ago)
Author:
Paul Kienzle <pkienzle@…>
Branches:
master, core_shell_microgels, costrafo411, magnetic_model, ticket-1257-vesicle-product, ticket_1156, ticket_1265_superball, ticket_822_more_unit_tests
Children:
5abbad7
Parents:
a80e64c
git-author:
Paul Kienzle <pkienzle@…> (10/14/16 01:39:19)
git-committer:
Paul Kienzle <pkienzle@…> (10/14/16 03:13:55)
Message:

add macros for theta-phi-psi orientation

File:
1 edited

Legend:

Unmodified
Added
Removed

  • sasmodels/models/barbell.c

    r0d6e865 r11ca2ab  
    3939    // translate dx in [-1,1] to dx in [lower,upper] 
    4040    const double integral = total*zm; 
    41     const double bell_Fq = 2*M_PI*cube(radius_bell)*integral; 
    42     return bell_Fq; 
     41    const double bell_fq = 2*M_PI*cube(radius_bell)*integral; 
     42    return bell_fq; 
    4343} 
     44 
     45static double 
     46_fq(double q, double h, 
     47    double radius_bell, double radius, double half_length, 
     48    double sin_alpha, double cos_alpha) 
     49{ 
     50    const double bell_fq = _bell_kernel(q, h, radius_bell, half_length, sin_alpha, cos_alpha); 
     51    const double bj = sas_J1c(q*radius*sin_alpha); 
     52    const double si = sinc(q*half_length*cos_alpha); 
     53    const double cyl_fq = 2.0*M_PI*radius*radius*half_length*bj*si; 
     54    const double Aq = bell_fq + cyl_fq; 
     55    return Aq; 
     56} 
     57 
    4458 
    4559double form_volume(double radius_bell, 
    … …  
    4761        double length) 
    4862{ 
    49  
    5063    // bell radius should never be less than radius when this is called 
    5164    const double hdist = sqrt(square(radius_bell) - square(radius)); 
    … …  
    7184        double sin_alpha, cos_alpha; // slots to hold sincos function output 
    7285        SINCOS(alpha, sin_alpha, cos_alpha); 
    73  
    74         const double bell_Fq = _bell_kernel(q, h, radius_bell, half_length, sin_alpha, cos_alpha); 
    75         const double bj = sas_J1c(q*radius*sin_alpha); 
    76         const double si = sinc(q*half_length*cos_alpha); 
    77         const double cyl_Fq = M_PI*radius*radius*length*bj*si; 
    78         const double Aq = bell_Fq + cyl_Fq; 
     86        const double Aq = _fq(q, h, radius_bell, radius, half_length, sin_alpha, cos_alpha); 
    7987        total += Gauss76Wt[i] * Aq * Aq * sin_alpha; 
    8088    } 
    … …  
    93101        double theta, double phi) 
    94102{ 
    95     // Compute angle alpha between q and the cylinder axis 
    96     double sn, cn; // slots to hold sincos function output 
    97     SINCOS(phi*M_PI_180, sn, cn); 
    98     const double q = sqrt(qx*qx + qy*qy); 
    99     const double cos_val = (q==0. ? 1.0 : (cn*qx + sn*qy)*sin(theta*M_PI_180)/q); 
    100     const double alpha = acos(cos_val); // rod angle relative to q 
     103    double q, sin_alpha, cos_alpha; 
     104    ORIENT_SYMMETRIC(qx, qy, theta, phi, q, sin_alpha, cos_alpha); 
    101105 
    102106    const double h = -sqrt(square(radius_bell) - square(radius)); 
    103     const double half_length = 0.5*length; 
    104  
    105     double sin_alpha, cos_alpha; // slots to hold sincos function output 
    106     SINCOS(alpha, sin_alpha, cos_alpha); 
    107     const double bell_Fq = _bell_kernel(q, h, radius_bell, half_length, sin_alpha, cos_alpha); 
    108     const double bj = sas_J1c(q*radius*sin_alpha); 
    109     const double si = sinc(q*half_length*cos_alpha); 
    110     const double cyl_Fq = M_PI*radius*radius*length*bj*si; 
    111     const double Aq = cyl_Fq + bell_Fq; 
     107    const double Aq = _fq(q, h, radius_bell, radius, 0.5*length, sin_alpha, cos_alpha); 
    112108 
    113109    // Multiply by contrast^2 and convert to cm-1 
Note: See TracChangeset for help on using the changeset viewer.