source: sasmodels/sasmodels/models/elliptical_cylinder.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: 6.4 KB
RevLine 
[c612162]1# pylint: disable=line-too-long
[a8b3cdb]2r"""
3
[5111921]4.. figure:: img/elliptical_cylinder_geometry.png
[a8b3cdb]5
[40a87fa]6   Elliptical cylinder geometry $a = r_\text{minor}$
7   and $\nu = r_\text{major} / r_\text{minor}$ is the *axis_ratio*.
[a8b3cdb]8
9The function calculated is
10
11.. math::
[fa8011eb]12
[40a87fa]13    I(\vec q)=\frac{1}{V_\text{cyl}}\int{d\psi}\int{d\phi}\int{
[fcb33e4]14        p(\theta,\phi,\psi)F^2(\vec q,\alpha,\psi)\sin(\alpha)d\alpha}
[a8b3cdb]15
16with the functions
17
18.. math::
[fa8011eb]19
[fcb33e4]20    F(q,\alpha,\psi) = 2\frac{J_1(a)\sin(b)}{ab}
[40a87fa]21
22where
23
24.. math::
25
[fcb33e4]26    a = qr'\sin(\alpha)
[404ebbd]27
[fcb33e4]28    b = q\frac{L}{2}\cos(\alpha)
[404ebbd]29
[fcb33e4]30    r'=\frac{r_{minor}}{\sqrt{2}}\sqrt{(1+\nu^{2}) + (1-\nu^{2})cos(\psi)}
[40a87fa]31
[a8b3cdb]32
[fcb33e4]33and the angle $\psi$ is defined as the orientation of the major axis of the
[40a87fa]34ellipse with respect to the vector $\vec q$. The angle $\alpha$ is the angle
35between the axis of the cylinder and $\vec q$.
[a8b3cdb]36
37
[eda8b30]38For 1D scattering, with no preferred orientation, the form factor is averaged over all possible orientations and normalized
[40a87fa]39by the particle volume
[a8b3cdb]40
41.. math::
42
[40a87fa]43    P(q) = \text{scale}  <F^2> / V
[a8b3cdb]44
[74768cb]45For 2d data the orientation of the particle is required, described using a different set
46of angles as in the diagrams below, for further details of the calculation and angular
[eda8b30]47dispersions  see :ref:`orientation` .
[40a87fa]48
[a8b3cdb]49
[15a90c1]50.. figure:: img/elliptical_cylinder_angle_definition.png
[a8b3cdb]51
[eda8b30]52    Note that the angles here are not the same as in the equations for the scattering function.
53    Rotation $\theta$, initially in the $xz$ plane, is carried out first, then
54    rotation $\phi$ about the $z$ axis, finally rotation $\Psi$ is now around the axis of the cylinder.
55    The neutron or X-ray beam is along the $z$ axis.
[a8b3cdb]56
[15a90c1]57.. figure:: img/elliptical_cylinder_angle_projection.png
[a8b3cdb]58
[40a87fa]59    Examples of the angles for oriented elliptical cylinders against the
[15a90c1]60    detector plane, with $\Psi$ = 0.
[a8b3cdb]61
[404ebbd]62The $\theta$ and $\phi$ parameters to orient the cylinder only appear in the model when fitting 2d data.
[eda8b30]63
[9802ab3]64
[40a87fa]65NB: The 2nd virial coefficient of the cylinder is calculated based on the
66averaged radius $(=\sqrt{r_\text{minor}^2 * \text{axis ratio}})$ and length
67values, and used as the effective radius for $S(Q)$ when $P(Q)*S(Q)$ is applied.
[a8b3cdb]68
69
70Validation
71----------
72
[40a87fa]73Validation of our code was done by comparing the output of the 1D calculation
74to the angular average of the output of the 2D calculation over all possible
75angles.
[a8b3cdb]76
[40a87fa]77In the 2D average, more binning in the angle $\phi$ is necessary to get the
78proper result. The following figure shows the results of the averaging by
79varying the number of angular bins.
[a8b3cdb]80
[5111921]81.. figure:: img/elliptical_cylinder_averaging.png
[a8b3cdb]82
[fa8011eb]83    The intensities averaged from 2D over different numbers of bins and angles.
[a8b3cdb]84
[e65a3e7]85References
86----------
[a8b3cdb]87
[0507e09]88.. [#] L A Feigin and D I Svergun, *Structure Analysis by Small-Angle X-Ray and Neutron Scattering*, Plenum, New York, (1987) [see table 3.4]
89.. [#] L. Onsager, *Ann. New York Acad. Sci.*, 51 (1949) 627-659
90
91Source
92------
93
94`elliptical_cylinder.py <https://github.com/SasView/sasmodels/blob/master/sasmodels/models/elliptical_cylinder.py>`_
95
96`elliptical_cylinder.c <https://github.com/SasView/sasmodels/blob/master/sasmodels/models/elliptical_cylinder.c>`_
[fcb33e4]97
98Authorship and Verification
99----------------------------
100
101* **Author:**
[404ebbd]102* **Last Modified by:**
[fcb33e4]103* **Last Reviewed by:**  Richard Heenan - corrected equation in docs **Date:** December 21, 2016
[0507e09]104* **Source added by :** Steve King **Date:** March 25, 2019
[a8b3cdb]105"""
106
[2d81cfe]107import numpy as np
[0b56f38]108from numpy import pi, inf, sqrt, sin, cos
[a8b3cdb]109
110name = "elliptical_cylinder"
111title = "Form factor for an elliptical cylinder."
112description = """
[c612162]113    Form factor for an elliptical cylinder.
114    See L A Feigin and D I Svergun, Structure Analysis by Small-Angle X-Ray and Neutron Scattering, Plenum, New York, (1987).
[a8b3cdb]115"""
116category = "shape:cylinder"
117
118# pylint: disable=bad-whitespace, line-too-long
119#             ["name", "units", default, [lower, upper], "type","description"],
[a807206]120parameters = [["radius_minor",     "Ang",        20.0,  [0, inf],    "volume",      "Ellipse minor radius"],
[42356c8]121              ["axis_ratio",   "",          1.5,   [1, inf],    "volume",      "Ratio of major radius over minor radius"],
[a8b3cdb]122              ["length",      "Ang",        400.0, [1, inf],    "volume",      "Length of the cylinder"],
[42356c8]123              ["sld",         "1e-6/Ang^2", 4.0,   [-inf, inf], "sld",         "Cylinder scattering length density"],
124              ["sld_solvent", "1e-6/Ang^2", 1.0,   [-inf, inf], "sld",         "Solvent scattering length density"],
[9b79f29]125              ["theta",       "degrees",    90.0,  [-360, 360], "orientation", "cylinder axis to beam angle"],
126              ["phi",         "degrees",    0,     [-360, 360], "orientation", "rotation about beam"],
127              ["psi",         "degrees",    0,     [-360, 360], "orientation", "rotation about cylinder axis"]]
[a8b3cdb]128
129# pylint: enable=bad-whitespace, line-too-long
130
[74768cb]131source = ["lib/polevl.c", "lib/sas_J1.c", "lib/gauss76.c", "elliptical_cylinder.c"]
[71b751d]132have_Fq = True
[a34b811]133radius_effective_modes = [
[b297ba9]134    "equivalent cylinder excluded volume",
135    "equivalent volume sphere", "average radius", "min radius", "max radius",
[ee60aa7]136    "equivalent circular cross-section",
137    "half length", "half min dimension", "half max dimension", "half diagonal",
138    ]
[a8b3cdb]139
[a807206]140demo = dict(scale=1, background=0, radius_minor=100, axis_ratio=1.5, length=400.0,
[40a87fa]141            sld=4.0, sld_solvent=1.0, theta=10.0, phi=20, psi=30,
142            theta_pd=10, phi_pd=2, psi_pd=3)
[a8b3cdb]143
[404ebbd]144def random():
[b297ba9]145    """Return a random parameter set for the model."""
[404ebbd]146    # V = pi * radius_major * radius_minor * length;
[2d81cfe]147    volume = 10**np.random.uniform(3, 9)
[404ebbd]148    length = 10**np.random.uniform(1, 3)
[d9ec8f9]149    axis_ratio = 10**np.random.uniform(0, 2)
[b297ba9]150    radius_minor = sqrt(volume/length/axis_ratio)
[2d81cfe]151    volfrac = 10**np.random.uniform(-4, -2)
[404ebbd]152    pars = dict(
153        #background=0, sld=0, sld_solvent=1,
[2d81cfe]154        scale=1e9*volfrac/volume,
[404ebbd]155        length=length,
156        radius_minor=radius_minor,
157        axis_ratio=axis_ratio,
158    )
159    return pars
160
[0b56f38]161q = 0.1
162# april 6 2017, rkh added a 2d unit test, NOT READY YET pull #890 branch assume correct!
163qx = q*cos(pi/6.0)
164qy = q*sin(pi/6.0)
[a8b3cdb]165
[40a87fa]166tests = [
[b297ba9]167    #[{'radius_minor': 20.0, 'axis_ratio': 1.5, 'length':400.0}, 'ER', 79.89245454155024],
168    #[{'radius_minor': 20.0, 'axis_ratio': 1.2, 'length':300.0}, 'VR', 1],
[40a87fa]169
170    # The SasView test result was 0.00169, with a background of 0.001
[a807206]171    [{'radius_minor': 20.0, 'axis_ratio': 1.5, 'sld': 4.0, 'length':400.0,
[40a87fa]172      'sld_solvent':1.0, 'background':0.0},
173     0.001, 675.504402],
[2d81cfe]174    #[{'theta':80., 'phi':10.}, (qx, qy), 7.88866563001 ],
[40a87fa]175]
Note: See TracBrowser for help on using the repository browser.