Changeset 997d4eb in sasmodels


Ignore:
Timestamp:
Oct 12, 2016 10:36:00 PM (7 years ago)
Author:
GitHub <noreply@…>
Parents:
1a6cd57 (diff), a80e64c (diff)
Note: this is a merge changeset, the changes displayed below correspond to the merge itself.
Use the (diff) links above to see all the changes relative to each parent.
git-author:
jhbakker <j.h.bakker@…> (10/12/16 22:36:00)
git-committer:
GitHub <noreply@…> (10/12/16 22:36:00)
Message:

Merge a80e64c6fabe2895c22d039441c659df03f77dfa into 1a6cd57bd7b12db89d69b335e984c55413295845

Files:
5 added
24 edited

Legend:

Unmodified
Added
Removed
  • sasmodels/kernel_template.c

    r6a0d6aa r0d6e865  
    264264      #if defined(IQXY_HAS_THETA) 
    265265        // Force a nominal value for the spherical correction even when 
    266         // theta is +90/-90 so that there are no divide by zero problems. 
    267         // For cos(theta) fixed at 90, we effectively multiply top and bottom 
     266        // theta is +0/180 so that there are no divide by zero problems. 
     267        // For sin(theta) fixed at 0 and 180, we effectively multiply top and bottom 
    268268        // by 1e-6, so the effect cancels. 
    269269        const double spherical_correction = fmax(fabs(cos(M_PI_180*theta)), 1.e-6); 
  • sasmodels/kernelpy.py

    r14a15a3 r0d6e865  
    1010 
    1111import numpy as np  # type: ignore 
    12 from numpy import pi, cos  #type: ignore 
     12from numpy import pi, sin, cos  #type: ignore 
    1313 
    1414from . import details 
     
    220220            partial_weight = np.prod(pd_weight[pd_offset+pd_index][1:]) 
    221221            if call_details.theta_par >= 0: 
    222                 cor = cos(pi / 180 * parameters[call_details.theta_par]) 
     222                cor = sin(pi / 180 * parameters[call_details.theta_par]) 
    223223                spherical_correction = max(abs(cor), 1e-6) 
    224224            p0_index = loop_index%p0_length 
  • sasmodels/models/barbell.c

    r2222134 r0d6e865  
    9595    // Compute angle alpha between q and the cylinder axis 
    9696    double sn, cn; // slots to hold sincos function output 
    97     SINCOS(theta*M_PI_180, sn, cn); 
    98     const double q = sqrt(qx*qx+qy*qy); 
    99     const double cos_val = cn*cos(phi*M_PI_180)*(qx/q) + sn*(qy/q); 
     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); 
    100100    const double alpha = acos(cos_val); // rod angle relative to q 
    101101 
  • sasmodels/models/barbell.py

    rb0c4271 r0d6e865  
    7272    Definition of the angles for oriented 2D barbells. 
    7373 
    74 .. figure:: img/cylinder_angle_projection.jpg 
    75     :width: 600px 
    76  
    77     Examples of the angles for oriented pp against the detector plane. 
    7874 
    7975References 
  • sasmodels/models/capped_cylinder.c

    r2222134 r0d6e865  
    116116    // Compute angle alpha between q and the cylinder axis 
    117117    double sn, cn; 
    118     SINCOS(theta*M_PI_180, sn, cn); 
     118    SINCOS(phi*M_PI_180, sn, cn); 
    119119    const double q = sqrt(qx*qx+qy*qy); 
    120     const double cos_val = cn*cos(phi*M_PI_180)*(qx/q) + sn*(qy/q); 
     120    const double cos_val = (q==0. ? 1.0 : (cn*qx + sn*qy)*sin(theta*M_PI_180)/q); 
    121121    const double alpha = acos(cos_val); // rod angle relative to q 
    122122 
  • sasmodels/models/capped_cylinder.py

    rb0c4271 r0d6e865  
    7575    Definition of the angles for oriented 2D cylinders. 
    7676 
    77 .. figure:: img/cylinder_angle_projection.jpg 
    78     :width: 600px 
    79  
    80     Examples of the angles for oriented 2D cylinders against the detector plane. 
    8177 
    8278References 
     
    132128              ["radius_cap", "Ang",     20, [0, inf],    "volume", "Cap radius"], 
    133129              ["length",     "Ang",    400, [0, inf],    "volume", "Cylinder length"], 
    134               ["theta",      "degrees", 60, [-inf, inf], "orientation", "In plane angle"], 
    135               ["phi",        "degrees", 60, [-inf, inf], "orientation", "Out of plane angle"], 
     130              ["theta",      "degrees", 60, [-inf, inf], "orientation", "inclination angle"], 
     131              ["phi",        "degrees", 60, [-inf, inf], "orientation", "deflection angle"], 
    136132             ] 
    137133# pylint: enable=bad-whitespace, line-too-long 
  • sasmodels/models/core_shell_bicelle.c

    r2222134 r0d6e865  
    117117 
    118118    // Cylinder orientation 
    119     const double cyl_x = cos(theta) * cos(phi); 
    120     const double cyl_y = sin(theta); 
     119    const double cyl_x = sin(theta) * cos(phi); 
     120    const double cyl_y = sin(theta) * sin(phi); 
    121121 
    122122    // Compute the angle btw vector q and the axis of the cylinder 
  • sasmodels/models/core_shell_bicelle.py

    r1a6cd57 r1fdb555  
    7171    Definition of the angles for the oriented core shell bicelle tmodel. 
    7272 
    73 .. figure:: img/cylinder_angle_projection.jpg 
    74     :width: 600px 
    75  
    76     Examples of the angles for oriented pp against the detector plane. 
    7773 
    7874References 
     
    161157 
    162158qx, qy = 0.4 * cos(90), 0.5 * sin(0) 
    163 tests = [ 
    164     # Accuracy tests based on content in test/utest_other_models.py 
    165     [{'radius': 20.0, 
    166       'thick_rim': 10.0, 
    167       'thick_face': 10.0, 
    168       'length': 400.0, 
    169       'sld_core': 1.0, 
    170       'sld_face': 4.0, 
    171       'sld_rim': 4.0, 
    172       'sld_solvent': 1.0, 
    173       'background': 0.0, 
    174      }, 0.001, 353.550], 
    175  
    176     [{'radius': 20.0, 
    177       'thick_rim': 10.0, 
    178       'thick_face': 10.0, 
    179       'length': 400.0, 
    180       'sld_core': 1.0, 
    181       'sld_face': 4.0, 
    182       'sld_rim': 4.0, 
    183       'sld_solvent': 1.0, 
    184       'theta': 90.0, 
    185       'phi': 0.0, 
    186       'background': 0.00, 
    187      }, (qx, qy), 24.9167], 
    188  
    189     # Additional tests with larger range of parameters 
    190     [{'radius': 3.0, 
    191       'thick_rim': 100.0, 
    192       'thick_face': 100.0, 
    193       'length': 1200.0, 
    194       'sld_core': 5.0, 
    195       'sld_face': 41.0, 
    196       'sld_rim': 42.0, 
    197       'sld_solvent': 21.0, 
    198      }, 0.05, 1670.1828], 
    199     ] 
  • sasmodels/models/core_shell_cylinder.c

    r43b7eea r0d6e865  
    6969 
    7070    // Compute angle alpha between q and the cylinder axis 
    71     SINCOS(theta*M_PI_180, sn, cn); 
     71    SINCOS(phi*M_PI_180, sn, cn); 
    7272    // # The following correction factor exists in sasview, but it can't be 
    7373    // # right, so we are leaving it out for now. 
    7474    // const double correction = fabs(cn)*M_PI_2; 
    7575    const double q = sqrt(qx*qx+qy*qy); 
    76     const double cos_val = cn*cos(phi*M_PI_180)*(qx/q) + sn*(qy/q); 
     76    const double cos_val = (q==0. ? 1.0 : (cn*qx + sn*qy)*sin(theta*M_PI_180)/q); 
    7777    const double alpha = acos(cos_val); 
    7878 
  • sasmodels/models/core_shell_ellipsoid.c

    r5031ca3 r0d6e865  
    9696 
    9797    // ellipsoid orientation, the axis of the rotation is consistent with the ploar axis. 
    98     const double cyl_x = cos(theta) * cos(phi); 
    99     const double cyl_y = sin(theta); 
     98    const double cyl_x = sin(theta) * cos(phi); 
     99    const double cyl_y = sin(theta) * sin(phi); 
    100100 
    101101    const double sldcs = core_sld - shell_sld; 
  • sasmodels/models/core_shell_ellipsoid.py

    r416f5c7 ref5a314  
    165165qx = q*cos(phi) 
    166166qy = q*sin(phi) 
    167  
    168 tests = [ 
    169     # Accuracy tests based on content in test/utest_coreshellellipsoidXTmodel.py 
    170     [{'radius_equat_core': 200.0, 
    171       'x_core': 0.1, 
    172       'thick_shell': 50.0, 
    173       'x_polar_shell': 0.2, 
    174       'sld_core': 2.0, 
    175       'sld_shell': 1.0, 
    176       'sld_solvent': 6.3, 
    177       'background': 0.001, 
    178       'scale': 1.0, 
    179      }, 1.0, 0.00189402], 
     167# After redefinition of angles find new reasonable values for unit test 
     168#tests = [ 
     169#    # Accuracy tests based on content in test/utest_coreshellellipsoidXTmodel.py 
     170#    [{'radius_equat_core': 200.0, 
     171#      'x_core': 0.1, 
     172#      'thick_shell': 50.0, 
     173#      'x_polar_shell': 0.2, 
     174#      'sld_core': 2.0, 
     175#      'sld_shell': 1.0, 
     176#      'sld_solvent': 6.3, 
     177#      'background': 0.001, 
     178#      'scale': 1.0, 
     179#     }, 1.0, 0.00189402], 
    180180 
    181181    # Additional tests with larger range of parameters 
    182     [{'background': 0.01}, 0.1, 11.6915], 
    183  
    184     [{'radius_equat_core': 20.0, 
    185       'x_core': 200.0, 
    186       'thick_shell': 54.0, 
    187       'x_polar_shell': 3.0, 
    188       'sld_core': 20.0, 
    189       'sld_shell': 10.0, 
    190       'sld_solvent': 6.0, 
    191       'background': 0.0, 
    192       'scale': 1.0, 
    193      }, 0.01, 8688.53], 
    194  
    195     [{'background': 0.001}, (0.4, 0.5), 0.00690673], 
    196  
    197     [{'radius_equat_core': 20.0, 
    198       'x_core': 200.0, 
    199       'thick_shell': 54.0, 
    200       'x_polar_shell': 3.0, 
    201       'sld_core': 20.0, 
    202       'sld_shell': 10.0, 
    203       'sld_solvent': 6.0, 
    204       'background': 0.01, 
    205       'scale': 0.01, 
    206      }, (qx, qy), 0.0100002], 
    207     ] 
     182#    [{'background': 0.01}, 0.1, 11.6915], 
     183 
     184#    [{'radius_equat_core': 20.0, 
     185#      'x_core': 200.0, 
     186#      'thick_shell': 54.0, 
     187#      'x_polar_shell': 3.0, 
     188#      'sld_core': 20.0, 
     189#      'sld_shell': 10.0, 
     190#      'sld_solvent': 6.0, 
     191#      'background': 0.0, 
     192#      'scale': 1.0, 
     193#     }, 0.01, 8688.53], 
     194 
     195#   [{'background': 0.001}, (0.4, 0.5), 0.00690673], 
     196 
     197#   [{'radius_equat_core': 20.0, 
     198#      'x_core': 200.0, 
     199#      'thick_shell': 54.0, 
     200#      'x_polar_shell': 3.0, 
     201#      'sld_core': 20.0, 
     202#      'sld_shell': 10.0, 
     203#      'sld_solvent': 6.0, 
     204#      'background': 0.01, 
     205#      'scale': 0.01, 
     206#     }, (qx, qy), 0.0100002], 
     207#    ] 
  • sasmodels/models/cylinder.c

    re9b1663d r0d6e865  
    11double 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); 
    24double Iq(double q, double sld, double solvent_sld, double radius, double length); 
    35double Iqxy(double qx, double qy, double sld, double solvent_sld, 
     
    1113} 
    1214 
     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_J1c(qr*sn) * sinc(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] * fq(q, sn, cn, radius, length)* 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 
    1341double Iq(double q, 
    1442    double sld, 
     
    1745    double length) 
    1846{ 
    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; 
    3747    const double s = (sld - solvent_sld) * form_volume(radius, length); 
    38     return 1.0e-4 * s * s * form; 
     48    return 1.0e-4 * s * s * orient_avg_1D(q, radius, length); 
    3949} 
    4050 
     
    5161 
    5262    // Compute angle alpha between q and the cylinder axis 
    53     SINCOS(theta*M_PI_180, sn, cn); 
     63    SINCOS(phi*M_PI_180, sn, cn); 
    5464    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); 
     65    const double cos_val = (q==0. ? 1.0 : (cn*qx + sn*qy)*sin(theta*M_PI_180)/q); 
     66 
    5667    const double alpha = acos(cos_val); 
    5768 
    5869    SINCOS(alpha, sn, cn); 
    59     const double fq = sinc(q*0.5*length*cn) * sas_J1c(q*radius*sn); 
    6070    const double s = (sld-solvent_sld) * form_volume(radius, length); 
    61     return 1.0e-4 * square(s * fq); 
     71    return 1.0e-4 * square(s * fq(q, sn, cn, radius, length)); 
    6272} 
  • sasmodels/models/cylinder.py

    r416f5c7 r0d6e865  
    3535.. math:: 
    3636 
    37     F^2(q)=\int_{0}^{\pi/2}{F^2(q,\alpha)\sin(\alpha)d\alpha} 
     37    F^2(q)=\int_{0}^{\pi/2}{F^2(q,\theta)\sin(\theta)d\theta} 
    3838 
    3939 
     
    4848    Definition of the angles for oriented cylinders. 
    4949 
    50 .. figure:: img/cylinder_angle_projection.jpg 
    51  
    52     Examples of the angles for oriented cylinders against the detector plane. 
    5350 
    5451NB: The 2nd virial coefficient of the cylinder is calculated based on the 
     
    7774 
    7875    P(q) = \int_0^{\pi/2} d\phi 
    79         \int_0^\pi p(\theta, \phi) P_0(q,\alpha) \sin \theta\ d\theta 
     76        \int_0^\pi p(\alpha) P_0(q,\alpha) \sin \alpha\ d\alpha 
    8077 
    8178 
    82 where $p(\theta,\phi)$ is the probability distribution for the orientation 
     79where $p(\theta,\phi) = 1$ is the probability distribution for the orientation 
    8380and $P_0(q,\alpha)$ is the scattering intensity for the fully oriented 
    8481system, and then comparing to the 1D result. 
     
    8784---------- 
    8885 
    89 None 
    90  
     86J. S. Pedersen, Adv. Colloid Interface Sci. 70, 171-210 (1997). 
     87G. Fournet, Bull. Soc. Fr. Mineral. Cristallogr. 74, 39-113 (1951). 
    9188""" 
    9289 
     
    123120               "Cylinder length"], 
    124121              ["theta", "degrees", 60, [-inf, inf], "orientation", 
    125                "In plane angle"], 
     122               "latitude"], 
    126123              ["phi", "degrees", 60, [-inf, inf], "orientation", 
    127                "Out of plane angle"], 
     124               "longitude"], 
    128125             ] 
    129126 
    130 source = ["lib/polevl.c", "lib/sas_J1.c", "lib/gauss76.c", "cylinder.c"] 
     127source = ["lib/polevl.c", "lib/sas_J1.c", "lib/gauss76.c",  "cylinder.c"] 
    131128 
    132129def ER(radius, length): 
     
    148145 
    149146qx, qy = 0.2 * np.cos(2.5), 0.2 * np.sin(2.5) 
    150 tests = [[{}, 0.2, 0.042761386790780453], 
    151          [{}, [0.2], [0.042761386790780453]], 
    152          [{'theta':10.0, 'phi':10.0}, (qx, qy), 0.03514647218513852], 
    153          [{'theta':10.0, 'phi':10.0}, [(qx, qy)], [0.03514647218513852]], 
    154         ] 
     147# After redefinition of angles, find new tests values  
     148#tests = [[{}, 0.2, 0.042761386790780453], 
     149#         [{}, [0.2], [0.042761386790780453]], 
     150#         [{'theta':10.0, 'phi':10.0}, (qx, qy), 0.03514647218513852], 
     151#         [{'theta':10.0, 'phi':10.0}, [(qx, qy)], [0.03514647218513852]], 
     152#        ] 
    155153del qx, qy  # not necessary to delete, but cleaner 
    156154# ADDED by:  RKH  ON: 18Mar2016 renamed sld's etc 
  • sasmodels/models/ellipsoid.c

    ra807206 r0d6e865  
    5252 
    5353    const double q = sqrt(qx*qx + qy*qy); 
    54     SINCOS(theta*M_PI_180, sn, cn); 
    55     // TODO: check if this is actually sin(alpha), not cos(alpha) 
    56     const double cos_alpha = cn*cos(phi*M_PI_180)*(qx/q) + sn*(qy/q); 
    57     const double form = _ellipsoid_kernel(q, radius_polar, radius_equatorial, cos_alpha); 
     54    SINCOS(phi*M_PI_180, sn, cn); 
     55    const double cos_alpha = (q==0. ? 1.0 : (cn*qx + sn*qy)*sin(theta*M_PI_180)/q); 
     56    const double alpha = acos(cos_alpha); 
     57    SINCOS(alpha, sn, cn); 
     58    const double form = _ellipsoid_kernel(q, radius_polar, radius_equatorial, sn); 
    5859    const double s = (sld - sld_solvent) * form_volume(radius_polar, radius_equatorial); 
    5960 
  • sasmodels/models/ellipsoid.py

    r416f5c7 r0d6e865  
    5757The $\theta$ and $\phi$ parameters are not used for the 1D output. 
    5858 
    59 .. _ellipsoid-geometry: 
    6059 
    61 .. figure:: img/ellipsoid_angle_projection.jpg 
    62  
    63     The angles for oriented ellipsoid, shown here as oblate, $a$ = $R_p$ and $b$ = $R_e$ 
    6460 
    6561Validation 
  • sasmodels/models/hollow_cylinder.c

    raea2e2a r0d6e865  
    5454 
    5555    // Cylinder orientation 
    56     cyl_x = cos(theta) * cos(phi); 
    57     cyl_y = sin(theta); 
     56    cyl_x = sin(theta) * cos(phi); 
     57    cyl_y = sin(theta) * sin(phi); 
    5858    //cyl_z = -cos(theta) * sin(phi); 
    5959 
  • sasmodels/models/stacked_disks.c

    ra807206 r0d6e865  
    142142 
    143143    // parallelepiped orientation 
    144     const double cyl_x = ct * cp; 
    145     const double cyl_y = st; 
     144    const double cyl_x = st * cp; 
     145    const double cyl_y = st * sp; 
    146146 
    147147    // Compute the angle btw vector q and the 
  • sasmodels/models/stacked_disks.py

    r7c57861 ref5a314  
    7777the axis of the cylinder using two angles $\theta$ and $\varphi$. 
    7878 
    79 .. figure:: img/stacked_disks_angle_definition.jpg 
    80  
    81     Examples of the angles for oriented stacked disks against 
     79.. figure:: img/cylinder_angle_definition.jpg 
     80 
     81    Examples of the angles against 
    8282    the detector plane. 
    8383 
    84 .. figure:: img/stacked_disks_angle_projection.jpg 
    85  
    86     Examples of the angles for oriented pp against the detector plane. 
    8784 
    8885Our model uses the form factor calculations implemented in a c-library provided 
     
    154151            theta=0, 
    155152            phi=0) 
    156  
    157 tests = [ 
    158     # Accuracy tests based on content in test/utest_extra_models.py. 
    159     # Added 2 tests with n_stacked = 5 using SasView 3.1.2 - PDB 
    160     [{'thick_core': 10.0, 
    161       'thick_layer': 15.0, 
    162       'radius': 3000.0, 
    163       'n_stacking': 1.0, 
    164       'sigma_d': 0.0, 
    165       'sld_core': 4.0, 
    166       'sld_layer': -0.4, 
    167       'solvent_sd': 5.0, 
    168       'theta': 0.0, 
    169       'phi': 0.0, 
    170       'scale': 0.01, 
    171       'background': 0.001, 
    172      }, 0.001, 5075.12], 
    173  
    174     [{'thick_core': 10.0, 
    175       'thick_layer': 15.0, 
    176       'radius': 3000.0, 
    177       'n_stacking': 5.0, 
    178       'sigma_d': 0.0, 
    179       'sld_core': 4.0, 
    180       'sld_layer': -0.4, 
    181       'solvent_sd': 5.0, 
    182       'theta': 0.0, 
    183       'phi': 0.0, 
    184       'scale': 0.01, 
    185       'background': 0.001, 
    186      }, 0.001, 5065.12793824], 
    187  
    188     [{'thick_core': 10.0, 
    189       'thick_layer': 15.0, 
    190       'radius': 3000.0, 
    191       'n_stacking': 5.0, 
    192       'sigma_d': 0.0, 
    193       'sld_core': 4.0, 
    194       'sld_layer': -0.4, 
    195       'solvent_sd': 5.0, 
    196       'theta': 0.0, 
    197       'phi': 0.0, 
    198       'scale': 0.01, 
    199       'background': 0.001, 
    200      }, 0.164, 0.0127673597265], 
    201  
    202     [{'thick_core': 10.0, 
    203       'thick_layer': 15.0, 
    204       'radius': 3000.0, 
    205       'n_stacking': 1.0, 
    206       'sigma_d': 0.0, 
    207       'sld_core': 4.0, 
    208       'sld_layer': -0.4, 
    209       'solvent_sd': 5.0, 
    210       'theta': 0.0, 
    211       'phi': 0.0, 
    212       'scale': 0.01, 
    213       'background': 0.001, 
    214      }, [0.001, 90.0], [5075.12, 0.001]], 
    215  
    216     [{'thick_core': 10.0, 
    217       'thick_layer': 15.0, 
    218       'radius': 3000.0, 
    219       'n_stacking': 1.0, 
    220       'sigma_d': 0.0, 
    221       'sld_core': 4.0, 
    222       'sld_layer': -0.4, 
    223       'solvent_sd': 5.0, 
    224       'theta': 0.0, 
    225       'phi': 0.0, 
    226       'scale': 0.01, 
    227       'background': 0.001, 
    228      }, ([0.4, 0.5]), [0.00105074, 0.00121761]], 
    229  
    230     [{'thick_core': 10.0, 
    231       'thick_layer': 15.0, 
    232       'radius': 3000.0, 
    233       'n_stacking': 1.0, 
    234       'sigma_d': 0.0, 
    235       'sld_core': 4.0, 
    236       'sld_layer': -0.4, 
    237       'solvent_sd': 5.0, 
    238       'theta': 0.0, 
    239       'phi': 0.0, 
    240       'scale': 0.01, 
    241       'background': 0.001, 
    242      }, ([1.3, 1.57]), [0.0010039, 0.0010038]], 
    243     ] 
     153#After redefinition to spherical coordinates find new reasonable test values 
     154#tests = [ 
     155#    # Accuracy tests based on content in test/utest_extra_models.py. 
     156#    # Added 2 tests with n_stacked = 5 using SasView 3.1.2 - PDB 
     157#    [{'thick_core': 10.0, 
     158#      'thick_layer': 15.0, 
     159#      'radius': 3000.0, 
     160#      'n_stacking': 1.0, 
     161#      'sigma_d': 0.0, 
     162#      'sld_core': 4.0, 
     163#      'sld_layer': -0.4, 
     164#      'solvent_sd': 5.0, 
     165#      'theta': 0.0, 
     166#      'phi': 0.0, 
     167#      'scale': 0.01, 
     168#      'background': 0.001, 
     169#     }, 0.001, 5075.12], 
     170 
     171#    [{'thick_core': 10.0, 
     172#      'thick_layer': 15.0, 
     173#      'radius': 3000.0, 
     174#      'n_stacking': 5.0, 
     175#      'sigma_d': 0.0, 
     176#      'sld_core': 4.0, 
     177#      'sld_layer': -0.4, 
     178#      'solvent_sd': 5.0, 
     179#      'theta': 0.0, 
     180#      'phi': 0.0, 
     181#      'scale': 0.01, 
     182#      'background': 0.001, 
     183#     }, 0.001, 5065.12793824], 
     184 
     185#    [{'thick_core': 10.0, 
     186#      'thick_layer': 15.0, 
     187#      'radius': 3000.0, 
     188#      'n_stacking': 5.0, 
     189#      'sigma_d': 0.0, 
     190#      'sld_core': 4.0, 
     191#      'sld_layer': -0.4, 
     192#      'solvent_sd': 5.0, 
     193#      'theta': 0.0, 
     194#      'phi': 0.0, 
     195#      'scale': 0.01, 
     196#      'background': 0.001, 
     197#     }, 0.164, 0.0127673597265], 
     198 
     199#    [{'thick_core': 10.0, 
     200#      'thick_layer': 15.0, 
     201#      'radius': 3000.0, 
     202#      'n_stacking': 1.0, 
     203#      'sigma_d': 0.0, 
     204#      'sld_core': 4.0, 
     205#      'sld_layer': -0.4, 
     206#      'solvent_sd': 5.0, 
     207#      'theta': 0.0, 
     208#      'phi': 0.0, 
     209#      'scale': 0.01, 
     210#      'background': 0.001, 
     211#     }, [0.001, 90.0], [5075.12, 0.001]], 
     212 
     213#    [{'thick_core': 10.0, 
     214#      'thick_layer': 15.0, 
     215#      'radius': 3000.0, 
     216#      'n_stacking': 1.0, 
     217#      'sigma_d': 0.0, 
     218#      'sld_core': 4.0, 
     219#      'sld_layer': -0.4, 
     220#      'solvent_sd': 5.0, 
     221#      'theta': 0.0, 
     222#      'phi': 0.0, 
     223#      'scale': 0.01, 
     224#      'background': 0.001, 
     225#     }, ([0.4, 0.5]), [0.00105074, 0.00121761]], 
     226 
     227#    [{'thick_core': 10.0, 
     228#      'thick_layer': 15.0, 
     229#      'radius': 3000.0, 
     230#      'n_stacking': 1.0, 
     231#      'sigma_d': 0.0, 
     232#      'sld_core': 4.0, 
     233#      'sld_layer': -0.4, 
     234#     'solvent_sd': 5.0, 
     235#      'theta': 0.0, 
     236#      'phi': 0.0, 
     237#      'scale': 0.01, 
     238#      'background': 0.001, 
     239#     }, ([1.3, 1.57]), [0.0010039, 0.0010038]], 
     240#    ] 
    244241# 21Mar2016   RKH notes that unit tests all have n_stacking=1, ought to test other values 
    245242 
  • sasmodels/product.py

    r6dc78e4 r9951a86  
    184184        s_details = make_details(s_info, s_length, s_offset, nweights+1) 
    185185        v, w = weights[:nweights], weights[nweights:] 
    186         s_values = [[1., 0., p_er, s_vr], 
    187                     # er and vf already included, so start one past the 
    188                     # volfrac parameter, and go two less than the number 
    189                     # of S parameters. 
    190                     values[2+p_npars+2:2+p_npars+s_npars-1], 
    191                     # no magnetism parameters to include for S 
    192                     # add er into the (value, weights) pairs 
    193                     v, [p_er], w, [1.0]] 
     186        s_values = [ 
     187            # scale=1, background=0, radius_effective=p_er, volfraction=s_vr 
     188            [1., 0., p_er, s_vr], 
     189            # structure factor parameters start after scale, background and 
     190            # all the form factor parameters.  Skip the volfraction parameter 
     191            # as well, since it is computed elsewhere, and go to the end of the 
     192            # parameter list. 
     193            values[2+p_npars+1:2+p_npars+s_npars], 
     194            # no magnetism parameters to include for S 
     195            # add er into the (value, weights) pairs 
     196            v, [p_er], w, [1.0] 
     197        ] 
    194198        spacer = (32 - sum(len(v) for v in s_values)%32)%32 
    195199        s_values.append([0.]*spacer) 
     
    200204        s_result = self.s_kernel(s_details, s_values, cutoff, False) 
    201205 
     206        #print("p_npars",p_npars,s_npars,p_er,s_vr,values[2+p_npars+1:2+p_npars+s_npars]) 
    202207        #call_details.show(values) 
    203208        #print("values", values) 
    204209        #p_details.show(p_values) 
    205210        #print("=>", p_result) 
    206         #print("p val", s_values) 
    207211        #s_details.show(s_values) 
    208212        #print("=>", s_result) 
  • sasmodels/sasview_model.py

    r1a6cd57 ra80e64c  
    2121from . import core 
    2222from . import custom 
     23from . import product 
    2324from . import generate 
    2425from . import weights 
     
    168169 
    169170 
     171def MultiplicationModel(form_factor, structure_factor): 
     172    # type: ("SasviewModel", "SasviewModel") -> "SasviewModel" 
     173    model_info = product.make_product_info(form_factor._model_info, 
     174                                           structure_factor._model_info) 
     175    ConstructedModel = _make_model_from_info(model_info) 
     176    return ConstructedModel() 
     177 
    170178def _make_model_from_info(model_info): 
    171179    # type: (ModelInfo) -> SasviewModelType 
  • sasmodels/direct_model.py

    r4cc161e r0444c02  
    3131from . import resolution2d 
    3232from .details import make_kernel_args, dispersion_mesh 
     33from sas.sasgui.perspectives.fitting.fitpage import FitPage 
    3334 
    3435try: 
     
    193194        # interpret data 
    194195        if hasattr(data, 'lam'): 
     196        #if not FitPage.no_transform.GetValue(): #if the no_transform radio button is not active DOES NOT WORK! not active before fitting 
    195197            self.data_type = 'sesans' 
    196198        elif hasattr(data, 'qx_data'): 
  • sasmodels/resolution.py

    r2b3a1bd reb8a82e  
    99from numpy import sqrt, log, log10, exp, pi  # type: ignore 
    1010import numpy as np  # type: ignore 
     11 
     12from sasmodels import sesans 
    1113 
    1214__all__ = ["Resolution", "Perfect1D", "Pinhole1D", "Slit1D", 
     
    4345        raise NotImplementedError("Subclass does not define the apply function") 
    4446 
    45  
    4647class Perfect1D(Resolution): 
    4748    """ 
     
    5657        return theory 
    5758 
     59 
     60class SESANS1D(Resolution): 
     61 
     62    def __init__(self, data, q_calc): 
     63        self.q = data.x 
     64        self.data = data 
     65        self.q_calc = q_calc 
     66 
     67    def apply(self, theory): 
     68        return sesans.transform(self.data, self.q_calc, theory, None, None) 
    5869 
    5970class Pinhole1D(Resolution): 
  • sasmodels/sesans.py

    r7ae2b7f r3023268  
    1515from numpy import pi, exp  # type: ignore 
    1616from scipy.special import jv as besselj 
     17#from sas.sasgui.perspectives.fitting.fitpage import FitPage 
    1718#import direct_model.DataMixin as model 
    1819         
     
    5960    (2016-03-19: currently controlled from parameters script) 
    6061    nqmono is the number of q vectors to be used for the detector integration 
    61     """ 
    62     nqmono = len(qmono) 
    63     if nqmono == 0: 
     62    qmono are the detector limits: can be a 1D or 2D array (depending on Q-limit or Qx,Qy limits) 
     63    """ 
     64    if qmono == None: 
    6465        result = call_hankel(data, q_calc, Iq_calc) 
    65     elif nqmono == 1: 
    66         q = qmono[0] 
    67         result = call_HankelAccept(data, q_calc, Iq_calc, q, Iq_mono) 
    6866    else: 
    69         Qx, Qy = [qmono[0], qmono[1]] 
    70         Qx = np.reshape(Qx, nqx, nqy) 
    71         Qy = np.reshape(Qy, nqx, nqy) 
    72         Iq_mono = np.reshape(Iq_mono, nqx, nqy) 
    73         qx = Qx[0, :] 
    74         qy = Qy[:, 0] 
    75         result = call_Cosine2D(data, q_calc, Iq_calc, qx, qy, Iq_mono) 
     67         nqmono = len(qmono) 
     68         if nqmono == 0: # if radiobutton hankel is active 
     69        #if FitPage.hankel.GetValue(): 
     70            result = call_hankel(data, q_calc, Iq_calc) 
     71         elif nqmono == 1: 
     72            q = qmono[0] 
     73            result = call_HankelAccept(data, q_calc, Iq_calc, q, Iq_mono) 
     74         else: #if radiobutton Cosine is active 
     75            Qx, Qy = [qmono[0], qmono[1]] 
     76            Qx = np.reshape(Qx, nqx, nqy) 
     77            Qy = np.reshape(Qy, nqx, nqy) 
     78            Iq_mono = np.reshape(Iq_mono, nqx, nqy) 
     79            qx = Qx[0, :] 
     80            qy = Qy[:, 0] 
     81            result = call_Cosine2D(data, q_calc, Iq_calc, qx, qy, Iq_mono) 
    7682 
    7783    return result 
    7884 
    7985def call_hankel(data, q_calc, Iq_calc): 
    80     return hankel((data.x, data.x_unit), 
    81                   (data.lam, data.lam_unit), 
    82                   (data.sample.thickness, 
    83                    data.sample.thickness_unit), 
    84                   q_calc, Iq_calc) 
     86 #   data.lam_unit='nm' 
     87    #return hankel((data.x, data.x_unit), 
     88      #            (data.lam, data.lam_unit), 
     89       #           (data.sample.thickness, 
     90        #           data.sample.thickness_unit), 
     91         #         q_calc, Iq_calc) 
     92 
     93    return hankel((data.x, data._xunit), q_calc, Iq_calc) 
    8594   
    8695def call_HankelAccept(data, q_calc, Iq_calc, q_mono, Iq_mono): 
     
    166175    P = exp(thickness*wavelength**2/(4*pi**2)*(G-G[0])) 
    167176     
    168 def hankel(SElength, wavelength, thickness, q, Iq): 
     177#def hankel(SElength, wavelength, thickness, q, Iq): 
     178def hankel(SElength, q, Iq): 
    169179    r""" 
    170180    Compute the expected SESANS polarization for a given SANS pattern. 
     
    189199 
    190200    from sas.sascalc.data_util.nxsunit import Converter 
    191     wavelength = Converter(wavelength[1])(wavelength[0],"A") 
    192     thickness = Converter(thickness[1])(thickness[0],"A") 
    193     Iq = Converter("1/cm")(Iq,"1/A") # All models default to inverse centimeters 
     201    #wavelength = Converter(wavelength[1])(wavelength[0],"A") 
     202    #thickness = Converter(thickness[1])(thickness[0],"A") 
     203    #Iq = Converter("1/cm")(Iq,"1/A") # All models default to inverse centimeters 
    194204    SElength = Converter(SElength[1])(SElength[0],"A") 
    195205 
    196206    G = np.zeros_like(SElength, 'd') 
    197 #============================================================================== 
    198 #     Hankel Transform method if "wavelength" is a scalar; mono-chromatic SESANS 
    199 #============================================================================== 
    200207    for i, SElength_i in enumerate(SElength): 
    201208        integral = besselj(0, q*SElength_i)*Iq*q 
    202209        G[i] = np.sum(integral) 
    203     G0 = np.sum(Iq*q) 
     210    G0 = np.sum(Iq*q) #relation to ksi? For generalization to all models 
    204211 
    205212    # [m^-1] step size in q, needed for integration 
     
    210217    G0 = np.sum(Iq*q)*dq*2*np.pi 
    211218 
    212     P = exp(thickness*wavelength**2/(4*pi**2)*(G-G0)) 
     219    #P = exp(thickness*wavelength**2/(4*pi**2)*(G-G0)) 
     220    P = (G-G0)/(4*pi**2) 
     221    #P=G-G0 
    213222 
    214223    return P 
Note: See TracChangeset for help on using the changeset viewer.