source: sasmodels/sasmodels/models/_spherepy.py @ a34b811

ticket-1257-vesicle-productticket_1156ticket_822_more_unit_tests
Last change on this file since a34b811 was a34b811, checked in by Paul Kienzle <pkienzle@…>, 5 years ago

use radius_effective/radius_effective_mode/radius_effective_modes consistently throughout the code

  • Property mode set to 100644
File size: 4.2 KB
Line 
1r"""
2For information about polarised and magnetic scattering, see
3the :ref:`magnetism` documentation.
4
5Definition
6----------
7
8The 1D scattering intensity is calculated in the following way (Guinier, 1955)
9
10.. math::
11
12    I(q) = \frac{\text{scale}}{V} \cdot \left[
13        3V(\Delta\rho) \cdot \frac{\sin(qr) - qr\cos(qr))}{(qr)^3}
14        \right]^2 + \text{background}
15
16where *scale* is a volume fraction, $V$ is the volume of the scatterer,
17$r$ is the radius of the sphere, *background* is the background level and
18*sld* and *sld_solvent* are the scattering length densities (SLDs) of the
19scatterer and the solvent respectively.
20
21Note that if your data is in absolute scale, the *scale* should represent
22the volume fraction (which is unitless) if you have a good fit. If not,
23it should represent the volume fraction times a factor (by which your data
24might need to be rescaled).
25
26The 2D scattering intensity is the same as above, regardless of the
27orientation of $\vec q$.
28
29Validation
30----------
31
32Validation of our code was done by comparing the output of the 1D model
33to the output of the software provided by the NIST (Kline, 2006).
34
35References
36----------
37
38A Guinier and G. Fournet, *Small-Angle Scattering of X-Rays*,
39John Wiley and Sons, New York, (1955)
40
41Source
42------
43
44`_spherepy.py <https://github.com/SasView/sasmodels/blob/master/sasmodels/models/_spherepy.py>`_
45`sphere.c <https://github.com/SasView/sasmodels/blob/master/sasmodels/models/sphere.c>`_
46
47Authorship and Verification
48----------------------------
49
50* **Author: P Kienzle**
51* **Last Modified by:**
52* **Last Reviewed by:** S King and P Parker **Date:** 2013/09/09 and 2014/01/06
53* **Source added by :** Steve King **Date:** March 25, 2019
54"""
55
56import numpy as np
57from numpy import pi, inf, sin, cos, sqrt, log
58
59name = " _sphere (python)"
60title = "PAK testing ideas for Spheres with uniform scattering length density"
61description = """\
62P(q)=(scale/V)*[3V(sld-sld_solvent)*(sin(qr)-qr cos(qr))
63                /(qr)^3]^2 + background
64    r: radius of sphere
65    V: The volume of the scatter
66    sld: the SLD of the sphere
67    sld_solvent: the SLD of the solvent
68"""
69category = "shape:sphere"
70
71#             ["name", "units", default, [lower, upper], "type","description"],
72parameters = [["sld", "1e-6/Ang^2", 1, [-inf, inf], "",
73               "Layer scattering length density"],
74              ["sld_solvent", "1e-6/Ang^2", 6, [-inf, inf], "",
75               "Solvent scattering length density"],
76              ["radius", "Ang", 50, [0, inf], "volume",
77               "Sphere radius"],
78             ]
79
80
81def form_volume(radius):
82    """Calculate volume for sphere"""
83    return 1.333333333333333 * pi * radius ** 3
84
85def radius_effective(mode, radius):
86    """Calculate R_eff for sphere"""
87    return radius if mode else 0.
88
89def Iq(q, sld, sld_solvent, radius):
90    """Calculate I(q) for sphere"""
91    #print "q",q
92    #print "sld,r",sld,sld_solvent,radius
93    qr = q * radius
94    sn, cn = sin(qr), cos(qr)
95    ## The natural expression for the bessel function is the following:
96    ##     bes = 3 * (sn-qr*cn)/qr**3 if qr>0 else 1
97    ## however, to support vector q values we need to handle the conditional
98    ## as a vector, which we do by first evaluating the full expression
99    ## everywhere, then fixing it up where it is broken.   We should probably
100    ## set numpy to ignore the 0/0 error before we do though...
101    bes = 3 * (sn - qr * cn) / qr ** 3 # may be 0/0 but we fix that next line
102    bes[qr == 0] = 1
103    fq = bes * (sld - sld_solvent) * form_volume(radius)
104    return 1.0e-4 * fq ** 2
105Iq.vectorized = True  # Iq accepts an array of q values
106
107def sesans(z, sld, sld_solvent, radius):
108    """
109    Calculate SESANS-correlation function for a solid sphere.
110
111    Wim Bouwman after formulae Timofei Kruglov J.Appl.Cryst. 2003 article
112    """
113    d = z / radius
114    g = np.zeros_like(z)
115    g[d == 0] = 1.
116    low = ((d > 0) & (d < 2))
117    dlow = d[low]
118    dlow2 = dlow ** 2
119    g[low] = (sqrt(1 - dlow2/4.) * (1 + dlow2/8.)
120              + dlow2/2.*(1 - dlow2/16.) * log(dlow / (2. + sqrt(4. - dlow2))))
121    return g
122sesans.vectorized = True  # sesans accepts an array of z values
123
124demo = dict(scale=1, background=0,
125            sld=6, sld_solvent=1,
126            radius=120,
127            radius_pd=.2, radius_pd_n=45)
Note: See TracBrowser for help on using the repository browser.