[f237c05] | 1 | .. _Writing_a_Plugin: |
---|
| 2 | |
---|
| 3 | Writing a Plugin |
---|
| 4 | ================ |
---|
| 5 | |
---|
| 6 | Users can write their own models and save it to the the SasView |
---|
| 7 | *plugin_models* folder |
---|
| 8 | |
---|
| 9 | *C:\\Users\\[username]\\.sasview\\plugin_models* - (on Windows) |
---|
| 10 | |
---|
| 11 | The next time SasView is started it will compile the plugin and add |
---|
| 12 | it to the list of *Customized Models*. It is recommended that an |
---|
| 13 | existing model be used as a template. |
---|
| 14 | |
---|
| 15 | This page was originally written by our MOST experienced developers, |
---|
| 16 | but has subsequently been edited by our LEAST experienced developer who felt |
---|
| 17 | some instructions could have been clearer, and learnt one or two things that |
---|
| 18 | were missing altogether! But they succeeded in converting a model that passed |
---|
| 19 | testing, so there is no reason why you should not be able to do the same. |
---|
| 20 | |
---|
| 21 | SasView has three ways of writing models: |
---|
| 22 | |
---|
| 23 | - As a pure python model : Example - |
---|
| 24 | `broadpeak.py <https://github.com/SasView/sasmodels/blob/master/sasmodels/models/broad_peak.py>`_ |
---|
| 25 | - As a python model with embedded C : Example - |
---|
| 26 | `sphere.py <https://github.com/SasView/sasmodels/blob/master/sasmodels/models/sphere.py>`_ |
---|
| 27 | - As a python wrapper with separate C code : Example - |
---|
| 28 | `cylinder.py <https://github.com/SasView/sasmodels/blob/master/sasmodels/models/cylinder.py>`_, |
---|
| 29 | `cylinder.c <https://github.com/SasView/sasmodels/blob/master/sasmodels/models/cylinder.c>`_ |
---|
| 30 | |
---|
| 31 | Many models are available for download from the |
---|
| 32 | `model marketplace <http://marketplace.sasview.org/>`_. |
---|
| 33 | |
---|
| 34 | The builtin modules are available in the *sasmodels-data/models* subdirectory |
---|
| 35 | of the sasview distribution. On Windows, this will be something like |
---|
| 36 | *C:\Program Files (x86)\SasView\models*. On MacOSX, these will be within |
---|
| 37 | the application bundle as |
---|
| 38 | */Applications/SasView 4.0.app/Contents/Resources/sasmodels-data/models*. |
---|
| 39 | |
---|
| 40 | Create New Model Files |
---|
| 41 | ^^^^^^^^^^^^^^^^^^^^^^ |
---|
| 42 | |
---|
| 43 | In the *~/.sasview/plugin_models* directory, copy the appropriate files |
---|
| 44 | (using the examples above as templates) to mymodel.py (and mymodel.c, etc) |
---|
| 45 | as required, where "mymodel" is the name for the model you are creating. |
---|
| 46 | |
---|
| 47 | *Please follow these naming rules:* |
---|
| 48 | |
---|
| 49 | - No capitalization and thus no CamelCase. |
---|
| 50 | - If necessary use underscore to separate (i.e. barbell not BarBell or |
---|
| 51 | broad_peak not BroadPeak) |
---|
| 52 | - Don't include âmodelâ in the name (i.e. barbell not BarBellModel) |
---|
| 53 | |
---|
| 54 | |
---|
| 55 | Edit New Model Files |
---|
| 56 | ^^^^^^^^^^^^^^^^^^^^ |
---|
| 57 | |
---|
| 58 | The model interface definition is in the .py file. This file contains: |
---|
| 59 | |
---|
| 60 | - a **model name**: |
---|
| 61 | - this is the **name** string in the *.py* file |
---|
| 62 | - titles should be: |
---|
| 63 | |
---|
| 64 | - all in *lower* case |
---|
| 65 | - without spaces (use underscores to separate words instead) |
---|
| 66 | - without any capitalization or CamelCase |
---|
| 67 | - without incorporating the word 'model' |
---|
| 68 | - examples: *barbell* **not** *BarBell*; *broad_peak* **not** *BroadPeak*; |
---|
| 69 | *barbell* **not** *BarBellModel* |
---|
| 70 | |
---|
| 71 | - a **model title**: |
---|
| 72 | - this is the **title** string in the *.py* file |
---|
| 73 | - this is a one or two line description of the model, which will appear |
---|
| 74 | at the start of the model documentation and as a tooltip in the GUI |
---|
| 75 | |
---|
| 76 | - a **short discription**: |
---|
| 77 | - this is the **description** string in the *.py* file |
---|
| 78 | - this is a medium length description which appears when you click |
---|
| 79 | *Description* on the model fit page |
---|
| 80 | |
---|
| 81 | - a **parameter table**: |
---|
| 82 | - this will be auto-generated from the *parameters* in the *.py* file |
---|
| 83 | |
---|
| 84 | - a **long description**: |
---|
| 85 | - this is ReStructuredText enclosed between the r""" and """ delimiters |
---|
| 86 | at the top of the *.py* file |
---|
| 87 | |
---|
| 88 | - a **definition** of the model: |
---|
| 89 | - as part of the **long description** |
---|
| 90 | |
---|
| 91 | - a **formula** defining the function the model calculates: |
---|
| 92 | - as part of the **long description** |
---|
| 93 | |
---|
| 94 | - an **explanation of the parameters**: |
---|
| 95 | - as part of the **long description** |
---|
| 96 | - explaining how the symbols in the formula map to the model parameters |
---|
| 97 | |
---|
| 98 | - a **plot** of the function, with a **figure caption**: |
---|
| 99 | - this is automatically generated from the default parameters |
---|
| 100 | |
---|
| 101 | - at least one **reference**: |
---|
| 102 | - as part of the **long description** |
---|
| 103 | - specifying where the reader can obtain more information about the model |
---|
| 104 | |
---|
| 105 | - the **name of the author** |
---|
| 106 | - as part of the **long description** |
---|
| 107 | - the *.py* file should also contain a comment identifying *who* |
---|
| 108 | converted/created the model file |
---|
| 109 | |
---|
| 110 | More complete documentation for the sasmodels package can be found at |
---|
| 111 | `<http://www.sasview.org/sasmodels>`_. In particular, |
---|
| 112 | `<http://www.sasview.org/sasmodels/api/generate.html#module-sasmodels.generate>`_ |
---|
| 113 | describes the structure of a model. |
---|
| 114 | |
---|
| 115 | |
---|
| 116 | Model Documentation |
---|
| 117 | ................... |
---|
| 118 | |
---|
| 119 | The *.py* file starts with an r (for raw) and three sets of quotes |
---|
| 120 | to start the doc string and ends with a second set of three quotes. |
---|
| 121 | For example:: |
---|
| 122 | |
---|
| 123 | r""" |
---|
| 124 | Definition |
---|
| 125 | ---------- |
---|
| 126 | |
---|
| 127 | The 1D scattering intensity of the sphere is calculated in the following |
---|
| 128 | way (Guinier, 1955) |
---|
| 129 | |
---|
| 130 | .. math:: |
---|
| 131 | |
---|
| 132 | I(q) = \frac{\text{scale}}{V} \cdot \left[ |
---|
| 133 | 3V(\Delta\rho) \cdot \frac{\sin(qr) - qr\cos(qr))}{(qr)^3} |
---|
| 134 | \right]^2 + \text{background} |
---|
| 135 | |
---|
| 136 | where *scale* is a volume fraction, $V$ is the volume of the scatterer, |
---|
| 137 | $r$ is the radius of the sphere and *background* is the background level. |
---|
| 138 | *sld* and *sld_solvent* are the scattering length densities (SLDs) of the |
---|
| 139 | scatterer and the solvent respectively, whose difference is $\Delta\rho$. |
---|
| 140 | |
---|
| 141 | You can included figures in your documentation, as in the following |
---|
| 142 | figure for the cylinder model. |
---|
| 143 | |
---|
| 144 | .. figure:: img/cylinder_angle_definition.jpg |
---|
| 145 | |
---|
| 146 | Definition of the angles for oriented cylinders. |
---|
| 147 | |
---|
| 148 | References |
---|
| 149 | ---------- |
---|
| 150 | |
---|
| 151 | A Guinier, G Fournet, *Small-Angle Scattering of X-Rays*, |
---|
| 152 | John Wiley and Sons, New York, (1955) |
---|
| 153 | """ |
---|
| 154 | |
---|
| 155 | This is where the FULL documentation for the model goes (to be picked up by |
---|
| 156 | the automatic documentation system). Although it feels odd, you |
---|
| 157 | should start the documentation immediately with the **definition**---the model |
---|
| 158 | name, a brief description and the parameter table are automatically inserted |
---|
| 159 | above the definition, and the a plot of the model is automatically inserted |
---|
| 160 | before the **reference**. |
---|
| 161 | |
---|
| 162 | Figures can be included using the *figure* command, with the name |
---|
| 163 | of the *.png* file containing the figure and a caption to appear below the |
---|
| 164 | figure. Figure numbers will be added automatically. |
---|
| 165 | |
---|
| 166 | See this `Sphinx cheat sheet <http://matplotlib.org/sampledoc/cheatsheet.html>`_ |
---|
| 167 | for a quick guide to the documentation layout commands, or the |
---|
| 168 | `Sphinx Documentation <http://www.sphinx-doc.org/en/stable/>`_ for |
---|
| 169 | complete details. |
---|
| 170 | |
---|
| 171 | The model should include a **formula** written using LaTeX markup. |
---|
| 172 | The above example uses the *math* command to make a displayed equation. You |
---|
| 173 | can also use *\$formula\$* for an inline formula. This is handy for defining |
---|
| 174 | the relationship between the model parameters and formula variables, such |
---|
| 175 | as the phrase "\$r\$ is the radius" used above. The live demo MathJax |
---|
| 176 | page `<http://www.mathjax.org/>`_ is handy for checking that the equations |
---|
| 177 | will look like you expect. |
---|
| 178 | |
---|
| 179 | Math layout uses the `amsmath <http://www.ams.org/publications/authors/tex/amslatex>`_ |
---|
| 180 | package for aligning equations (see amsldoc.pdf on that page for complete documentation). |
---|
| 181 | You will automatically be in an aligned environment, with blank lines separating |
---|
| 182 | the lines of the equation. Place an ampersand before the operator on which to |
---|
| 183 | align. For example:: |
---|
| 184 | |
---|
| 185 | .. math:: |
---|
| 186 | |
---|
| 187 | x + y &= 1 \\ |
---|
| 188 | y &= x - 1 |
---|
| 189 | |
---|
| 190 | produces |
---|
| 191 | |
---|
| 192 | .. math:: |
---|
| 193 | |
---|
| 194 | x + y &= 1 \\ |
---|
| 195 | y &= x - 1 |
---|
| 196 | |
---|
| 197 | If you need more control, use:: |
---|
| 198 | |
---|
| 199 | .. math:: |
---|
| 200 | :nowrap: |
---|
| 201 | |
---|
| 202 | |
---|
| 203 | Model Definition |
---|
| 204 | ................ |
---|
| 205 | |
---|
| 206 | Following the documentation string, there are a series of definitions:: |
---|
| 207 | |
---|
| 208 | name = "sphere" # optional: defaults to the filename without .py |
---|
| 209 | title = "Spheres with uniform scattering length density" |
---|
| 210 | description = """\ |
---|
| 211 | P(q)=(scale/V)*[3V(sld-sld_solvent)*(sin(qr)-qr cos(qr)) |
---|
| 212 | /(qr)^3]^2 + background |
---|
| 213 | r: radius of sphere |
---|
| 214 | V: The volume of the scatter |
---|
| 215 | sld: the SLD of the sphere |
---|
| 216 | sld_solvent: the SLD of the solvent |
---|
| 217 | """ |
---|
| 218 | category = "shape:sphere" |
---|
| 219 | single = True # optional: defaults to True |
---|
| 220 | opencl = False # optional: defaults to False |
---|
| 221 | structure_factor = False # optional: defaults to False |
---|
| 222 | |
---|
| 223 | **name = "mymodel"** defines the name of the model that is shown to the user. |
---|
| 224 | If it is not provided, it will use the name of the model file, with '_' |
---|
| 225 | replaced by spaces and the parts capitalized. So *adsorbed_layer.py* will |
---|
| 226 | become *Adsorbed Layer*. The predefined models all use the name of the |
---|
| 227 | model file as the name of the model, so the default may be changed. |
---|
| 228 | |
---|
| 229 | **title = "short description"** is short description of the model which |
---|
| 230 | is included after the model name in the automatically generated documentation. |
---|
| 231 | The title can also be used for a tooltip, for example. |
---|
| 232 | |
---|
| 233 | **description = """doc string"""** is a longer description of the model. It |
---|
| 234 | shows up when you press the "Description" button of the SasView fit page. |
---|
| 235 | It should give a brief description of the equation and the parameters |
---|
| 236 | without the need to read the entire model documentation. The triple quotes |
---|
| 237 | allow you to write the description over multiple lines. Keep the lines |
---|
| 238 | short since the GUI will wrap each one separately if they are too long. |
---|
| 239 | **Make sure the parameter names in the description match the model definition.** |
---|
| 240 | |
---|
| 241 | **category = "shape:sphere"** defines where the model will appear in the |
---|
| 242 | model documentation. In this example, the model will appear alphabetically |
---|
| 243 | in the list of spheroid models. |
---|
| 244 | |
---|
| 245 | **single = True** indicates that the model can be run using single |
---|
| 246 | precision floating point values. Set it to False if the numerical |
---|
| 247 | calculation for the model is unstable, which is the case for about 20 of |
---|
| 248 | the built in models. It is worthwhile modifying the calculation to support |
---|
| 249 | single precision, allowing models to run up to 10 times faster. The |
---|
| 250 | section `Test_Your_New_Model`_ describes how to compare model values for |
---|
| 251 | single vs. double precision so you can decide if you need to set |
---|
| 252 | single to False. |
---|
| 253 | |
---|
| 254 | **opencl = False** indicates that the model should not be run using OpenCL. |
---|
| 255 | This may be because the model definition includes code that cannot be |
---|
| 256 | compiled for the GPU (for example, goto statements). It can also be used |
---|
| 257 | for large models which can't run on most GPUs. This flag has not been |
---|
| 258 | used on any of the built in models; models which were failing were |
---|
| 259 | streamlined so this flag was not necessary. |
---|
| 260 | |
---|
| 261 | **structure_factor = True** indicates that the model can be used as a |
---|
| 262 | structure factor to account for interactions between particles. See |
---|
| 263 | `Form_Factors`_ for more details. |
---|
| 264 | |
---|
| 265 | Model Parameters |
---|
| 266 | ................ |
---|
| 267 | |
---|
| 268 | Next comes the parameter table. For example:: |
---|
| 269 | |
---|
| 270 | # pylint: disable=bad-whitespace, line-too-long |
---|
| 271 | # ["name", "units", default, [min, max], "type", "description"], |
---|
| 272 | parameters = [ |
---|
| 273 | ["sld", "1e-6/Ang^2", 1, [-inf, inf], "sld", "Layer scattering length density"], |
---|
| 274 | ["sld_solvent", "1e-6/Ang^2", 6, [-inf, inf], "sld", "Solvent scattering length density"], |
---|
| 275 | ["radius", "Ang", 50, [0, inf], "volume", "Sphere radius"], |
---|
| 276 | ] |
---|
| 277 | # pylint: disable=bad-whitespace, line-too-long |
---|
| 278 | |
---|
| 279 | **parameters = [["name", "units", default, [min,max], "type", "tooltip"],...]** |
---|
| 280 | defines the parameters form the model. |
---|
| 281 | |
---|
| 282 | - **the order of the parameters in the definition will be the order of |
---|
| 283 | the parameters in the user interface and the order of the parameters |
---|
| 284 | in Iq(), Iqxy() and form_volume().** |
---|
| 285 | |
---|
| 286 | - **scale and background parameters are implicit to all models, so |
---|
| 287 | they do not need to be included in the parameter table** |
---|
| 288 | |
---|
| 289 | - **"name"** is the name of the parameter shown on the fit screen |
---|
| 290 | |
---|
| 291 | - parameter names should follow the mathematical convention; e.g., |
---|
| 292 | *radius_core* not *core_radius*, or *sld_solvent* not *solvent_sld* |
---|
| 293 | - model parameter names should be consistent between different models, |
---|
| 294 | so *sld_solvent*, for example, should have exactly the same name |
---|
| 295 | in every model |
---|
| 296 | - to see all the parameter names currently in use, type the following in the |
---|
| 297 | python shell/editor under the Tools menu:: |
---|
| 298 | |
---|
| 299 | import sasmodels.list_pars |
---|
| 300 | sasmodels.list_pars.list_pars() |
---|
| 301 | |
---|
| 302 | *re-use* as many as possible!!! |
---|
| 303 | - use "name[n]" for multiplicity parameters, where *n* is the name of |
---|
| 304 | the parameter defining the number of shells/layers/segments, etc. |
---|
| 305 | |
---|
| 306 | - **"units"** are displayed along with the parameter name |
---|
| 307 | |
---|
| 308 | - every parameter should have units; use "None" if there are no units |
---|
| 309 | - **sld's should be given in units of 1e-6/Ang^2, and not simply |
---|
| 310 | 1/Ang^2 to be consistent with the builtin models. Adjust your formulas |
---|
| 311 | appropriately.** |
---|
| 312 | - fancy units markup is available for some units, including:: |
---|
| 313 | |
---|
| 314 | Ang, 1/Ang, 1/Ang^2, 1e-6/Ang^2, degrees, 1/cm, Ang/cm, g/cm^3, mg/m^2 |
---|
| 315 | |
---|
| 316 | - the list of units is defined in the variable *RST_UNITS* within |
---|
| 317 | `sasmodels/generate.py <https://github.com/SasView/sasmodels/tree/master/sasmodels/generate.py>`_ |
---|
| 318 | |
---|
| 319 | - new units can be added using the macros defined in *doc/rst_prolog* |
---|
| 320 | in the sasmodels source. |
---|
| 321 | - units should be properly formatted using sub-/super-scripts |
---|
| 322 | and using negative exponents instead of the / operator, though |
---|
| 323 | the unit name should use the / operator for consistency. |
---|
| 324 | - post a message to the sasview developers list with the changes |
---|
| 325 | |
---|
| 326 | - **default** is the initial value for the parameter |
---|
| 327 | |
---|
| 328 | - **the parameter default values are used to auto-generate a plot of |
---|
| 329 | the model function in the documentation.** |
---|
| 330 | |
---|
| 331 | - **[min, max]** are the lower and upper limits on the parameter |
---|
| 332 | |
---|
| 333 | - lower and upper limits can be any number, or -inf or inf. |
---|
| 334 | - the limits will show up as the default limits for the fit making it easy, |
---|
| 335 | for example, to force the radius to always be greater than zero. |
---|
| 336 | |
---|
| 337 | - **"type"** can be one of: "", "sld", "volume", or "orientation" |
---|
| 338 | |
---|
| 339 | - "sld" parameters can have magnetic moments when fitting magnetic models; |
---|
| 340 | depending on the spin polarization of the beam and the $q$ value being |
---|
| 341 | examined, the effective sld for that material will be used to compute the |
---|
| 342 | scattered intensity |
---|
| 343 | - "volume" parameters are passed to Iq(), Iqxy(), and form_volume(), and |
---|
| 344 | have polydispersity loops generated automatically. |
---|
| 345 | - "orientation" parameters are only passed to Iqxy(), and have angular |
---|
| 346 | dispersion. |
---|
| 347 | |
---|
| 348 | |
---|
| 349 | Model Computation |
---|
| 350 | ................. |
---|
| 351 | |
---|
| 352 | Models can be defined as pure python models, or they can be a mixture of |
---|
| 353 | python and C models. C models are run on the GPU if it is available, |
---|
| 354 | otherwise they are compiled and run on the CPU. |
---|
| 355 | |
---|
| 356 | Models are defined by the scattering kernel, which takes a set of parameter |
---|
| 357 | values defining the shape, orientation and material, and returns the |
---|
| 358 | expected scattering. Polydispersity and angular dispersion are defined |
---|
| 359 | by the computational infrastructure. Any parameters defined as "volume" |
---|
| 360 | parameters are polydisperse, with polydispersity defined in proportion |
---|
| 361 | to their value. "orientation" parameters use angular dispersion defined |
---|
| 362 | in degrees, and are not relative to the current angle. |
---|
| 363 | |
---|
| 364 | Based on a weighting function $G(x)$ and a number of points $n$, the |
---|
| 365 | computed value is |
---|
| 366 | |
---|
| 367 | .. math:: |
---|
| 368 | |
---|
| 369 | \hat I(q) |
---|
| 370 | = \frac{\int G(x) I(q, x)\,dx}{\int G(x) V(x)\,dx} |
---|
| 371 | \approx \frac{\sum_{i=1}^n G(x_i) I(q,x_i)}{\sum_{i=1}^n G(x_i) V(x_i)} |
---|
| 372 | |
---|
| 373 | That is, the indivdual models do not need to include polydispersity |
---|
| 374 | calculations, but instead rely on numerical integration to compute the |
---|
| 375 | appropriately smeared pattern. Angular dispersion values over polar angle |
---|
| 376 | $\theta$ requires an additional $\cos \theta$ weighting due to decreased |
---|
| 377 | arc length for the equatorial angle $\phi$ with increasing latitude. |
---|
| 378 | |
---|
| 379 | Python Models |
---|
| 380 | ............. |
---|
| 381 | |
---|
| 382 | For pure python models, define the Iq funtion:: |
---|
| 383 | |
---|
| 384 | import numpy as np |
---|
| 385 | from numpy import cos, sin, ... |
---|
| 386 | |
---|
| 387 | def Iq(q, par1, par2, ...): |
---|
| 388 | return I(q, par1, par2, ...) |
---|
| 389 | Iq.vectorized = True |
---|
| 390 | |
---|
| 391 | The parameters *par1, par2, ...* are the list of non-orientation parameters |
---|
| 392 | to the model in the order that they appear in the parameter table. |
---|
| 393 | **Note that the autogenerated model file uses *x* rather than *q*.** |
---|
| 394 | |
---|
| 395 | The *.py* file should import trigonometric and exponential functions from |
---|
| 396 | numpy rather that from math. This lets us evaluate the model for the whole |
---|
| 397 | range of $q$ values at once rather than looping over each $q$ separately in |
---|
| 398 | python. With $q$ as a vector, you cannot use if statements, but must instead |
---|
| 399 | do tricks like |
---|
| 400 | |
---|
| 401 | :: |
---|
| 402 | |
---|
| 403 | a = x*q*(q>0) + y*q*(q<=0) |
---|
| 404 | |
---|
| 405 | or |
---|
| 406 | |
---|
| 407 | :: |
---|
| 408 | |
---|
| 409 | a = np.empty_like(q) |
---|
| 410 | index = q>0 |
---|
| 411 | a[index] = x*q[index] |
---|
| 412 | a[~index] = y*q[~index] |
---|
| 413 | |
---|
| 414 | which sets $a$ to $q \cdot x$ if $q$ is positive or $q \cdot y$ if $q$ |
---|
| 415 | is zero or negative. If you have not converted your function to use $q$ |
---|
| 416 | vectors, you can set the following and it will only receive one $q$ |
---|
| 417 | value at a time:: |
---|
| 418 | |
---|
| 419 | Iq.vectorized = False |
---|
| 420 | |
---|
| 421 | Return np.NaN if the parameters are not valid (e.g., cap_radius < radius in |
---|
| 422 | barbell). If I(q; pars) is NaN for any $q$, then those parameters will be |
---|
| 423 | ignored, and not included in the calculation of the weighted polydispersity. |
---|
| 424 | |
---|
| 425 | Similar to *Iq*, you can define *Iqxy(qx, qy, par1, par2, ...)* where the |
---|
| 426 | parameter list includes any orientation parameters. If *Iqxy* is not defined, |
---|
| 427 | then it will default to *Iqxy = Iq(sqrt(qx**2+qy**2), par1, par2, ...)*. |
---|
| 428 | |
---|
| 429 | Models should define *form_volume(par1, par2, ...)* where the parameter |
---|
| 430 | list includes the *volume* parameters in order. This is used for a weighted |
---|
| 431 | volume normalization so that scattering is on an absolute scale. If |
---|
| 432 | *form_volume* is not definded, then the default *form_volume = 1.0* will be |
---|
| 433 | used. |
---|
| 434 | |
---|
| 435 | Embedded C Models |
---|
| 436 | ................. |
---|
| 437 | |
---|
| 438 | Like pure python models, inline C models need define an *Iq* function:: |
---|
| 439 | |
---|
| 440 | Iq = """ |
---|
| 441 | return I(q, par1, par2, ...); |
---|
| 442 | """ |
---|
| 443 | |
---|
| 444 | This expands into the equivalent C code:: |
---|
| 445 | |
---|
| 446 | #include <math.h> |
---|
| 447 | double Iq(double q, double par1, double par2, ...); |
---|
| 448 | double Iq(double q, double par1, double par2, ...) |
---|
| 449 | { |
---|
| 450 | return I(q, par1, par2, ...); |
---|
| 451 | } |
---|
| 452 | |
---|
| 453 | The C model operates on a single $q$ value at a time. The code will be |
---|
| 454 | run in parallel across different $q$ values, either on the graphics card |
---|
| 455 | or the processor. |
---|
| 456 | |
---|
| 457 | Rather than returning NAN from Iq, you must define the *INVALID(v)*. The |
---|
| 458 | *v* parameter lets you access all the parameters in the model using |
---|
| 459 | *v.par1*, *v.par2*, etc. For example:: |
---|
| 460 | |
---|
| 461 | #define INVALID(v) (v.bell_radius < v.radius) |
---|
| 462 | |
---|
| 463 | *Iqxy* is similar to *Iq*, except it uses parameters *qx, qy* instead of *q*, |
---|
| 464 | and it includes orientation parameters. As in python models, *form_volume* |
---|
| 465 | includes only the volume parameters. *Iqxy* will default to |
---|
| 466 | *Iq(sqrt(qx**2 + qy**2), par1, ...)* and *form_volume* will default to 1.0. |
---|
| 467 | |
---|
| 468 | The C code follows the C99 standard, including the usual math functions, |
---|
| 469 | as defined in |
---|
| 470 | `OpenCL <https://www.khronos.org/registry/cl/sdk/1.1/docs/man/xhtml/mathFunctions.html>`_. |
---|
| 471 | |
---|
| 472 | The standard constants and functions include the following:: |
---|
| 473 | |
---|
| 474 | M_PI = pi |
---|
| 475 | M_PI_2 = pi/2 |
---|
| 476 | M_PI_4 = pi/4 |
---|
| 477 | M_E = e |
---|
| 478 | M_SQRT1_2 = 1/sqrt(2) |
---|
| 479 | NAN = NaN |
---|
| 480 | INFINITY = 1/0 |
---|
| 481 | erf(x) = error function |
---|
| 482 | erfc(x) = 1-erf(x) |
---|
| 483 | expm1(x) = exp(x) - 1 |
---|
| 484 | tgamma(x) = gamma function |
---|
| 485 | |
---|
| 486 | Some non-standard constants and functions are also provided:: |
---|
| 487 | |
---|
| 488 | M_PI_180 = pi/180 |
---|
| 489 | M_4PI_3 = 4pi/3 |
---|
| 490 | square(x) = x*x |
---|
| 491 | cube(x) = x*x*x |
---|
| 492 | sinc(x) = sin(x)/x, with sin(0)/0 -> 1 |
---|
| 493 | SINCOS(x, s, c) sets s=sin(angle) and c=cos(angle) |
---|
| 494 | powr(x, y) = x^y for x >= 0 |
---|
| 495 | pown(x, n) = x^n for n integer |
---|
| 496 | |
---|
| 497 | **source=['lib/fn.c', ...]** includes the listed C source files in the |
---|
| 498 | program before *Iq* and *Iqxy* are defined. This allows you to extend the |
---|
| 499 | library of available C functions. Additional special functions and |
---|
| 500 | scattering calculations are defined in |
---|
| 501 | `sasmodels/models/lib <https://github.com/SasView/sasmodels/tree/master/sasmodels/models/lib>`_, |
---|
| 502 | including:: |
---|
| 503 | |
---|
| 504 | sph_j1c(x) = 3 j1(x)/x = 3 (sin(x) - x cos(x))/x^3 [spherical bessel function] |
---|
| 505 | sas_J1c(x) = 2 J1(x)/x [bessel function of the first kind] |
---|
| 506 | sas_gamma(x) = gamma function [tgamma is unstable below 1] |
---|
| 507 | sas_erf(x) = error function [erf is broken on some Intel OpenCL drivers] |
---|
| 508 | sas_erfc(x) = 1-erf(x) |
---|
| 509 | sas_J0(x) = J0(x) |
---|
| 510 | sas_J1(x) = J1(x) |
---|
| 511 | sas_JN(x) = JN(x) |
---|
| 512 | Si(x) = integral sin(z)/z from 0 to x |
---|
| 513 | Gauss76Wt = gaussian quadrature weights for 76 point integral |
---|
| 514 | Gauss76Z = gaussian quadrature values for 76 point integral |
---|
| 515 | |
---|
| 516 | These functions have been tuned to be fast and numerically stable down |
---|
| 517 | to $q=0$ even in single precision. In some cases they work around bugs |
---|
| 518 | which appear on some platforms but not others. |
---|
| 519 | |
---|
| 520 | Models are defined using double precision declarations for the |
---|
| 521 | parameters and return values. Declarations and constants will be converted |
---|
| 522 | to float or long double depending on the precision requested. |
---|
| 523 | **Floating point constants must include the decimal point.** This allows us |
---|
| 524 | to convert values such as 1.0 (double precision) to 1.0f (single precision) |
---|
| 525 | so that expressions that use these values are not promoted to double precision |
---|
| 526 | expressions. Some graphics card drivers are confused when functions |
---|
| 527 | that expect floating point values are passed integers, such as 4*atan(1); it |
---|
| 528 | is safest to not use integers in floating point expressions. Even better, |
---|
| 529 | use the builtin constant M_PI rather than 4*atan(1); it is faster and smaller! |
---|
| 530 | |
---|
| 531 | FLOAT_SIZE is the number of bytes in the converted variables. If your |
---|
| 532 | algorithm depends on precision (which is not uncommon for numerical |
---|
| 533 | algorithms), use the following:: |
---|
| 534 | |
---|
| 535 | #if FLOAT_SIZE>4 |
---|
| 536 | ... code for double precision ... |
---|
| 537 | #else |
---|
| 538 | ... code for single precision ... |
---|
| 539 | #endif |
---|
| 540 | |
---|
| 541 | A value defined as SAS_DOUBLE will stay double precision; this should |
---|
| 542 | not be used since some graphics card don't support double precision. |
---|
| 543 | |
---|
| 544 | |
---|
| 545 | External C Models |
---|
| 546 | ................. |
---|
| 547 | |
---|
| 548 | External C models are very much like embedded C models, except that |
---|
| 549 | *Iq*, *Iqxy* and *form_volume* are defined in an external source file |
---|
| 550 | loaded using the *source=[...]* method. You need to supply the function |
---|
| 551 | declarations for each of these that you need instead of building them |
---|
| 552 | automatically from the parameter table. |
---|
| 553 | |
---|
| 554 | |
---|
| 555 | .. _Form_Factors: |
---|
| 556 | |
---|
| 557 | Form Factors |
---|
| 558 | ............ |
---|
| 559 | |
---|
| 560 | :: |
---|
| 561 | |
---|
| 562 | def ER(radius, thickness): |
---|
| 563 | """Effective radius of a core-shell sphere.""" |
---|
| 564 | return radius + thickness |
---|
| 565 | |
---|
| 566 | Away from the dilute limit you can estimate scattering including |
---|
| 567 | particle-particle interactions using $I(q) = P(q)*S(q)$ where $P(q)$ |
---|
| 568 | is the form factor and $S(q)$ is the structure factor. The simplest |
---|
| 569 | structure factor is the *hardsphere* interaction, which |
---|
| 570 | uses the effective radius of the form factor as an input to the structure |
---|
| 571 | factor model. The effective radius is the average radius of the |
---|
| 572 | form averaged over all the polydispersity values. |
---|
| 573 | |
---|
| 574 | Consider the *core_shell_sphere*, which has a simple effective radius |
---|
| 575 | equal to the radius of the core plus the thickness of the shell, as |
---|
| 576 | shown above. Given polydispersity over *(r1, r2, ..., rm)* in radius and |
---|
| 577 | *(t1, t2, ..., tn)* in thickness, *ER* is called with a mesh |
---|
| 578 | grid covering all possible combinations of radius and thickness. |
---|
| 579 | That is, *radius* is *(r1, r2, ..., rm, r1, r2, ..., rm, ...)* |
---|
| 580 | and *thickness* is *(t1, t1, ... t1, t2, t2, ..., t2, ...)*. |
---|
| 581 | The *ER* function returns one effective radius for each combination. |
---|
| 582 | The effective radius calculator weights each of these according to |
---|
| 583 | the polydispersity distributions and calls the structure factor |
---|
| 584 | with the average *ER*. |
---|
| 585 | |
---|
| 586 | :: |
---|
| 587 | |
---|
| 588 | def VR(radius, thickness): |
---|
| 589 | """Sphere and shell volumes for a core-shell sphere.""" |
---|
| 590 | whole = 4.0/3.0 * pi * (radius + thickness)**3 |
---|
| 591 | core = 4.0/3.0 * pi * radius**3 |
---|
| 592 | return whole, whole - core |
---|
| 593 | |
---|
| 594 | Core-shell type models have an additional volume ratio which scales |
---|
| 595 | the structure factor. The *VR* function returns the volume of |
---|
| 596 | the whole sphere and the volume of the shell. Like *ER*, there is |
---|
| 597 | one return value for each point in the mesh grid. |
---|
| 598 | |
---|
| 599 | *NOTE: we may be removing or modifying this feature soon.* As of this |
---|
| 600 | writing, core-shell sphere returns (1., 1.) for *VR*, giving a volume |
---|
| 601 | ratio of 1.0. |
---|
| 602 | |
---|
| 603 | Unit Tests |
---|
| 604 | .......... |
---|
| 605 | |
---|
| 606 | THESE ARE VERY IMPORTANT. Include at least one test for each model and |
---|
| 607 | PLEASE make sure that the answer value is correct (i.e. not a random number). |
---|
| 608 | |
---|
| 609 | :: |
---|
| 610 | |
---|
| 611 | tests = [ |
---|
| 612 | [{}, 0.2, 0.726362], |
---|
| 613 | [{"scale": 1., "background": 0., "sld": 6., "sld_solvent": 1., |
---|
| 614 | "radius": 120., "radius_pd": 0.2, "radius_pd_n":45}, |
---|
| 615 | 0.2, 0.228843], |
---|
| 616 | [{"radius": 120., "radius_pd": 0.2, "radius_pd_n":45}, "ER", 120.], |
---|
| 617 | [{"radius": 120., "radius_pd": 0.2, "radius_pd_n":45}, "VR", 1.], |
---|
| 618 | ] |
---|
| 619 | |
---|
| 620 | |
---|
| 621 | **tests=[[{parameters}, q, result], ...]** is a list of lists. |
---|
| 622 | Each list is one test and contains, in order: |
---|
| 623 | |
---|
| 624 | - a dictionary of parameter values. This can be {} using the default |
---|
| 625 | parameters, or filled with some parameters that will be different |
---|
| 626 | from the default, such as {âradiusâ:10.0, âsldâ:4}. Unlisted parameters |
---|
| 627 | will be given the default values. |
---|
| 628 | - the input $q$ value or tuple of $(q_x, q_y)$ values. |
---|
| 629 | - the output $I(q)$ or $I(q_x,q_y)$ expected of the model for the parameters |
---|
| 630 | and input value given. |
---|
| 631 | - input and output values can themselves be lists if you have several |
---|
| 632 | $q$ values to test for the same model parameters. |
---|
| 633 | - for testing *ER* and *VR*, give the inputs as "ER" and "VR" respectively; |
---|
| 634 | the output for *VR* should be the sphere/shell ratio, not the individual |
---|
| 635 | sphere and shell values. |
---|
| 636 | |
---|
| 637 | .. _Test_Your_New_Model: |
---|
| 638 | |
---|
| 639 | Test Your New Model |
---|
| 640 | ^^^^^^^^^^^^^^^^^^^ |
---|
| 641 | |
---|
| 642 | If you are editing your model from the SasView GUI, you can test it |
---|
| 643 | by selecting *Run -> Compile* from the *Model Editor* menu bar. An |
---|
| 644 | *Info* box will appear with the results of the compilation and a |
---|
| 645 | check that the model runs. |
---|
| 646 | |
---|
| 647 | If the model compiles and runs, you can next run the unit tests that |
---|
| 648 | you have added using the **test=** values. Switch to the *Shell* tab |
---|
| 649 | and type the following:: |
---|
| 650 | |
---|
| 651 | from sasmodels.model_test import run_one |
---|
| 652 | run_one("~/.sasview/plugin_models/model.py") |
---|
| 653 | |
---|
| 654 | This should print:: |
---|
| 655 | |
---|
| 656 | test_model_python (sasmodels.model_test.ModelTestCase) ... ok |
---|
| 657 | |
---|
| 658 | To check whether single precision is good enough, type the following:: |
---|
| 659 | |
---|
| 660 | from sasmodels.compare import main |
---|
| 661 | main("~/.sasview/plugin_models/model.py") |
---|
| 662 | |
---|
| 663 | This will pop up a plot showing the difference between single precision |
---|
| 664 | and double precision on a range of $q$ values. |
---|
| 665 | |
---|
| 666 | :: |
---|
| 667 | |
---|
| 668 | demo = dict(scale=1, background=0, |
---|
| 669 | sld=6, sld_solvent=1, |
---|
| 670 | radius=120, |
---|
| 671 | radius_pd=.2, radius_pd_n=45) |
---|
| 672 | |
---|
| 673 | **demo={'par': value, ...}** in the model file sets the default values for |
---|
| 674 | the comparison. You can include polydispersity parameters such as |
---|
| 675 | *radius_pd=0.2, radius_pd_n=45* which would otherwise be zero. |
---|
| 676 | |
---|
| 677 | The options to compare are quite extensive; type the following for help:: |
---|
| 678 | |
---|
| 679 | main() |
---|
| 680 | |
---|
| 681 | Options will need to be passed as separate strings. |
---|
| 682 | For example to run your model with a random set of parameters:: |
---|
| 683 | |
---|
| 684 | main("-random", "-pars", "~/.sasview/plugin_models/model.py") |
---|
| 685 | |
---|
| 686 | For the random models, |
---|
| 687 | |
---|
| 688 | - sld will be in(-0.5,10.5), |
---|
| 689 | - angles (theta, phi, psi) will be in (-180,180), |
---|
| 690 | - angular dispersion will be in (0,45), |
---|
| 691 | - polydispersity will be in (0,1) |
---|
| 692 | - other values will be in (0, 2*v) where v is the value of the parameter in demo. |
---|
| 693 | |
---|
| 694 | Dispersion parameters n, sigma and type will be unchanged from demo so that |
---|
| 695 | run times are predictable. |
---|
| 696 | |
---|
| 697 | If your model has 2D orientational calculation, then you should also |
---|
| 698 | test with:: |
---|
| 699 | |
---|
| 700 | main("-2d", "~/.sasview/plugin_models/model.py") |
---|
| 701 | |
---|
| 702 | |
---|
| 703 | Clean Lint |
---|
| 704 | ^^^^^^^^^^ |
---|
| 705 | |
---|
| 706 | **NB: For now we are not providing pylint with SasView; unless you have a |
---|
| 707 | SasView development environment available, you can ignore this section.** |
---|
| 708 | |
---|
| 709 | Run the lint check with:: |
---|
| 710 | |
---|
| 711 | python -m pylint --rcfile=extra/pylint.rc ~/.sasview/plugin_models/model.py |
---|
| 712 | |
---|
| 713 | We are not aiming for zero lint just yet, only keeping it to a minimum. |
---|
| 714 | For now, don't worry too much about *invalid-name*. If you really want a |
---|
| 715 | variable name *Rg* for example because $R_g$ is the right name for the model |
---|
| 716 | parameter then ignore the lint errors. Also, ignore *missing-docstring* |
---|
| 717 | for standard model functions *Iq*, *Iqxy*, etc. |
---|
| 718 | |
---|
| 719 | We will have delinting sessions at the SasView code camps, where we can |
---|
| 720 | decide on standards for model files, parameter names, etc. |
---|
| 721 | |
---|
| 722 | For now, you can tell pylint to ignore things. For example, to align you |
---|
| 723 | parameters in blocks:: |
---|
| 724 | |
---|
| 725 | # pylint: disable=bad-whitespace,line-too-long |
---|
| 726 | # ["name", "units", default, [lower, upper], "type", "description"], |
---|
| 727 | parameters = [ |
---|
| 728 | ["contrast_factor", "barns", 10.0, [-inf, inf], "", "Contrast factor of the polymer"], |
---|
| 729 | ["bjerrum_length", "Ang", 7.1, [0, inf], "", "Bjerrum length"], |
---|
| 730 | ["virial_param", "1/Ang^2", 12.0, [-inf, inf], "", "Virial parameter"], |
---|
| 731 | ["monomer_length", "Ang", 10.0, [0, inf], "", "Monomer length"], |
---|
| 732 | ["salt_concentration", "mol/L", 0.0, [-inf, inf], "", "Concentration of monovalent salt"], |
---|
| 733 | ["ionization_degree", "", 0.05, [0, inf], "", "Degree of ionization"], |
---|
| 734 | ["polymer_concentration", "mol/L", 0.7, [0, inf], "", "Polymer molar concentration"], |
---|
| 735 | ] |
---|
| 736 | # pylint: enable=bad-whitespace,line-too-long |
---|
| 737 | |
---|
| 738 | Don't put in too many pylint statements, though, since they make the code ugly. |
---|
| 739 | |
---|
| 740 | Check The Docs |
---|
| 741 | ^^^^^^^^^^^^^^ |
---|
| 742 | |
---|
| 743 | You can get a rough idea of how the documentation will look using the |
---|
| 744 | following:: |
---|
| 745 | |
---|
| 746 | from sasmodels.generate import view_html |
---|
| 747 | view_html('~/.sasview/plugin_models/model.py') |
---|
| 748 | |
---|
| 749 | This does not use the same styling as the SasView docs, but it will allow |
---|
| 750 | you to check that your ReStructuredText and LaTeX formatting. Here are |
---|
| 751 | some tools to help with the inevitable syntax errors: |
---|
| 752 | |
---|
| 753 | - `Sphinx cheat sheet <http://matplotlib.org/sampledoc/cheatsheet.html>`_ |
---|
| 754 | - `Sphinx Documentation <http://www.sphinx-doc.org/en/stable/>`_ |
---|
| 755 | - `MathJax <http://www.mathjax.org/>`_ |
---|
| 756 | - `amsmath <http://www.ams.org/publications/authors/tex/amslatex>`_ |
---|
| 757 | |
---|
| 758 | Finally |
---|
| 759 | ^^^^^^^ |
---|
| 760 | |
---|
| 761 | Once compare and the unit test(s) pass properly and everything is done, |
---|
| 762 | consider adding your model to the |
---|
| 763 | `model marketplace <http://marketplace.sasview.org/>`_. |
---|
| 764 | |
---|