Changes in sasmodels/models/capped_cylinder.c [0d6e865:3a48772] in sasmodels
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
sasmodels/models/capped_cylinder.c
r0d6e865 r3a48772 45 45 // translate dx in [-1,1] to dx in [lower,upper] 46 46 const double integral = total*zm; 47 const double cap_Fq = 2 *M_PI*cube(radius_cap)*integral;47 const double cap_Fq = 2.0*M_PI*cube(radius_cap)*integral; 48 48 return cap_Fq; 49 } 50 51 static double 52 _fq(double q, double h, double radius_cap, double radius, double half_length, 53 double sin_alpha, double cos_alpha) 54 { 55 const double cap_Fq = _cap_kernel(q, h, radius_cap, half_length, sin_alpha, cos_alpha); 56 const double bj = sas_J1c(q*radius*sin_alpha); 57 const double si = sinc(q*half_length*cos_alpha); 58 const double cyl_Fq = 2.0*M_PI*radius*radius*half_length*bj*si; 59 const double Aq = cap_Fq + cyl_Fq; 60 return Aq; 49 61 } 50 62 … … 93 105 SINCOS(alpha, sin_alpha, cos_alpha); 94 106 95 const double cap_Fq = _cap_kernel(q, h, radius_cap, half_length, sin_alpha, cos_alpha); 96 const double bj = sas_J1c(q*radius*sin_alpha); 97 const double si = sinc(q*half_length*cos_alpha); 98 const double cyl_Fq = M_PI*radius*radius*length*bj*si; 99 const double Aq = cap_Fq + cyl_Fq; 100 total += Gauss76Wt[i] * Aq * Aq * sin_alpha; // sin_alpha for spherical coord integration 107 const double Aq = _fq(q, h, radius_cap, radius, half_length, sin_alpha, cos_alpha); 108 // sin_alpha for spherical coord integration 109 total += Gauss76Wt[i] * Aq * Aq * sin_alpha; 101 110 } 102 111 // translate dx in [-1,1] to dx in [lower,upper] … … 114 123 double theta, double phi) 115 124 { 116 // Compute angle alpha between q and the cylinder axis 117 double sn, cn; 118 SINCOS(phi*M_PI_180, sn, cn); 119 const double q = sqrt(qx*qx+qy*qy); 120 const double cos_val = (q==0. ? 1.0 : (cn*qx + sn*qy)*sin(theta*M_PI_180)/q); 121 const double alpha = acos(cos_val); // rod angle relative to q 125 double q, sin_alpha, cos_alpha; 126 ORIENT_SYMMETRIC(qx, qy, theta, phi, q, sin_alpha, cos_alpha); 122 127 123 128 const double h = sqrt(radius_cap*radius_cap - radius*radius); 124 const double half_length = 0.5*length; 125 126 double sin_alpha, cos_alpha; // slots to hold sincos function output 127 SINCOS(alpha, sin_alpha, cos_alpha); 128 const double cap_Fq = _cap_kernel(q, h, radius_cap, half_length, sin_alpha, cos_alpha); 129 const double bj = sas_J1c(q*radius*sin_alpha); 130 const double si = sinc(q*half_length*cos_alpha); 131 const double cyl_Fq = M_PI*radius*radius*length*bj*si; 132 const double Aq = cap_Fq + cyl_Fq; 129 const double Aq = _fq(q, h, radius_cap, radius, 0.5*length, sin_alpha, cos_alpha); 133 130 134 131 // Multiply by contrast^2 and convert to cm-1
Note: See TracChangeset
for help on using the changeset viewer.