← Previous Change
Next Change →

Changeset d0b0f5d in sasmodels for doc

Timestamp:

Mar 31, 2019 10:20:32 AM (6 years ago)

Author:

GitHub <noreply@…>

Branches:

master

Children:

d827c5e

Parents:

7050455 (diff), e2da671 (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:

Andrew Jackson <andrew.jackson@…> (03/31/19 10:20:32)

git-committer:

GitHub <noreply@…> (03/31/19 10:20:32)

Message:

Merge pull request #67 from SasView?/ticket-608-user-defined-weights

Ticket 608 user defined weights

Location:

doc/guide

Files:

: 1 added
: 4 edited

pd/polydispersity.rst (modified) (11 diffs)
fitting_sq.rst (added)
index.rst (modified) (1 diff)
plugin.rst (modified) (15 diffs)
resolution.rst (modified) (3 diffs)

Legend:

: Unmodified
: Added
: Removed

doc/guide/pd/polydispersity.rst

-                      rd089a00
+                      ra5cb9bc
 --------------------------------------------
 For some models we can calculate the average intensity for a population of
 particles that possess size and/or orientational (ie, angular) distributions.
 In SasView we call the former *polydispersity* but use the parameter *PD* to
 parameterise both. In other words, the meaning of *PD* in a model depends on
+For some models we can calculate the average intensity for a population of
+particles that possess size and/or orientational (ie, angular) distributions.
+In SasView we call the former *polydispersity* but use the parameter *PD* to
+parameterise both. In other words, the meaning of *PD* in a model depends on
 the actual parameter it is being applied too.
 The resultant intensity is then normalized by the average particle volume such
+The resultant intensity is then normalized by the average particle volume such
 that
 …
   P(q) = \text{scale} \langle F^* F \rangle / V + \text{background}
 where $F$ is the scattering amplitude and $\langle\cdot\rangle$ denotes an
+where $F$ is the scattering amplitude and $\langle\cdot\rangle$ denotes an
 average over the distribution $f(x; \bar x, \sigma)$, giving
 .. math::
   P(q) = \frac{\text{scale}}{V} \int_\mathbb{R}
+  P(q) = \frac{\text{scale}}{V} \int_\mathbb{R}
   f(x; \bar x, \sigma) F^2(q, x)\, dx + \text{background}
 Each distribution is characterized by a center value $\bar x$ or
 $x_\text{med}$, a width parameter $\sigma$ (note this is *not necessarily*
+the standard deviation, so read the description of the distribution carefully),
+the number of sigmas $N_\sigma$ to include from the tails of the distribution,
+and the number of points used to compute the average. The center of the
+distribution is set by the value of the model parameter.
+The distribution width applied to *volume* (ie, shape-describing) parameters
+is relative to the center value such that $\sigma = \mathrm{PD} \cdot \bar x$.
+However, the distribution width applied to *orientation* parameters is just
+$\sigma = \mathrm{PD}$.
+the standard deviation, so read the description carefully), the number of
+sigmas $N_\sigma$ to include from the tails of the distribution, and the
+number of points used to compute the average. The center of the distribution
+is set by the value of the model parameter. The meaning of a polydispersity
+parameter *PD* (not to be confused with a molecular weight distributions
+in polymer science) in a model depends on the type of parameter it is being
+applied too.
+The distribution width applied to *volume* (ie, shape-describing) parameters
+is relative to the center value such that $\sigma = \mathrm{PD} \cdot \bar x$.
+However, the distribution width applied to *orientation* (ie, angle-describing)
+parameters is just $\sigma = \mathrm{PD}$.
 $N_\sigma$ determines how far into the tails to evaluate the distribution,
 …
 Users should note that the averaging computation is very intensive. Applying
 polydispersion and/or orientational distributions to multiple parameters at
 the same time, or increasing the number of points in the distribution, will
 require patience! However, the calculations are generally more robust with
+polydispersion and/or orientational distributions to multiple parameters at
+the same time, or increasing the number of points in the distribution, will
+require patience! However, the calculations are generally more robust with
 more data points or more angles.
 …
 *  *Schulz Distribution*
 *  *Array Distribution*
+*  *User-defined Distributions*
 These are all implemented as *number-average* distributions.
-Additional distributions are under consideration.
 **Beware: when the Polydispersity & Orientational Distribution panel in SasView is**
 …
 **This may not be suitable. See Suggested Applications below.**
 .. note:: In 2009 IUPAC decided to introduce the new term 'dispersity' to replace
            the term 'polydispersity' (see `Pure Appl. Chem., (2009), 81(2),
 -353 <http://media.iupac.org/publications/pac/2009/pdf/8102x0351.pdf>`_
            in order to make the terminology describing distributions of chemical
            properties unambiguous. However, these terms are unrelated to the
            proportional size distributions and orientational distributions used in
+.. note:: In 2009 IUPAC decided to introduce the new term 'dispersity' to replace
+           the term 'polydispersity' (see `Pure Appl. Chem., (2009), 81(2),
+-353 <http://media.iupac.org/publications/pac/2009/pdf/8102x0351.pdf>`_
+           in order to make the terminology describing distributions of chemical
+           properties unambiguous. However, these terms are unrelated to the
+           proportional size distributions and orientational distributions used in
            SasView models.
 …
 or angular orientations, consider using the Gaussian or Boltzmann distributions.
 If applying polydispersion to parameters describing angles, use the Uniform
 distribution. Beware of using distributions that are always positive (eg, the
+If applying polydispersion to parameters describing angles, use the Uniform
+distribution. Beware of using distributions that are always positive (eg, the
 Lognormal) because angles can be negative!
+The array distribution allows a user-defined distribution to be applied.
+The array distribution provides a very simple means of implementing a user-
+defined distribution, but without any fittable parameters. Greater flexibility
+is conferred by the user-defined distribution.
 .. ZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZ
 …
 .. ZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZ
+User-defined Distributions
+^^^^^^^^^^^^^^^^^^^^^^^^^^
+You can also define your own distribution by creating a python file defining a
+*Distribution* object with a *_weights* method.  The *_weights* method takes
+*center*, *sigma*, *lb* and *ub* as arguments, and can access *self.npts*
+and *self.nsigmas* from the distribution.  They are interpreted as follows:
+* *center* the value of the shape parameter (for size dispersity) or zero
+  if it is an angular dispersity.  This parameter may be fitted.
+* *sigma* the width of the distribution, which is the polydispersity parameter
+  times the center for size dispersity, or the polydispersity parameter alone
+  for angular dispersity.  This parameter may be fitted.
+* *lb*, *ub* are the parameter limits (lower & upper bounds) given in the model
+  definition file.  For example, a radius parameter has *lb* equal to zero.  A
+  volume fraction parameter would have *lb* equal to zero and *ub* equal to one.
+* *self.nsigmas* the distance to go into the tails when evaluating the
+  distribution.  For a two parameter distribution, this value could be
+  co-opted to use for the second parameter, though it will not be available
+  for fitting.
+* *self.npts* the number of points to use when evaluating the distribution.
+  The user will adjust this to trade calculation time for accuracy, but the
+  distribution code is free to return more or fewer, or use it for the third
+  parameter in a three parameter distribution.
+As an example, the code following wraps the Laplace distribution from scipy stats::
+    import numpy as np
+    from scipy.stats import laplace
+    from sasmodels import weights
+    class Dispersion(weights.Dispersion):
+        r"""
+        Laplace distribution
+        .. math::
+            w(x) = e^{-\sigma |x - \mu|}
+        """
+        type = "laplace"
+        default = dict(npts=35, width=0, nsigmas=3)  # default values
+        def _weights(self, center, sigma, lb, ub):
+            x = self._linspace(center, sigma, lb, ub)
+            wx = laplace.pdf(x, center, sigma)
+            return x, wx
+You can plot the weights for a given value and width using the following::
+    from numpy import inf
+    from matplotlib import pyplot as plt
+    from sasmodels import weights
+    # reload the user-defined weights
+    weights.load_weights()
+    x, wx = weights.get_weights('laplace', n=35, width=0.1, nsigmas=3, value=50,
+                                limits=[0, inf], relative=True)
+    # plot the weights
+    plt.interactive(True)
+    plt.plot(x, wx, 'x')
+The *self.nsigmas* and *self.npts* parameters are normally used to control
+the accuracy of the distribution integral. The *self._linspace* function
+uses them to define the *x* values (along with the *center*, *sigma*,
+*lb*, and *ub* which are passed as parameters).  If you repurpose npts or
+nsigmas you will need to generate your own *x*.  Be sure to honour the
+limits *lb* and *ub*, for example to disallow a negative radius or constrain
+the volume fraction to lie between zero and one.
+To activate a user-defined distribution, put it in a file such as *distname.py*
+in the *SAS_WEIGHTS_PATH* folder.  This is defined with an environment
+variable, defaulting to::
+    SAS_WEIGHTS_PATH=~/.sasview/weights
+The weights path is loaded on startup.  To update the distribution definition
+in a running application you will need to enter the following python commands::
+    import sasmodels.weights
+    sasmodels.weights.load_weights('path/to/distname.py')
+.. ZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZ
 Note about DLS polydispersity
 ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
 Several measures of polydispersity abound in Dynamic Light Scattering (DLS) and
 it should not be assumed that any of the following can be simply equated with
+Several measures of polydispersity abound in Dynamic Light Scattering (DLS) and
+it should not be assumed that any of the following can be simply equated with
 the polydispersity *PD* parameter used in SasView.
 The dimensionless **Polydispersity Index (PI)** is a measure of the width of the
 distribution of autocorrelation function decay rates (*not* the distribution of
 particle sizes itself, though the two are inversely related) and is defined by
+The dimensionless **Polydispersity Index (PI)** is a measure of the width of the
+distribution of autocorrelation function decay rates (*not* the distribution of
+particle sizes itself, though the two are inversely related) and is defined by
 ISO 22412:2017 as
 …
     PI = \mu_{2} / \bar \Gamma^2
 where $\mu_\text{2}$ is the second cumulant, and $\bar \Gamma^2$ is the
+where $\mu_\text{2}$ is the second cumulant, and $\bar \Gamma^2$ is the
 intensity-weighted average value, of the distribution of decay rates.
 …
     PI = \sigma^2 / 2\bar \Gamma^2
 where $\sigma$ is the standard deviation, allowing a **Relative Polydispersity (RP)**
+where $\sigma$ is the standard deviation, allowing a **Relative Polydispersity (RP)**
 to be defined as
 …
     RP = \sigma / \bar \Gamma = \sqrt{2 \cdot PI}
 PI values smaller than 0.05 indicate a highly monodisperse system. Values
+PI values smaller than 0.05 indicate a highly monodisperse system. Values
 greater than 0.7 indicate significant polydispersity.
 The **size polydispersity P-parameter** is defined as the relative standard
 deviation coefficient of variation
+The **size polydispersity P-parameter** is defined as the relative standard
+deviation coefficient of variation
 .. math::
 …
 where $\nu$ is the variance of the distribution and $\bar R$ is the mean
 value of $R$. Here, the product $P \bar R$ is *equal* to the standard
+value of $R$. Here, the product $P \bar R$ is *equal* to the standard
 deviation of the Lognormal distribution.

doc/guide/index.rst

-                      rda5536f
+                      rbc69321
    pd/polydispersity.rst
    resolution.rst
+   plugin.rst
+   fitting_sq.rst
    magnetism/magnetism.rst
    orientation/orientation.rst
    sesans/sans_to_sesans.rst
    sesans/sesans_fitting.rst
-   plugin.rst
    scripting.rst
    refs.rst

doc/guide/plugin.rst

-                      r9150036
+                      re15a822
 **Note: The order of the parameters in the definition will be the order of the
 parameters in the user interface and the order of the parameters in Fq(), Iq(),
 Iqac(), Iqabc(), form_volume() and shell_volume().
+Iqac(), Iqabc(), radius_effective(), form_volume() and shell_volume().
 And** *scale* **and** *background* **parameters are implicit to all models,
 so they do not need to be included in the parameter table.**
 …
 can take arbitrary values, even for integer parameters, so your model should
 round the incoming parameter value to the nearest integer inside your model
+you should round to the nearest integer.  In C code, you can do this using::
+you should round to the nearest integer.  In C code, you can do this using:
+.. code-block:: c
     static double
 …
 .............
+.. note::
+   Pure python models do not yet support direct computation of $<F(Q)^2>$ or
+   $<F(Q)>^2$. Neither do they support orientational distributions or magnetism
+   (use C models if these are required).
 For pure python models, define the *Iq* function::
 …
 Models should define *form_volume(par1, par2, ...)* where the parameter
 list includes the *volume* parameters in order.  This is used for a weighted
+volume normalization so that scattering is on an absolute scale.  If
+*form_volume* is not defined, then the default *form_volume = 1.0* will be
+used.
+volume normalization so that scattering is on an absolute scale.  For
+solid shapes, the *I(q)* function should use *form_volume* squared
+as its scale factor.  If *form_volume* is not defined, then the default
+*form_volume = 1.0* will be used.
 Hollow shapes, where the volume fraction of particle corresponds to the
 material in the shell rather than the volume enclosed by the shape, must
 also define a *shell_volume(par1, par2, ...)* function.  The parameters
+are the same as for *form_volume*.  The *I(q)* calculation should use
+*shell_volume* squared as its scale factor for the volume normalization.
+The structure factor calculation needs *form_volume* in order to properly
+scale the volume fraction parameter, so both functions are required for
+hollow shapes.
+Note: Pure python models do not yet support direct computation of the
+average of $F(q)$ and $F^2(q)$.
+are the same as for *form_volume*.  Here the *I(q)* function should use
+*shell_volume* squared instead of *form_volume* squared so that the scale
+parameter corresponds to the volume fraction of material in the sample.
+The structure factor calculation needs the volume fraction of the filled
+shapes for its calculation, so the volume fraction parameter in the model
+is automatically scaled by *form_volume/shell_volume* prior to calling the
+structure factor.
 Embedded C Models
 …
     """
+This expands into the equivalent C code::
+This expands into the equivalent C code:
+.. code-block:: c
     double Iq(double q, double par1, double par2, ...);
 …
 includes only the volume parameters.
 *form_volume* defines the volume of the shell for hollow shapes. As in
+*shell_volume* defines the volume of the shell for hollow shapes. As in
 python models, it includes only the volume parameters.
 …
 Rather than returning NAN from Iq, you must define the *INVALID(v)*.  The
 *v* parameter lets you access all the parameters in the model using
+*v.par1*, *v.par2*, etc. For example::
+*v.par1*, *v.par2*, etc. For example:
+.. code-block:: c
     #define INVALID(v) (v.bell_radius < v.radius)
 …
 Instead of defining the *Iq* function, models can define *Fq* as
+something like::
+something like:
+.. code-block:: c
     double Fq(double q, double *F1, double *F2, double par1, double par2, ...);
 …
 laboratory frame and beam travelling along $-z$.
+The oriented C model is called using *Iqabc(qa, qb, qc, par1, par2, ...)* where
+The oriented C model (oriented pure Python models are not supported)
+is called using *Iqabc(qa, qb, qc, par1, par2, ...)* where
 *par1*, etc. are the parameters to the model.  If the shape is rotationally
 symmetric about *c* then *psi* is not needed, and the model is called
 …
 Using the $z, w$ values for Gauss-Legendre integration in "lib/gauss76.c", the
+numerical integration is then::
+numerical integration is then:
+.. code-block:: c
     double outer_sum = 0.0;
 …
 to compute the proper magnetism and orientation, which you can implement
 using *Iqxy(qx, qy, par1, par2, ...)*.
+**Note: Magnetism is not supported in pure Python models.**
 Special Functions
 …
   memory, and wrong answers computed. The conclusion from a very long and
   strange debugging session was that any arrays that you declare in your
+  model should be a multiple of four. For example::
+  model should be a multiple of four. For example:
+  .. code-block:: c
       double Iq(q, p1, p2, ...)
 …
 structure factor is the *hardsphere* interaction, which
 uses the effective radius of the form factor as an input to the structure
+factor model.  The effective radius is the average radius of the
+form averaged over all the polydispersity values.
+::
+    def ER(radius, thickness):
+        """Effective radius of a core-shell sphere."""
+        return radius + thickness
+Now consider the *core_shell_sphere*, which has a simple effective radius
+equal to the radius of the core plus the thickness of the shell, as
+shown above. Given polydispersity over *(r1, r2, ..., rm)* in radius and
+*(t1, t2, ..., tn)* in thickness, *ER* is called with a mesh
+grid covering all possible combinations of radius and thickness.
+That is, *radius* is *(r1, r2, ..., rm, r1, r2, ..., rm, ...)*
+and *thickness* is *(t1, t1, ... t1, t2, t2, ..., t2, ...)*.
+The *ER* function returns one effective radius for each combination.
+The effective radius calculator weights each of these according to
+the polydispersity distributions and calls the structure factor
+with the average *ER*.
+::
+    def VR(radius, thickness):
+        """Sphere and shell volumes for a core-shell sphere."""
+        whole = 4.0/3.0 * pi * (radius + thickness)**3
+        core = 4.0/3.0 * pi * radius**3
+        return whole, whole - core
+Core-shell type models have an additional volume ratio which scales
+the structure factor.  The *VR* function returns the volume of
+the whole sphere and the volume of the shell. Like *ER*, there is
+one return value for each point in the mesh grid.
+*NOTE: we may be removing or modifying this feature soon. As of the
+time of writing, core-shell sphere returns (1., 1.) for VR, giving a volume
+ratio of 1.0.*
+factor model.  The effective radius is the weighted average over all
+values of the shape in polydisperse systems.
+There can be many notions of effective radius, depending on the shape.  For
+a sphere it is clearly just the radius, but for an ellipsoid of revolution
+we provide average curvature, equivalent sphere radius, minimum radius and
+maximum radius.  These options are listed as *radius_effective_modes* in
+the python model defintion, and must be computed by the *radius_effective*
+function which takes the *radius_effective_mode* parameter as an integer,
+along with the various model parameters.  Unlike normal C/Python arrays,
+the first mode is 1, the second is 2, etc.  Mode 0 indicates that the
+effective radius from the shape is to be ignored in favour of the the
+effective radius parameter in the structure factor model.
+Consider the core-shell sphere, which defines the following effective radius
+modes in the python model::
+    radius_effective_modes = [
+        "outer radius",
+        "core radius",
+    ]
+and the following function in the C-file for the model:
+.. code-block:: c
+    static double
+    radius_effective(int mode, double radius, double thickness)
+    {
+        switch (mode) {
+            case 0: return radius + thickness;
+            case 1: return radius;
+            default: return 0.;
+        }
+    }
+    static double
+    form_volume(double radius, double thickness)
+    {
+        return M_4PI_3 * cube(radius + thickness);
+    }
+Given polydispersity over *(r1, r2, ..., rm)* in radius and *(t1, t2, ..., tn)*
+in thickness, *radius_effective* is called over a mesh grid covering all
+possible combinations of radius and thickness, with a single *(ri, tj)* pair
+in each call. The weights each of these results according to the
+polydispersity distributions and calls the structure factor with the average
+effective radius.  Similarly, for *form_volume*.
+Hollow models have an additional volume ratio which is needed to scale the
+structure factor.  The structure factor uses the volume fraction of the filled
+particles as part of its density estimate, but the scale factor for the
+scattering intensity (as non-solvent volume fraction / volume) is determined
+by the shell volume only.  Therefore the *shell_volume* function is
+needed to compute the form:shell volume ratio, which then scales the
+*volfraction* parameter prior to calling the structure factor calculator.
+In the case of a hollow sphere, this would be:
+.. code-block:: c
+    static double
+    shell_volume(double radius, double thickness)
+    {
+        double whole = M_4PI_3 * cube(radius + thickness);
+        double core = M_4PI_3 * cube(radius);
+        return whole - core;
+    }
+If *shell_volume* is not present, then *form_volume* and *shell_volume* are
+assumed to be equal, and the shape is considered solid.
 Unit Tests
 …
 and a check that the model runs.
-If you are not using sasmodels from SasView, skip this step.
 Recommended Testing
 ...................
+**NB: For now, this more detailed testing is only possible if you have a
+SasView build environment available!**
 If the model compiles and runs, you can next run the unit tests that
 …
 | 2016-10-25 Steve King
 | 2017-05-07 Paul Kienzle - Moved from sasview to sasmodels docs
+| 2019-03-28 Paul Kienzle - Update docs for radius_effective and shell_volume

doc/guide/resolution.rst

-                      r0db85af
+                      rdb1d9d5
 .. sm_help.rst
 .. This is a port of the original SasView html help file to ReSTructured text
 .. by S King, ISIS, during SasView CodeCamp-III in Feb 2015.
+.. resolution.rst
+.. This is a port of the original SasView html help file sm_help to ReSTructured
+.. text by S King, ISIS, during SasView CodeCamp-III in Feb 2015.
 …
 resolution contribution into a model calculation/simulation (which by definition
 will be exact) to make it more representative of what has been measured
+experimentally - a process called *smearing*. Sasmodels does the latter.
+experimentally - a process called *smearing*. The Sasmodels component of SasView
+does the latter.
 Both smearing and desmearing rely on functions to describe the resolution
 …
 for the instrument and stored with the data file.  If not, they will need to
 be set manually before fitting.
+.. note::
+    Problems may be encountered if the data set loaded by SasView is a
+    concatenation of SANS data from several detector distances where, of
+    course, the worst Q resolution is next to the beam stop at each detector
+    distance. (This will also be noticeable in the residuals plot where
+    there will be poor overlap). SasView sensibly orders all the input
+    data points by increasing Q for nicer-looking plots, however, the dQ
+    data can then vary considerably from point to point. If 'Use dQ data'
+    smearing is selected then spikes may appear in the model fits, whereas
+    if 'None' or 'Custom Pinhole Smear' are selected the fits look normal.
+    In such instances, possible solutions are to simply remove the data
+    with poor Q resolution from the shorter detector distances, or to fit
+    the data from different detector distances simultaneously.

Note: See TracChangeset for help on using the changeset viewer.

SasView