Fix limits on angles and angular dispersity

doc/guide/orientation/orientation.rst

@@ -103,6 +103,14 @@ to understand the 2d patterns fully. A number of different shapes of
 distribution are available, as described for size dispersity, see
 :ref:`polydispersityhelp`.
 
+The parameter limits for the orientation parameters are generous, with
+$\theta \in [-180, 180]$ and $\phi, \Psi \in [-360, 360]$. Although
+there will be multiple equivalent solutions in the space, these limits
+allow you to set bounds around the solution even when it is near a pole.
+For example, if the solution is within 4 degrees of $\phi = 179 \pm 4$,
+having strict bounds would force you to search both in $[175, 180]$ and
+in $[-180, -177]$ instead of searching within $[175, 183]$.
+
 Given that the angular dispersion distribution is defined in cartesian space,
 over a cube defined by
 
@@ -130,9 +138,9 @@ will be the case for $\theta=90$ and $\Psi=90$) the scattering pattern
 should be circularly symmetric, but it will go to zero at $q_x = 0$ due to the
 $\sin(\Delta\theta)$ correction. This problem does not appear for a shape
 that is tumbling freely around the $a$ axis, with $\Delta\phi$ uniform in
-$[-180, 180]$, so swap the $a$ and $b$ axes so $\Delta\theta < \Delta\phi$
-and adjust $\Psi$ by 90. This works with the current sasmodels shapes due to
-symmetry.
+$[-180, 180]$. For shapes such as parallelepiped or triaxial ellipsoid where
+the choice of a, b and c is arbitrary, you can minimize the problem by
+relabelling the axes and adjusting the orientation.
 
 Alternative projections were considered.
 The `sinusoidal projection <https://en.wikipedia.org/wiki/Sinusoidal_projection>`_

sasmodels/TwoYukawa/CalTYSk.py

@@ -1,8 +1,8 @@
 import numpy as np
-from numpy import pi, mean
+from numpy import mean, pi
 
-from .Ecoefficient import TYCoeff
 from .CalcRealRoot import CalRealRoot
+from .Ecoefficient import TYCoeff
 from .TInvFourier import TInvFourier
 
 # Supplied Q vector must go out to at least this value to calculate g(r).

sasmodels/TwoYukawa/CalcRealRoot.py

@@ -3,6 +3,7 @@
 from .Ecoefficient import TYCoeff
 from .Epoly import make_poly
 
+
 def CalRealRoot(coeff: TYCoeff):
 
     Ecoefficient = coeff.Ecoefficient

sasmodels/TwoYukawa/Ecoefficient.py

@@ -1,7 +1,7 @@
 from typing import Callable, Tuple
 
 import numpy as np
-from numpy import exp, pi, cos, sin, cosh
+from numpy import cos, cosh, exp, pi, sin
 from numpy.typing import NDArray
 
 # CalCoeff.m
@@ -194,17 +194,17 @@ def ABC12(d1, d2):
               Ccd2_22*Cdd1_12*d1*d2**2 -
               Ccd1_21*Cdd2_22*d1*d2**2 -
               Ccd2_22*d2*k1 + Ccd1_21*d1*k2)))/
-        ((d1*((
+        (d1*(
               (-Ccd1_21)*Ccd2_12*d2 +
-                Ccd1_11*Ccd2_22*d2)))))
+                Ccd1_11*Ccd2_22*d2)))
         C2 = ((-((Ccd2_12*d2*
                 (((-Cd1_1)*d1 - Cdd1_11*d1**2 -
                     Cdd1_12*d1*d2 + k1)) -
               Ccd1_11*d1*
                 (((-Cd2_2)*d2 - Cdd2_12*d1*d2 -
                     Cdd2_22*d2**2 + k2))))) /
-        (((-Ccd1_21)*Ccd2_12*d1*d2 +
-            Ccd1_11*Ccd2_22*d1*d2)))
+        ((-Ccd1_21)*Ccd2_12*d1*d2 +
+            Ccd1_11*Ccd2_22*d1*d2))
         return A, B, C1, C2
     self.ABC12 = ABC12
 
@@ -370,10 +370,10 @@ def tau_d21(s):
 
     E1d02 = 12*c1F01*phi*sigma_d01(z1) - 12*c1F01*exp(-z1)*phi*tau_d01(z1)
 
-    E1d11 = (((12*c1F10*phi*sigma_d01(z1) + \
+    E1d11 = (12*c1F10*phi*sigma_d01(z1) + \
           12*c1F01*phi*sigma_d10(z1) - \
           12*c1F10*exp(-z1)*phi*tau_d01(z1) - \
-          12*c1F01*exp(-z1)*phi*tau_d10(z1))))
+          12*c1F01*exp(-z1)*phi*tau_d10(z1))
     E1d12 = (-c1F01 + 12*c1F11*phi*sigma_d01(z1) +
              12*c1F01*phi*sigma_d11(z1) -
              12*c1F11*exp(-z1)*phi*tau_d01(z1) -

sasmodels/TwoYukawa/TFourier.py

@@ -1,6 +1,7 @@
-from numpy import exp, pi, arange, linspace, abs, ceil, log2, interp
+from numpy import abs, arange, ceil, exp, interp, linspace, log2, pi
 from numpy.fft import fft
 
+
 def TFourier(x, y, deltaX):
     """
     Compute the Fourier transform of a function y(x) with sampling interval deltaX.

sasmodels/TwoYukawa/TInvFourier.py

@@ -1,6 +1,7 @@
-from numpy import exp, ceil, log2, pi, arange, interp
+from numpy import arange, ceil, exp, interp, log2, pi
 from numpy.fft import fft
 
+
 def TInvFourier(x, y, deltaX):
     """
     Inverse Fourier transform implementation

sasmodels/bumps_model.py

@@ -16,6 +16,17 @@
 
 import numpy as np  # type: ignore
 
+try:
+    # Optional import. This allows the doc builder and tests to run even
+    # when bumps is not on the path.
+    from bumps.parameter import Parameter as BumpsParameter  # type: ignore
+    from bumps.parameter import Reference
+except ImportError:
+    class BumpsParameter:
+        """See parameter.Parameter in the bumps documentation."""
+    class Reference:
+        """See parameter.Reference in the bumps documentation."""
+
 from .data import plot_theory
 from .direct_model import DataMixin
 
@@ -61,22 +72,34 @@ def create_parameters(model_info, **kwargs):
     parameters for each model parameter, and a dictionary of
     *{name: str}* containing the polydispersity distribution types.
     """
-    pars = {}     # type: dict[str, BumpsParameter]
-    pd_types = {} # type: dict[str, str]
+    def addpar(name: str, default: float, limits: tuple[float]) -> None:
+        value = kwargs.pop(name, default)
+        pars[name] = BumpsParameter.default(value, name=name, limits=limits)
+
+    pars: dict[str, BumpsParameter] = {}
+    pd_types: dict[str, str] = {}
     for p in model_info.parameters.call_parameters:
-        value = kwargs.pop(p.name, p.default)
-        pars[p.name] = BumpsParameter.default(value, name=p.name, limits=p.limits)
+        # To check the limits on parameters in sasmodels, the following
+        # regex captures most (all?) parameter definition lines:
+        #    grep "^\(parameters *= *\[\)\? *\[['\"]" sasmodels/models/*.py
+        # Skip those which have inf, 360 or 180 in the definition:
+        #    ... | grep -v inf | grep -v 360 | grep -v 180
+        # Looking at the remaining values they all represent hard limits on
+        # the model parameter that should not be exceeded by the fitter.
+
+        # Bumps respects hard limits, so make sure they are big enough. We probably don't
+        # need to do this since the values in the models should be set right.
+        #   angular = (p.type == "orientation") or (p.type == "magnetic" and p.units == "degrees")
+        #   limits = p.limits if not angular else (-180, 180) if "theta" in p.name else (-360, 360)
+        addpar(p.name, p.default, p.limits)
         if p.polydisperse:
-            pd_limits = (-360.0, 360.0) if p.type == "orientation" else (0., 1.)
-            for part, default, limits in [
-                    ('_pd', 0., pd_limits),
-                    ('_pd_n', 35., (0, 1000)),
-                    ('_pd_nsigma', 3., (0, 10)),
-                ]:
-                name = p.name + part
-                value = kwargs.pop(name, default)
-                pars[name] = BumpsParameter.default(value, name=name, limits=limits)
+            addpar(f"{p.name}_pd", 0., limits=(0.0, np.inf))
+            addpar(f"{p.name}_pd_n", 35, limits=(0, 1000))
+            addpar(f"{p.name}_pd_nsigma", 3., limits=(0.0, 10.0))
             name = p.name + '_pd_type'
+            # If angular we should be defaulting to a cyclic gaussian, but this is not
+            # yet in sasmodels. There is an implementation in the example/weights directory.
+            # When we do add it be sure to update orientation.rst describing how to use it.
             pd_types[name] = str(kwargs.pop(name, 'gaussian'))
 
     if kwargs:  # args not corresponding to parameters

sasmodels/modelinfo.py

@@ -641,19 +641,19 @@ def _get_call_parameters(self):
                           'magnetic', 'fraction of spin up incident'),
                 Parameter('up_frac_f', '', 0., [0., 1.],
                           'magnetic', 'fraction of spin up final'),
-                Parameter('up_theta', 'degrees', 90., [0., 360.],
+                Parameter('up_theta', 'degrees', 90., [-180., 180.],
                           'magnetic', 'polarization axis rotation angle'),
-                Parameter('up_phi', 'degrees', 0., [0., 180.],
+                Parameter('up_phi', 'degrees', 0., [-360., 360.],
                           'magnetic', 'polarization axis inclination angle'),
             ])
             slds = [p for p in full_list if p.type == 'sld']
             for p in slds:
                 full_list.extend([
                     Parameter(p.id+'_M0', '1e-6/Ang^2', 0., [-np.inf, np.inf],
                               'magnetic', 'magnetic amplitude for '+p.description),
-                    Parameter(p.id+'_mtheta', 'degrees', 0., [-90., 90.],
+                    Parameter(p.id+'_mtheta', 'degrees', 0., [-180., 180.],
                               'magnetic', 'magnetic latitude for '+p.description),
-                    Parameter(p.id+'_mphi', 'degrees', 0., [-180., 180.],
+                    Parameter(p.id+'_mphi', 'degrees', 0., [-360., 360.],
                               'magnetic', 'magnetic longitude for '+p.description),
                 ])
         #print("call parameters", full_list)

sasmodels/models/barbell.py

@@ -112,7 +112,7 @@
     ["radius_bell", "Ang",         40, [0, inf],    "volume",      "Spherical bell radius"],
     ["radius",      "Ang",         20, [0, inf],    "volume",      "Cylindrical bar radius"],
     ["length",      "Ang",        400, [0, inf],    "volume",      "Cylinder bar length"],
-    ["theta",       "degrees",     60, [-360, 360], "orientation", "Barbell axis to beam angle"],
+    ["theta",       "degrees",     60, [-180, 180], "orientation", "Barbell axis to beam angle"],
     ["phi",         "degrees",     60, [-360, 360], "orientation", "Rotation about beam"],
     ]
 # pylint: enable=bad-whitespace, line-too-long

sasmodels/models/bcc_paracrystal.py

@@ -188,7 +188,7 @@
               ["radius",      "Ang",        40,    [0, inf],    "volume",      "Particle radius"],
               ["sld",         "1e-6/Ang^2",  4,    [-inf, inf], "sld",         "Particle scattering length density"],
               ["sld_solvent", "1e-6/Ang^2",  1,    [-inf, inf], "sld",         "Solvent scattering length density"],
-              ["theta",       "degrees",    60,    [-360, 360], "orientation", "c axis to beam angle"],
+              ["theta",       "degrees",    60,    [-180, 180], "orientation", "c axis to beam angle"],
               ["phi",         "degrees",    60,    [-360, 360], "orientation", "rotation about beam"],
               ["psi",         "degrees",    60,    [-360, 360], "orientation", "rotation about c axis"]
              ]

sasmodels/models/capped_cylinder.py

@@ -132,7 +132,7 @@
               # both models, one would be a pill.
               ["radius_cap", "Ang",     20, [0, inf],    "volume", "Cap radius"],
               ["length",     "Ang",    400, [0, inf],    "volume", "Cylinder length"],
-              ["theta",      "degrees", 60, [-360, 360], "orientation", "cylinder axis to beam angle"],
+              ["theta",      "degrees", 60, [-180, 180], "orientation", "cylinder axis to beam angle"],
               ["phi",        "degrees", 60, [-360, 360], "orientation", "rotation about beam"],
              ]
 # pylint: enable=bad-whitespace, line-too-long

sasmodels/models/core_shell_bicelle.py

@@ -146,7 +146,7 @@
     ["sld_face",       "1e-6/Ang^2", 4, [-inf, inf], "sld",         "Cylinder face scattering length density"],
     ["sld_rim",        "1e-6/Ang^2", 4, [-inf, inf], "sld",         "Cylinder rim scattering length density"],
     ["sld_solvent",    "1e-6/Ang^2", 1, [-inf, inf], "sld",         "Solvent scattering length density"],
-    ["theta",          "degrees",   90, [-360, 360], "orientation", "cylinder axis to beam angle"],
+    ["theta",          "degrees",   90, [-180, 180], "orientation", "cylinder axis to beam angle"],
     ["phi",            "degrees",    0, [-360, 360], "orientation", "rotation about beam"]
     ]
 

sasmodels/models/core_shell_bicelle_elliptical.py

@@ -137,7 +137,7 @@
     ["sld_face",       "1e-6/Ang^2", 7, [-inf, inf], "sld",         "Cylinder face scattering length density"],
     ["sld_rim",        "1e-6/Ang^2", 1, [-inf, inf], "sld",         "Cylinder rim scattering length density"],
     ["sld_solvent",    "1e-6/Ang^2", 6, [-inf, inf], "sld",         "Solvent scattering length density"],
-    ["theta",       "degrees",    90.0, [-360, 360], "orientation", "Cylinder axis to beam angle"],
+    ["theta",       "degrees",    90.0, [-180, 180], "orientation", "Cylinder axis to beam angle"],
     ["phi",         "degrees",    0,    [-360, 360], "orientation", "Rotation about beam"],
     ["psi",         "degrees",    0,    [-360, 360], "orientation", "Rotation about cylinder axis"]
     ]

sasmodels/models/core_shell_bicelle_elliptical_belt_rough.py

@@ -150,7 +150,7 @@
     ["sld_rim",        "1e-6/Ang^2", 1, [-inf, inf], "sld",         "Cylinder rim scattering length density"],
     ["sld_solvent",    "1e-6/Ang^2", 6, [-inf, inf], "sld",         "Solvent scattering length density"],
     ["sigma",       "Ang",        0,    [0, inf],    "",            "Interfacial roughness"],
-    ["theta",       "degrees",    90.0, [-360, 360], "orientation", "Cylinder axis to beam angle"],
+    ["theta",       "degrees",    90.0, [-180, 180], "orientation", "Cylinder axis to beam angle"],
     ["phi",         "degrees",    0,    [-360, 360], "orientation", "Rotation about beam"],
     ["psi",         "degrees",    0,    [-360, 360], "orientation", "Rotation about cylinder axis"],
     ]

sasmodels/models/core_shell_cylinder.py

@@ -127,7 +127,7 @@
                "Cylinder shell thickness"],
               ["length", "Ang", 400, [0, inf], "volume",
                "Cylinder length"],
-              ["theta", "degrees", 60, [-360, 360], "orientation",
+              ["theta", "degrees", 60, [-180, 180], "orientation",
                "cylinder axis to beam angle"],
               ["phi", "degrees", 60, [-360, 360], "orientation",
                "rotation about beam"],

sasmodels/models/core_shell_ellipsoid.py

@@ -144,7 +144,7 @@
     ["sld_core",      "1e-6/Ang^2", 2,   [-inf, inf], "sld",         "Core scattering length density"],
     ["sld_shell",     "1e-6/Ang^2", 1,   [-inf, inf], "sld",         "Shell scattering length density"],
     ["sld_solvent",   "1e-6/Ang^2", 6.3, [-inf, inf], "sld",         "Solvent scattering length density"],
-    ["theta",         "degrees",    0,   [-360, 360], "orientation", "elipsoid axis to beam angle"],
+    ["theta",         "degrees",    0,   [-180, 180], "orientation", "elipsoid axis to beam angle"],
     ["phi",           "degrees",    0,   [-360, 360], "orientation", "rotation about beam"],
     ]
 # pylint: enable=bad-whitespace, line-too-long

sasmodels/models/core_shell_parallelepiped.py

@@ -216,7 +216,7 @@
                "Thickness of B rim"],
               ["thick_rim_c", "Ang", 10, [0, inf], "volume",
                "Thickness of C rim"],
-              ["theta", "degrees", 0, [-360, 360], "orientation",
+              ["theta", "degrees", 0, [-180, 180], "orientation",
                "c axis to beam angle"],
               ["phi", "degrees", 0, [-360, 360], "orientation",
                "rotation about beam"],

sasmodels/models/cylinder.py

@@ -137,7 +137,7 @@
     ["sld_solvent", "1e-6/Ang^2",   1, [-inf, inf], "sld",         "Solvent scattering length density"],
     ["radius",      "Ang",         20, [0, inf],    "volume",      "Cylinder radius"],
     ["length",      "Ang",        400, [0, inf],    "volume",      "Cylinder length"],
-    ["theta",       "degrees",     60, [-360, 360], "orientation", "cylinder axis to beam angle"],
+    ["theta",       "degrees",     60, [-180, 180], "orientation", "cylinder axis to beam angle"],
     ["phi",         "degrees",     60, [-360, 360], "orientation", "rotation about beam"],
     ]
 # pylint: enable=bad-whitespace, line-too-long

sasmodels/models/ellipsoid.py

@@ -147,7 +147,7 @@
     ["sld_solvent", "1e-6/Ang^2",   1, [-inf, inf], "sld",         "Solvent scattering length density"],
     ["radius_polar",       "Ang",  20, [0, inf],    "volume",      "Polar radius"],
     ["radius_equatorial",  "Ang", 400, [0, inf],    "volume",      "Equatorial radius"],
-    ["theta",          "degrees",  60, [-360, 360], "orientation", "ellipsoid axis to beam angle"],
+    ["theta",          "degrees",  60, [-180, 180], "orientation", "ellipsoid axis to beam angle"],
     ["phi",            "degrees",  60, [-360, 360], "orientation", "rotation about beam"],
     ]
 # pylint: enable=bad-whitespace, line-too-long

sasmodels/models/elliptical_cylinder.py

@@ -117,7 +117,7 @@
               ["length",      "Ang",        400.0, [1, inf],    "volume",      "Length of the cylinder"],
               ["sld",         "1e-6/Ang^2", 4.0,   [-inf, inf], "sld",         "Cylinder scattering length density"],
               ["sld_solvent", "1e-6/Ang^2", 1.0,   [-inf, inf], "sld",         "Solvent scattering length density"],
-              ["theta",       "degrees",    90.0,  [-360, 360], "orientation", "cylinder axis to beam angle"],
+              ["theta",       "degrees",    90.0,  [-180, 180], "orientation", "cylinder axis to beam angle"],
               ["phi",         "degrees",    0,     [-360, 360], "orientation", "rotation about beam"],
               ["psi",         "degrees",    0,     [-360, 360], "orientation", "rotation about cylinder axis"]]
 

sasmodels/models/fcc_paracrystal.py

@@ -186,7 +186,7 @@
               ["radius", "Ang", 40, [0, inf], "volume", "Particle radius"],
               ["sld", "1e-6/Ang^2", 4, [-inf, inf], "sld", "Particle scattering length density"],
               ["sld_solvent", "1e-6/Ang^2", 1, [-inf, inf], "sld", "Solvent scattering length density"],
-              ["theta",       "degrees",    60,    [-360, 360], "orientation", "c axis to beam angle"],
+              ["theta",       "degrees",    60,    [-180, 180], "orientation", "c axis to beam angle"],
               ["phi",         "degrees",    60,    [-360, 360], "orientation", "rotation about beam"],
               ["psi",         "degrees",    60,    [-360, 360], "orientation", "rotation about c axis"]
              ]

sasmodels/models/hollow_cylinder.py

@@ -93,7 +93,7 @@
     ["length",      "Ang",    400.0, [0, inf],    "volume",      "Cylinder total length"],
     ["sld",         "1e-6/Ang^2",  6.3, [-inf, inf], "sld",      "Cylinder sld"],
     ["sld_solvent", "1e-6/Ang^2",  1,   [-inf, inf], "sld",      "Solvent sld"],
-    ["theta",       "degrees", 90,   [-360, 360], "orientation", "Cylinder axis to beam angle"],
+    ["theta",       "degrees", 90,   [-180, 180], "orientation", "Cylinder axis to beam angle"],
     ["phi",         "degrees",  0,   [-360, 360], "orientation", "Rotation about beam"],
     ]
 # pylint: enable=bad-whitespace, line-too-long

sasmodels/models/hollow_rectangular_prism.py

@@ -138,7 +138,7 @@
                "Ratio sides c/a"],
               ["thickness", "Ang", 1, [0, inf], "volume",
                "Thickness of parallelepiped"],
-              ["theta", "degrees", 0, [-360, 360], "orientation",
+              ["theta", "degrees", 0, [-180, 180], "orientation",
                "c axis to beam angle"],
               ["phi", "degrees", 0, [-360, 360], "orientation",
                "rotation about beam"],

sasmodels/models/parallelepiped.py

@@ -223,7 +223,7 @@
                "Second side of the parallelepiped"],
               ["length_c", "Ang", 400, [0, inf], "volume",
                "Larger side of the parallelepiped"],
-              ["theta", "degrees", 60, [-360, 360], "orientation",
+              ["theta", "degrees", 60, [-180, 180], "orientation",
                "c axis to beam angle"],
               ["phi", "degrees", 60, [-360, 360], "orientation",
                "rotation about beam"],

sasmodels/models/rectangular_prism.py

@@ -136,7 +136,7 @@
                "Ratio sides b/a"],
               ["c2a_ratio", "", 1, [0, inf], "volume",
                "Ratio sides c/a"],
-              ["theta", "degrees", 0, [-360, 360], "orientation",
+              ["theta", "degrees", 0, [-180, 180], "orientation",
                "c axis to beam angle"],
               ["phi", "degrees", 0, [-360, 360], "orientation",
                "rotation about beam"],

sasmodels/models/sc_paracrystal.py

@@ -183,7 +183,7 @@
               ["radius",      "Ang",        40.0, [0.0, inf],  "volume",      "Radius of sphere"],
               ["sld",  "1e-6/Ang^2",         3.0, [0.0, inf],  "sld",         "Sphere scattering length density"],
               ["sld_solvent", "1e-6/Ang^2",  6.3, [0.0, inf],  "sld",         "Solvent scattering length density"],
-              ["theta",       "degrees",    0,    [-360, 360], "orientation", "c axis to beam angle"],
+              ["theta",       "degrees",    0,    [-180, 180], "orientation", "c axis to beam angle"],
               ["phi",         "degrees",    0,    [-360, 360], "orientation", "rotation about beam"],
               ["psi",         "degrees",    0,    [-360, 360], "orientation", "rotation about c axis"]
              ]