Changes in / [69ef533:e0de72f] in sasmodels


Ignore:
Location:
sasmodels/models
Files:
2 deleted
2 edited

Legend:

Unmodified
Added
Removed

  • sasmodels/models/core_shell_ellipsoid.c

    r2222134 r5031ca3  
    11double form_volume(double radius_equat_core, 
    2                    double radius_polar_core, 
    3                    double radius_equat_shell, 
    4                    double radius_polar_shell); 
     2                   double polar_core, 
     3                   double equat_shell, 
     4                   double polar_shell); 
    55double Iq(double q, 
    66          double radius_equat_core, 
    7           double radius_polar_core, 
    8           double radius_equat_shell, 
    9           double radius_polar_shell, 
    10           double sld_core, 
    11           double sld_shell, 
    12           double sld_solvent); 
     7          double x_core, 
     8          double thick_shell, 
     9          double x_polar_shell, 
     10          double core_sld, 
     11          double shell_sld, 
     12          double solvent_sld); 
    1313 
    1414 
    1515double Iqxy(double qx, double qy, 
    1616          double radius_equat_core, 
    17           double radius_polar_core, 
    18           double radius_equat_shell, 
    19           double radius_polar_shell, 
    20           double sld_core, 
    21           double sld_shell, 
    22           double sld_solvent, 
     17          double x_core, 
     18          double thick_shell, 
     19          double x_polar_shell, 
     20          double core_sld, 
     21          double shell_sld, 
     22          double solvent_sld, 
    2323          double theta, 
    2424          double phi); 
    … …  
    2626 
    2727double form_volume(double radius_equat_core, 
    28                    double radius_polar_core, 
    29                    double radius_equat_shell, 
    30                    double radius_polar_shell) 
     28                   double x_core, 
     29                   double thick_shell, 
     30                   double x_polar_shell) 
    3131{ 
    32     double vol = 4.0*M_PI/3.0*radius_equat_shell*radius_equat_shell*radius_polar_shell; 
     32    const double equat_shell = radius_equat_core + thick_shell; 
     33    const double polar_shell = radius_equat_core*x_core + thick_shell*x_polar_shell; 
     34    double vol = 4.0*M_PI/3.0*equat_shell*equat_shell*polar_shell; 
    3335    return vol; 
    3436} 
    3537 
    3638static double 
    37 core_shell_ellipsoid_kernel(double q, 
     39core_shell_ellipsoid_xt_kernel(double q, 
    3840          double radius_equat_core, 
    39           double radius_polar_core, 
    40           double radius_equat_shell, 
    41           double radius_polar_shell, 
    42           double sld_core, 
    43           double sld_shell, 
    44           double sld_solvent) 
     41          double x_core, 
     42          double thick_shell, 
     43          double x_polar_shell, 
     44          double core_sld, 
     45          double shell_sld, 
     46          double solvent_sld) 
    4547{ 
    46  
    47     //upper and lower integration limits 
    4848    const double lolim = 0.0; 
    4949    const double uplim = 1.0; 
    … …  
    5151    double summ = 0.0;   //initialize intergral 
    5252 
    53     const double delpc = sld_core - sld_shell;    //core - shell 
    54     const double delps = sld_shell - sld_solvent; //shell - solvent 
     53    const double delpc = core_sld - shell_sld; //core - shell 
     54    const double delps = shell_sld - solvent_sld; //shell - solvent 
     55 
     56 
     57    const double polar_core = radius_equat_core*x_core; 
     58    const double equat_shell = radius_equat_core + thick_shell; 
     59    const double polar_shell = radius_equat_core*x_core + thick_shell*x_polar_shell; 
    5560 
    5661    for(int i=0;i<N_POINTS_76;i++) { 
    … …  
    5863        double yyy = Gauss76Wt[i] * gfn4(zi, 
    5964                                  radius_equat_core, 
    60                                   radius_polar_core, 
    61                                   radius_equat_shell, 
    62                                   radius_polar_shell, 
     65                                  polar_core, 
     66                                  equat_shell, 
     67                                  polar_shell, 
    6368                                  delpc, 
    6469                                  delps, 
    … …  
    6873 
    6974    double answer = (uplim-lolim)/2.0*summ; 
    70  
    7175    //convert to [cm-1] 
    7276    answer *= 1.0e-4; 
    … …  
    7680 
    7781static double 
    78 core_shell_ellipsoid_kernel_2d(double q, double q_x, double q_y, 
     82core_shell_ellipsoid_xt_kernel_2d(double q, double q_x, double q_y, 
    7983          double radius_equat_core, 
    80           double radius_polar_core, 
    81           double radius_equat_shell, 
    82           double radius_polar_shell, 
    83           double sld_core, 
    84           double sld_shell, 
    85           double sld_solvent, 
     84          double x_core, 
     85          double thick_shell, 
     86          double x_polar_shell, 
     87          double core_sld, 
     88          double shell_sld, 
     89          double solvent_sld, 
    8690          double theta, 
    8791          double phi) 
    … …  
    9195    phi = phi * M_PI_180; 
    9296 
    93  
    9497    // ellipsoid orientation, the axis of the rotation is consistent with the ploar axis. 
    9598    const double cyl_x = cos(theta) * cos(phi); 
    9699    const double cyl_y = sin(theta); 
    97100 
    98     const double sldcs = sld_core - sld_shell; 
    99     const double sldss = sld_shell- sld_solvent; 
     101    const double sldcs = core_sld - shell_sld; 
     102    const double sldss = shell_sld- solvent_sld; 
    100103 
    101104    // Compute the angle btw vector q and the 
    … …  
    103106    const double cos_val = cyl_x*q_x + cyl_y*q_y; 
    104107 
    105     // Call the IGOR library function to get the kernel: MUST use gfn4 not gf2 because of the def of params. 
     108    const double polar_core = radius_equat_core*x_core; 
     109    const double equat_shell = radius_equat_core + thick_shell; 
     110    const double polar_shell = radius_equat_core*x_core + thick_shell*x_polar_shell; 
     111 
     112    // Call the IGOR library function to get the kernel: 
     113    // MUST use gfn4 not gf2 because of the def of params. 
    106114    double answer = gfn4(cos_val, 
    107115                  radius_equat_core, 
    108                   radius_polar_core, 
    109                   radius_equat_shell, 
    110                   radius_polar_shell, 
     116                  polar_core, 
     117                  equat_shell, 
     118                  polar_shell, 
    111119                  sldcs, 
    112120                  sldss, 
    … …  
    121129double Iq(double q, 
    122130          double radius_equat_core, 
    123           double radius_polar_core, 
    124           double radius_equat_shell, 
    125           double radius_polar_shell, 
    126           double sld_core, 
    127           double sld_shell, 
    128           double sld_solvent) 
     131          double x_core, 
     132          double thick_shell, 
     133          double x_polar_shell, 
     134          double core_sld, 
     135          double shell_sld, 
     136          double solvent_sld) 
    129137{ 
    130     double intensity = core_shell_ellipsoid_kernel(q, 
     138    double intensity = core_shell_ellipsoid_xt_kernel(q, 
    131139           radius_equat_core, 
    132            radius_polar_core, 
    133            radius_equat_shell, 
    134            radius_polar_shell, 
    135            sld_core, 
    136            sld_shell, 
    137            sld_solvent); 
     140           x_core, 
     141           thick_shell, 
     142           x_polar_shell, 
     143           core_sld, 
     144           shell_sld, 
     145           solvent_sld); 
    138146 
    139147    return intensity; 
    … …  
    143151double Iqxy(double qx, double qy, 
    144152          double radius_equat_core, 
    145           double radius_polar_core, 
    146           double radius_equat_shell, 
    147           double radius_polar_shell, 
    148           double sld_core, 
    149           double sld_shell, 
    150           double sld_solvent, 
     153          double x_core, 
     154          double thick_shell, 
     155          double x_polar_shell, 
     156          double core_sld, 
     157          double shell_sld, 
     158          double solvent_sld, 
    151159          double theta, 
    152160          double phi) 
    … …  
    154162    double q; 
    155163    q = sqrt(qx*qx+qy*qy); 
    156     double intensity = core_shell_ellipsoid_kernel_2d(q, qx/q, qy/q, 
     164    double intensity = core_shell_ellipsoid_xt_kernel_2d(q, qx/q, qy/q, 
    157165                       radius_equat_core, 
    158                        radius_polar_core, 
    159                        radius_equat_shell, 
    160                        radius_polar_shell, 
    161                        sld_core, 
    162                        sld_shell, 
    163                        sld_solvent, 
     166                       x_core, 
     167                       thick_shell, 
     168                       x_polar_shell, 
     169                       core_sld, 
     170                       shell_sld, 
     171                       solvent_sld, 
    164172                       theta, 
    165173                       phi); 

  • sasmodels/models/core_shell_ellipsoid.py

    r2222134 r5031ca3  
    11r""" 
    2 This model provides the form factor, $P(q)$, for a core shell ellipsoid (below) 
    3 where the form factor is normalized by the volume of the outer [CHECK]. 
     2An alternative version of $P(q)$ for the core_shell_ellipsoid 
     3having as parameters the core axial ratio X and a shell thickness, 
     4which are more often what we would like to determine. 
    45 
    5 .. math:: 
    6  
    7     P(q) = \text{scale} * \left<f^2\right>/V + \text{background} 
    8  
    9 where the volume $V = (4/3)\pi(r_\text{major outer} r_\text{minor outer}^2)$ 
    10 and the averaging $< >$ is applied over all orientations for 1D. 
    11  
    12 .. figure:: img/core_shell_ellipsoid_geometry.png 
    13  
    14     The returned value is in units of |cm^-1|, on absolute scale. 
     6This model is also better behaved when polydispersity is applied than the four 
     7independent radii in core_shell_ellipsoid model. 
    158 
    169Definition 
    1710---------- 
    1811 
    19 The form factor calculated is 
     12.. figure:: img/core_shell_ellipsoid_geometry.png 
    2013 
    21 .. math:: 
     14The geometric parameters of this model are 
    2215 
    23     P(q) &= \frac{\text{scale}}{V}\int_0^1 
    24         \left|F(q,r_\text{minor core},r_\text{major core},\alpha) 
    25         + F(q,r_\text{minor outer},r_\text{major outer},\alpha)\right|^2 
    26         d\alpha 
    27         + \text{background} 
     16*radius_equat_core =* equatorial core radius *= R_minor_core* 
    2817 
    29     \left|F(q,r_\text{minor},r_\text{major},\alpha)\right| 
    30         &=(4\pi/3)r_\text{major}r_\text{minor}^2 \Delta \rho \cdot (3j_1(u)/u) 
     18*X_core = polar_core / radius_equat_core = Rmajor_core / Rminor_core* 
    3119 
    32     u &= q\left[ r_\text{major}^2\alpha ^2 
    33                   + r_\text{minor}^2(1-\alpha ^2)\right]^{1/2} 
     20*Thick_shell = equat_outer - radius_equat_core = Rminor_outer - Rminor_core* 
    3421 
    35 where 
     22*XpolarShell = Tpolar_shell / Thick_shell = (Rmajor_outer - Rmajor_core)/ 
     23(Rminor_outer - Rminor_core)* 
    3624 
    37 .. math:: 
     25In terms of the original radii 
    3826 
    39     j_1(u)=(\sin x - x \cos x)/x^2 
     27*polar_core = radius_equat_core * X_core* 
    4028 
    41 To provide easy access to the orientation of the core-shell ellipsoid, 
    42 we define the axis of the solid ellipsoid using two angles $\theta$ and $\phi$. 
    43 These angles are defined as for 
    44 :ref:`cylinder orientation <cylinder-angle-definition>`. 
    45 The contrast is defined as SLD(core) - SLD(shell) and SLD(shell) - SLD(solvent). 
     29*equat_shell = radius_equat_core + Thick_shell* 
    4630 
    47 Note: It is the users' responsibility to ensure that shell radii are larger than 
    48 the core radii, especially if both are polydisperse, in which case the 
    49 core_shell_ellipsoid_xt model may be much better. 
     31*polar_shell = radius_equat_core * X_core + Thick_shell * XpolarShell* 
    5032 
     33(where we note that "shell" perhaps confusingly, relates to the outer radius) 
     34When *X_core < 1* the core is oblate; when *X_core > 1* it is prolate. 
     35*X_core = 1* is a spherical core. 
    5136 
    52 .. note:: 
    53     The 2nd virial coefficient of the solid ellipsoid is calculated based on 
    54     the *radius_a* (= *radius_polar_shell)* and *radius_b (= radius_equat_shell)* values, 
    55     and used as the effective radius for *S(Q)* when $P(Q) * S(Q)$ is applied. 
     37For a fixed shell thickness *XpolarShell = 1*, to scale the shell thickness 
     38pro-rata with the radius *XpolarShell = X_core*. 
    5639 
    57 .. figure:: img/core_shell_ellipsoid_angle_projection.jpg 
     40When including an $S(q)$, the radius in $S(q)$ is calculated to be that of 
     41a sphere with the same 2nd virial coefficient of the outer surface of the 
     42ellipsoid. This may have some undesirable effects if the aspect ratio of the 
     43ellipsoid is large (ie, if $X << 1$ or $X >> 1$ ), when the $S(q)$ 
     44- which assumes spheres - will not in any case be valid. 
    5845 
    59     The angles for oriented core_shell_ellipsoid. 
    60  
    61 Our model uses the form factor calculations implemented in a c-library provided 
    62 by the NIST Center for Neutron Research (Kline, 2006). 
     46If SAS data are in absolute units, and the SLDs are correct, then scale should 
     47be the total volume fraction of the "outer particle". When $S(q)$ is introduced 
     48this moves to the $S(q)$ volume fraction, and scale should then be 1.0, 
     49or contain some other units conversion factor (for example, if you have SAXS data). 
    6350 
    6451References 
    6552---------- 
    6653 
    67 M Kotlarchyk, S H Chen, *J. Chem. Phys.*, 79 (1983) 2461 
     54R K Heenan, 2015, reparametrised the core_shell_ellipsoid model 
    6855 
    69 S J Berr, *Phys. Chem.*, 91 (1987) 4760 
    7056""" 
    7157 
    7258from numpy import inf, sin, cos, pi 
    7359 
    74 name = "core_shell_ellipsoid" 
     60name = "core_shell_ellipsoid_xt" 
    7561title = "Form factor for an spheroid ellipsoid particle with a core shell structure." 
    7662description = """ 
    77     [SpheroidCoreShellModel] Calculates the form factor for an spheroid 
    78     ellipsoid particle with a core_shell structure. 
    79     The form factor is averaged over all possible 
    80     orientations of the ellipsoid such that P(q) 
    81     = scale*<f^2>/Vol + bkg, where f is the 
    82     single particle scattering amplitude. 
    83     [Parameters]: 
    84     radius_equat_core = equatorial radius of core, Rminor_core, 
    85     radius_polar_core = polar radius of core, Rmajor_core, 
    86     radius_equat_shell = equatorial radius of shell, Rminor_outer, 
    87     radius_polar_shell = polar radius of shell, Rmajor_outer, 
    88     sld_core = scattering length density of core, 
    89     sld_shell = scattering length density of shell, 
    90     sld_solvent = scattering length density of solvent, 
    91     background = Incoherent bkg 
    92     scale =scale 
    93     Note:It is the users' responsibility to ensure 
    94     that shell radii are larger than core radii, 
    95     especially if both are polydisperse. 
    96     oblate: polar radius < equatorial radius 
    97     prolate :  polar radius > equatorial radius 
     63        [core_shell_ellipsoid_xt] Calculates the form factor for an spheroid 
     64        ellipsoid particle with a core_shell structure. 
     65        The form factor is averaged over all possible 
     66        orientations of the ellipsoid such that P(q) 
     67        = scale*<f^2>/Vol + bkg, where f is the 
     68        single particle scattering amplitude. 
     69        [Parameters]: 
     70        radius_equat_core = equatorial radius of core, 
     71        x_core = ratio of core polar/equatorial radii, 
     72        thick_shell = equatorial radius of outer surface, 
     73        x_polar_shell = ratio of polar shell thickness to equatorial shell thickness, 
     74        sld_core = SLD_core 
     75        sld_shell = SLD_shell 
     76        sld_solvent = SLD_solvent 
     77        background = Incoherent bkg 
     78        scale =scale 
     79        Note:It is the users' responsibility to ensure 
     80        that shell radii are larger than core radii. 
     81        oblate: polar radius < equatorial radius 
     82        prolate :  polar radius > equatorial radius - this new model will make this easier 
     83        and polydispersity integrals more logical (as previously the shell could disappear). 
    9884    """ 
    9985category = "shape:ellipsoid" 
    10086 
    10187# pylint: disable=bad-whitespace, line-too-long 
    102 #   ["name", "units", default, [lower, upper], "type", "description"], 
     88#             ["name", "units", default, [lower, upper], "type", "description"], 
    10389parameters = [ 
    104     ["radius_equat_core",  "Ang",      200,   [0, inf],    "volume",      "Equatorial radius of core, r minor core"], 
    105     ["radius_polar_core",  "Ang",       10,   [0, inf],    "volume",      "Polar radius of core, r major core"], 
    106     ["radius_equat_shell", "Ang",      250,   [0, inf],    "volume",      "Equatorial radius of shell, r minor outer"], 
    107     ["radius_polar_shell", "Ang",       30,   [0, inf],    "volume",      "Polar radius of shell, r major outer"], 
    108     ["sld_core",    "1e-6/Ang^2", 2,   [-inf, inf], "sld",         "Core scattering length density"], 
    109     ["sld_shell",   "1e-6/Ang^2", 1,   [-inf, inf], "sld",         "Shell scattering length density"], 
    110     ["sld_solvent", "1e-6/Ang^2", 6.3, [-inf, inf], "sld",         "Solvent scattering length density"], 
    111     ["theta",       "degrees",    0,   [-inf, inf], "orientation", "Oblate orientation wrt incoming beam"], 
    112     ["phi",         "degrees",    0,   [-inf, inf], "orientation", "Oblate orientation in the plane of the detector"], 
     90    ["radius_equat_core","Ang",     20,   [0, inf],    "volume",      "Equatorial radius of core"], 
     91    ["x_core",        "None",       3,   [0, inf],    "volume",      "axial ratio of core, X = r_polar/r_equatorial"], 
     92    ["thick_shell",   "Ang",       30,   [0, inf],    "volume",      "thickness of shell at equator"], 
     93    ["x_polar_shell", "",           1,   [0, inf],    "volume",      "ratio of thickness of shell at pole to that at equator"], 
     94    ["sld_core",      "1e-6/Ang^2", 2,   [-inf, inf], "sld",         "Core scattering length density"], 
     95    ["sld_shell",     "1e-6/Ang^2", 1,   [-inf, inf], "sld",         "Shell scattering length density"], 
     96    ["sld_solvent",   "1e-6/Ang^2", 6.3, [-inf, inf], "sld",         "Solvent scattering length density"], 
     97    ["theta",         "degrees",    0,   [-inf, inf], "orientation", "Oblate orientation wrt incoming beam"], 
     98    ["phi",           "degrees",    0,   [-inf, inf], "orientation", "Oblate orientation in the plane of the detector"], 
    11399    ] 
    114100# pylint: enable=bad-whitespace, line-too-long 
    115101 
    116 source = ["lib/sph_j1c.c", "lib/gfn.c", "lib/gauss76.c", "core_shell_ellipsoid.c"] 
     102source = ["lib/sph_j1c.c", "lib/gfn.c", "lib/gauss76.c", 
     103          "core_shell_ellipsoid_xt.c"] 
    117104 
    118 def ER(radius_equat_core, radius_polar_core, radius_equat_shell, radius_polar_shell): 
     105def ER(radius_equat_core, x_core, thick_shell, x_polar_shell): 
    119106    """ 
    120107        Returns the effective radius used in the S*P calculation 
    121108    """ 
    122109    from .ellipsoid import ER as ellipsoid_ER 
    123     return ellipsoid_ER(radius_polar_shell, radius_equat_shell) 
     110    polar_outer = radius_equat_core*x_core + thick_shell*x_polar_shell 
     111    equat_outer = radius_equat_core + thick_shell 
     112    return ellipsoid_ER(polar_outer, equat_outer) 
    124113 
    125114 
    126 demo = dict(scale=1, background=0.001, 
    127             radius_equat_core=200.0, 
    128             radius_polar_core=10.0, 
    129             radius_equat_shell=250.0, 
    130             radius_polar_shell=30.0, 
     115demo = dict(scale=0.05, background=0.001, 
     116            radius_equat_core=20.0, 
     117            x_core=3.0, 
     118            thick_shell=30.0, 
     119            x_polar_shell=1.0, 
    131120            sld_core=2.0, 
    132121            sld_shell=1.0, 
    … …  
    141130 
    142131tests = [ 
    143     # Accuracy tests based on content in test/utest_other_models.py 
     132    # Accuracy tests based on content in test/utest_coreshellellipsoidXTmodel.py 
    144133    [{'radius_equat_core': 200.0, 
    145       'radius_polar_core': 20.0, 
    146       'radius_equat_shell': 250.0, 
    147       'radius_polar_shell': 30.0, 
     134      'x_core': 0.1, 
     135      'thick_shell': 50.0, 
     136      'x_polar_shell': 0.2, 
    148137      'sld_core': 2.0, 
    149138      'sld_shell': 1.0, 
    … …  
    154143 
    155144    # Additional tests with larger range of parameters 
    156     [{'background': 0.01}, 0.1, 8.86741], 
     145    [{'background': 0.01}, 0.1, 11.6915], 
    157146 
    158147    [{'radius_equat_core': 20.0, 
    159       'radius_polar_core': 200.0, 
    160       'radius_equat_shell': 54.0, 
    161       'radius_polar_shell': 3.0, 
     148      'x_core': 200.0, 
     149      'thick_shell': 54.0, 
     150      'x_polar_shell': 3.0, 
    162151      'sld_core': 20.0, 
    163152      'sld_shell': 10.0, 
    … …  
    165154      'background': 0.0, 
    166155      'scale': 1.0, 
    167      }, 0.01, 26150.4], 
     156     }, 0.01, 8688.53], 
    168157 
    169     [{'background': 0.001}, (0.4, 0.5), 0.00170471], 
     158    [{'background': 0.001}, (0.4, 0.5), 0.00690673], 
    170159 
    171160    [{'radius_equat_core': 20.0, 
    172       'radius_polar_core': 200.0, 
    173       'radius_equat_shell': 54.0, 
    174       'radius_polar_shell': 3.0, 
     161      'x_core': 200.0, 
     162      'thick_shell': 54.0, 
     163      'x_polar_shell': 3.0, 
    175164      'sld_core': 20.0, 
    176165      'sld_shell': 10.0, 
    … …  
    178167      'background': 0.01, 
    179168      'scale': 0.01, 
    180      }, (qx, qy), 0.105764], 
     169     }, (qx, qy), 0.0100002], 
    181170    ] 
Note: See TracChangeset for help on using the changeset viewer.