source: sasmodels/sasmodels/modelinfo.py @ c1799d3

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

Merge branch 'beta_approx' into ticket-1157

  • Property mode set to 100644
File size: 46.2 KB
RevLine 
[0d99a6a]1"""
2Model Info and Parameter Tables
3===============================
4
5Defines :class:`ModelInfo` and :class:`ParameterTable` and the routines for
6manipulating them.  In particular, :func:`make_model_info` converts a kernel
7module into the model info block as seen by the rest of the sasmodels library.
8"""
[04045f4]9from __future__ import print_function
10
[6d6508e]11from copy import copy
12from os.path import abspath, basename, splitext
[da63656]13import inspect
[69aa451]14
[7ae2b7f]15import numpy as np  # type: ignore
[69aa451]16
[a86f6c0]17# Optional typing
[2d81cfe]18# pylint: disable=unused-import
[a86f6c0]19try:
[04045f4]20    from typing import Tuple, List, Union, Dict, Optional, Any, Callable, Sequence, Set
[e65c3ba]21    from types import ModuleType
[a86f6c0]22except ImportError:
23    pass
24else:
25    Limits = Tuple[float, float]
[04dc697]26    #LimitsOrChoice = Union[Limits, Tuple[Sequence[str]]]
[9a943d0]27    ParameterDef = Tuple[str, str, float, Limits, str, str]
[a86f6c0]28    ParameterSetUser = Dict[str, Union[float, List[float]]]
29    ParameterSet = Dict[str, float]
30    TestInput = Union[str, float, List[float], Tuple[float, float], List[Tuple[float, float]]]
31    TestValue = Union[float, List[float]]
32    TestCondition = Tuple[ParameterSetUser, TestInput, TestValue]
[2d81cfe]33# pylint: enable=unused-import
[a86f6c0]34
[6aee3ab]35# If MAX_PD changes, need to change the loop macros in kernel_iq.c
36MAX_PD = 5 #: Maximum number of simultaneously polydisperse parameters
[d19962c]37
[69aa451]38# assumptions about common parameters exist throughout the code, such as:
[108e70e]39# (1) kernel functions Iq, Iqxy, Iqac, Iqabc, form_volume, ... don't see them
[69aa451]40# (2) kernel drivers assume scale is par[0] and background is par[1]
41# (3) mixture models drop the background on components and replace the scale
42#     with a scale that varies from [-inf, inf]
43# (4) product models drop the background and reassign scale
44# and maybe other places.
45# Note that scale and background cannot be coordinated parameters whose value
46# depends on the some polydisperse parameter with the current implementation
[7b9e4dd]47DEFAULT_BACKGROUND = 1e-3
[7edce22]48COMMON_PARAMETERS = [
[9a943d0]49    ("scale", "", 1, (0.0, np.inf), "", "Source intensity"),
[7b9e4dd]50    ("background", "1/cm", DEFAULT_BACKGROUND, (-np.inf, np.inf), "", "Source background"),
[7edce22]51]
52assert (len(COMMON_PARAMETERS) == 2
[40a87fa]53        and COMMON_PARAMETERS[0][0] == "scale"
54        and COMMON_PARAMETERS[1][0] == "background"), "don't change common parameters"
[69aa451]55
[a86f6c0]56
[69aa451]57def make_parameter_table(pars):
[f6029fd]58    # type: (List[ParameterDef]) -> ParameterTable
[7edce22]59    """
60    Construct a parameter table from a list of parameter definitions.
61
62    This is used by the module processor to convert the parameter block into
63    the parameter table seen in the :class:`ModelInfo` for the module.
64    """
[69aa451]65    processed = []
66    for p in pars:
[a86f6c0]67        if not isinstance(p, (list, tuple)) or len(p) != 6:
[69aa451]68            raise ValueError("Parameter should be [name, units, default, limits, type, desc], but got %r"
69                             %str(p))
70        processed.append(parse_parameter(*p))
71    partable = ParameterTable(processed)
[95498a3]72    partable.check_angles()
[69aa451]73    return partable
74
[a86f6c0]75def parse_parameter(name, units='', default=np.NaN,
[04045f4]76                    user_limits=None, ptype='', description=''):
[60f03de]77    # type: (str, str, float, Sequence[Any], str, str) -> Parameter
[7edce22]78    """
79    Parse an individual parameter from the parameter definition block.
80
81    This does type and value checking on the definition, leading
82    to early failure in the model loading process and easier debugging.
83    """
[69aa451]84    # Parameter is a user facing class.  Do robust type checking.
85    if not isstr(name):
86        raise ValueError("expected string for parameter name %r"%name)
87    if not isstr(units):
88        raise ValueError("expected units to be a string for %s"%name)
[04045f4]89
90    # Process limits as [float, float] or [[str, str, ...]]
[c1a888b]91    choices = []  # type: List[str]
[04045f4]92    if user_limits is None:
93        limits = (-np.inf, np.inf)
94    elif not isinstance(user_limits, (tuple, list)):
95        raise ValueError("invalid limits for %s"%name)
96    else:
97        # if limits is [[str,...]], then this is a choice list field,
98        # and limits are 1 to length of string list
99        if isinstance(user_limits[0], (tuple, list)):
100            choices = user_limits[0]
101            limits = (0., len(choices)-1.)
102            if not all(isstr(k) for k in choices):
103                raise ValueError("choices must be strings for %s"%name)
104        else:
105            try:
106                low, high = user_limits
107                limits = (float(low), float(high))
108            except Exception:
[bb4b509]109                raise ValueError("invalid limits for %s: %r"%(name, user_limits))
[50ec515]110            if low >= high:
111                raise ValueError("require lower limit < upper limit")
[69aa451]112
[04045f4]113    # Process default value as float, making sure it is in range
[69aa451]114    if not isinstance(default, (int, float)):
115        raise ValueError("expected default %r to be a number for %s"
116                         % (default, name))
[04045f4]117    if default < limits[0] or default > limits[1]:
[69aa451]118        raise ValueError("default value %r not in range for %s"
119                         % (default, name))
120
[04045f4]121    # Check for valid parameter type
[a86f6c0]122    if ptype not in ("volume", "orientation", "sld", "magnetic", ""):
123        raise ValueError("unexpected type %r for %s" % (ptype, name))
[69aa451]124
[04045f4]125    # Check for valid parameter description
[69aa451]126    if not isstr(description):
127        raise ValueError("expected description to be a string")
128
129    # Parameter id for name[n] does not include [n]
130    if "[" in name:
131        if not name.endswith(']'):
132            raise ValueError("Expected name[len] for vector parameter %s"%name)
133        pid, ref = name[:-1].split('[', 1)
134        ref = ref.strip()
135    else:
136        pid, ref = name, None
137
138    # automatically identify sld types
[40a87fa]139    if ptype == '' and (pid.startswith('sld') or pid.endswith('sld')):
[a86f6c0]140        ptype = 'sld'
[69aa451]141
142    # Check if using a vector definition, name[k], as the parameter name
143    if ref:
144        if ref == '':
145            raise ValueError("Need to specify vector length for %s"%name)
146        try:
147            length = int(ref)
148            control = None
[a86f6c0]149        except ValueError:
[69aa451]150            length = None
151            control = ref
152    else:
153        length = 1
154        control = None
155
156    # Build the parameter
157    parameter = Parameter(name=name, units=units, default=default,
[a86f6c0]158                          limits=limits, ptype=ptype, description=description)
[69aa451]159
160    # TODO: need better control over whether a parameter is polydisperse
[a86f6c0]161    parameter.polydisperse = ptype in ('orientation', 'volume')
162    parameter.relative_pd = ptype == 'volume'
[69aa451]163    parameter.choices = choices
164    parameter.length = length
165    parameter.length_control = control
166    return parameter
167
[ce896fd]168
[4bfd277]169def expand_pars(partable, pars):
[a86f6c0]170    # type: (ParameterTable, ParameterSetUser) ->  ParameterSet
[ce896fd]171    """
[6d6508e]172    Create demo parameter set from key-value pairs.
[ce896fd]173
[4bfd277]174    *pars* are the key-value pairs to use for the parameters.  Any
175    parameters not specified in *pars* are set from the *partable* defaults.
[ce896fd]176
[4bfd277]177    If *pars* references vector fields, such as thickness[n], then support
[ce896fd]178    different ways of assigning the demo values, including assigning a
179    specific value (e.g., thickness3=50.0), assigning a new value to all
180    (e.g., thickness=50.0) or assigning values using list notation.
181    """
[4bfd277]182    if pars is None:
[ce896fd]183        result = partable.defaults
184    else:
[4bfd277]185        lookup = dict((p.id, p) for p in partable.kernel_parameters)
[ce896fd]186        result = partable.defaults.copy()
[9a943d0]187        scalars = dict((name, value) for name, value in pars.items()
188                       if name not in lookup or lookup[name].length == 1)
[40a87fa]189        vectors = dict((name, value) for name, value in pars.items()
[4bfd277]190                       if name in lookup and lookup[name].length > 1)
[256dfe1]191        #print("lookup", lookup)
192        #print("scalars", scalars)
193        #print("vectors", vectors)
[ce896fd]194        if vectors:
195            for name, value in vectors.items():
196                if np.isscalar(value):
197                    # support for the form
[4bfd277]198                    #    dict(thickness=0, thickness2=50)
199                    for k in range(1, lookup[name].length+1):
[ce896fd]200                        key = name+str(k)
201                        if key not in scalars:
[256dfe1]202                            scalars[key] = value
[ce896fd]203                else:
204                    # supoprt for the form
[4bfd277]205                    #    dict(thickness=[20,10,3])
[40a87fa]206                    for (k, v) in enumerate(value):
[256dfe1]207                        scalars[name+str(k+1)] = v
[9a943d0]208        result.update(scalars)
[256dfe1]209        #print("expanded", result)
[ce896fd]210
[6d6508e]211    return result
212
213def prefix_parameter(par, prefix):
[a86f6c0]214    # type: (Parameter, str) -> Parameter
[6d6508e]215    """
216    Return a copy of the parameter with its name prefixed.
217    """
218    new_par = copy(par)
219    new_par.name = prefix + par.name
220    new_par.id = prefix + par.id
221
222def suffix_parameter(par, suffix):
[a86f6c0]223    # type: (Parameter, str) -> Parameter
[6d6508e]224    """
225    Return a copy of the parameter with its name prefixed.
226    """
227    new_par = copy(par)
228    # If name has the form x[n], replace with x_suffix[n]
229    new_par.name = par.id + suffix + par.name[len(par.id):]
230    new_par.id = par.id + suffix
[ce896fd]231
[69aa451]232class Parameter(object):
233    """
234    The available kernel parameters are defined as a list, with each parameter
235    defined as a sublist with the following elements:
236
[f88e248]237    *name* is the name that will be displayed to the user.  Names
[69aa451]238    should be lower case, with words separated by underscore.  If
[f88e248]239    acronyms are used, the whole acronym should be upper case. For vector
240    parameters, the name will be followed by *[len]* where *len* is an
241    integer length of the vector, or the name of the parameter which
242    controls the length.  The attribute *id* will be created from name
243    without the length.
[69aa451]244
245    *units* should be one of *degrees* for angles, *Ang* for lengths,
246    *1e-6/Ang^2* for SLDs.
247
248    *default value* will be the initial value for  the model when it
249    is selected, or when an initial value is not otherwise specified.
250
251    *limits = [lb, ub]* are the hard limits on the parameter value, used to
252    limit the polydispersity density function.  In the fit, the parameter limits
253    given to the fit are the limits  on the central value of the parameter.
254    If there is polydispersity, it will evaluate parameter values outside
255    the fit limits, but not outside the hard limits specified in the model.
256    If there are no limits, use +/-inf imported from numpy.
257
258    *type* indicates how the parameter will be used.  "volume" parameters
[108e70e]259    will be used in all functions.  "orientation" parameters are not passed,
260    but will be used to convert from *qx*, *qy* to *qa*, *qb*, *qc* in calls to
261    *Iqxy* and *Imagnetic*.  If *type* is the empty string, the parameter will
[69aa451]262    be used in all of *Iq*, *Iqxy* and *Imagnetic*.  "sld" parameters
263    can automatically be promoted to magnetic parameters, each of which
264    will have a magnitude and a direction, which may be different from
[6e7ff6d]265    other sld parameters. The volume parameters are used for calls
[39a06c9]266    to form_volume within the kernel (required for volume normalization),
267    to shell_volume (for hollow shapes), and to effective_radius (for
268    structure factor interactions) respectively.
[69aa451]269
270    *description* is a short description of the parameter.  This will
271    be displayed in the parameter table and used as a tool tip for the
272    parameter value in the user interface.
273
274    Additional values can be set after the parameter is created:
275
[6e7ff6d]276    * *length* is the length of the field if it is a vector field
277
278    * *length_control* is the parameter which sets the vector length
279
280    * *is_control* is True if the parameter is a control parameter for a vector
281
282    * *polydisperse* is true if the parameter accepts a polydispersity
283
284    * *relative_pd* is true if that polydispersity is a portion of the
[40a87fa]285      value (so a 10% length dipsersity would use a polydispersity value
286      of 0.1) rather than absolute dispersisity (such as an angle plus or
287      minus 15 degrees).
[69aa451]288
[50ec515]289    *choices* is the option names for a drop down list of options, as for
290    example, might be used to set the value of a shape parameter.
291
[a5b8477]292    These values are set by :func:`make_parameter_table` and
293    :func:`parse_parameter` therein.
[69aa451]294    """
295    def __init__(self, name, units='', default=None, limits=(-np.inf, np.inf),
[a86f6c0]296                 ptype='', description=''):
[c1a888b]297        # type: (str, str, float, Limits, str, str) -> None
[a86f6c0]298        self.id = name.split('[')[0].strip() # type: str
299        self.name = name                     # type: str
300        self.units = units                   # type: str
301        self.default = default               # type: float
302        self.limits = limits                 # type: Limits
303        self.type = ptype                    # type: str
304        self.description = description       # type: str
[69aa451]305
[d19962c]306        # Length and length_control will be filled in once the complete
[69aa451]307        # parameter table is available.
[a86f6c0]308        self.length = 1                      # type: int
309        self.length_control = None           # type: Optional[str]
310        self.is_control = False              # type: bool
[69aa451]311
312        # TODO: need better control over whether a parameter is polydisperse
[a86f6c0]313        self.polydisperse = False            # type: bool
314        self.relative_pd = False             # type: bool
[69aa451]315
[ce896fd]316        # choices are also set externally.
[a86f6c0]317        self.choices = []                    # type: List[str]
[ce896fd]318
[69aa451]319    def as_definition(self):
[a86f6c0]320        # type: () -> str
[69aa451]321        """
322        Declare space for the variable in a parameter structure.
323
324        For example, the parameter thickness with length 3 will
325        return "double thickness[3];", with no spaces before and
326        no new line character afterward.
327        """
328        if self.length == 1:
329            return "double %s;"%self.id
330        else:
331            return "double %s[%d];"%(self.id, self.length)
332
333    def as_function_argument(self):
[a86f6c0]334        # type: () -> str
[40a87fa]335        r"""
[69aa451]336        Declare the variable as a function argument.
337
338        For example, the parameter thickness with length 3 will
[40a87fa]339        return "double \*thickness", with no spaces before and
[69aa451]340        no comma afterward.
341        """
342        if self.length == 1:
343            return "double %s"%self.id
344        else:
345            return "double *%s"%self.id
346
347    def as_call_reference(self, prefix=""):
[c1a888b]348        # type: (str) -> str
[bb4b509]349        """
350        Return *prefix* + parameter name.  For struct references, use "v."
351        as the prefix.
352        """
[69aa451]353        # Note: if the parameter is a struct type, then we will need to use
354        # &prefix+id.  For scalars and vectors we can just use prefix+id.
355        return prefix + self.id
356
357    def __str__(self):
[a86f6c0]358        # type: () -> str
[69aa451]359        return "<%s>"%self.name
360
361    def __repr__(self):
[a86f6c0]362        # type: () -> str
[69aa451]363        return "P<%s>"%self.name
364
[d19962c]365
[69aa451]366class ParameterTable(object):
[d19962c]367    """
368    ParameterTable manages the list of available parameters.
369
370    There are a couple of complications which mean that the list of parameters
371    for the kernel differs from the list of parameters that the user sees.
372
373    (1) Common parameters.  Scale and background are implicit to every model,
374    but are not passed to the kernel.
375
376    (2) Vector parameters.  Vector parameters are passed to the kernel as a
377    pointer to an array, e.g., thick[], but they are seen by the user as n
378    separate parameters thick1, thick2, ...
379
380    Therefore, the parameter table is organized by how it is expected to be
381    used. The following information is needed to set up the kernel functions:
382
383    * *kernel_parameters* is the list of parameters in the kernel parameter
[40a87fa]384      table, with vector parameter p declared as p[].
[d19962c]385
386    * *iq_parameters* is the list of parameters to the Iq(q, ...) function,
[40a87fa]387      with vector parameter p sent as p[].
[d19962c]388
389    * *form_volume_parameters* is the list of parameters to the form_volume(...)
[40a87fa]390      function, with vector parameter p sent as p[].
[d19962c]391
392    Problem details, which sets up the polydispersity loops, requires the
393    following:
394
395    * *theta_offset* is the offset of the theta parameter in the kernel parameter
[40a87fa]396      table, with vector parameters counted as n individual parameters
397      p1, p2, ..., or offset is -1 if there is no theta parameter.
[d19962c]398
399    * *max_pd* is the maximum number of polydisperse parameters, with vector
[40a87fa]400      parameters counted as n individual parameters p1, p2, ...  Note that
401      this number is limited to sasmodels.modelinfo.MAX_PD.
[d19962c]402
403    * *npars* is the total number of parameters to the kernel, with vector
[40a87fa]404      parameters counted as n individual parameters p1, p2, ...
[d19962c]405
[b6dab59]406    * *common_parameters* is the list of common parameters, with a unique
407      copy for each model so that structure factors can have a default
408      background of 0.0.
409
[d19962c]410    * *call_parameters* is the complete list of parameters to the kernel,
[40a87fa]411      including scale and background, with vector parameters recorded as
412      individual parameters p1, p2, ...
[d19962c]413
414    * *active_1d* is the set of names that may be polydisperse for 1d data
415
416    * *active_2d* is the set of names that may be polydisperse for 2d data
417
418    User parameters are the set of parameters visible to the user, including
419    the scale and background parameters that the kernel does not see.  User
420    parameters don't use vector notation, and instead use p1, p2, ...
421    """
[69aa451]422    def __init__(self, parameters):
[a86f6c0]423        # type: (List[Parameter]) -> None
[b6dab59]424
425        # scale and background are implicit parameters
426        # Need them to be unique to each model in case they have different
427        # properties, such as default=0.0 for structure factor backgrounds.
428        self.common_parameters = [Parameter(*p) for p in COMMON_PARAMETERS]
429
[d19962c]430        self.kernel_parameters = parameters
[60eab2a]431        self._set_vector_lengths()
[32e3c9b]432        self.npars = sum(p.length for p in self.kernel_parameters)
[9eb3632]433        self.nmagnetic = sum(p.length for p in self.kernel_parameters
[bb4b509]434                             if p.type == 'sld')
[9eb3632]435        self.nvalues = 2 + self.npars
436        if self.nmagnetic:
437            self.nvalues += 3 + 3*self.nmagnetic
[a86f6c0]438        self.call_parameters = self._get_call_parameters()
439        self.defaults = self._get_defaults()
[d19962c]440        #self._name_table= dict((p.id, p) for p in parameters)
[60eab2a]441
[d19962c]442        # Set the kernel parameters.  Assumes background and scale are the
443        # first two parameters in the parameter list, but these are not sent
444        # to the underlying kernel functions.
445        self.iq_parameters = [p for p in self.kernel_parameters
446                              if p.type not in ('orientation', 'magnetic')]
[108e70e]447        self.orientation_parameters = [p for p in self.kernel_parameters
448                                       if p.type == 'orientation']
[d19962c]449        self.form_volume_parameters = [p for p in self.kernel_parameters
450                                       if p.type == 'volume']
451
452        # Theta offset
453        offset = 0
454        for p in self.kernel_parameters:
455            if p.name == 'theta':
456                self.theta_offset = offset
457                break
458            offset += p.length
459        else:
460            self.theta_offset = -1
[69aa451]461
[d19962c]462        # number of polydisperse parameters
463        num_pd = sum(p.length for p in self.kernel_parameters if p.polydisperse)
464        # Don't use more polydisperse parameters than are available in the model
465        self.max_pd = min(num_pd, MAX_PD)
[69aa451]466
[d19962c]467        # true if has 2D parameters
468        self.has_2d = any(p.type in ('orientation', 'magnetic')
469                          for p in self.kernel_parameters)
[8698a0d]470        self.is_asymmetric = any(p.name == 'psi' for p in self.kernel_parameters)
[bb4b509]471        self.magnetism_index = [k for k, p in enumerate(self.call_parameters)
[610ef23]472                                if p.id.endswith('_M0')]
[69aa451]473
[d19962c]474        self.pd_1d = set(p.name for p in self.call_parameters
[a86f6c0]475                         if p.polydisperse and p.type not in ('orientation', 'magnetic'))
[32e3c9b]476        self.pd_2d = set(p.name for p in self.call_parameters if p.polydisperse)
[69aa451]477
[0535624]478    def set_zero_background(self):
479        """
480        Set the default background to zero for this model.  This is done for
481        structure factor models.
482        """
483        # type: () -> None
484        # Make sure background is the second common parameter.
485        assert self.common_parameters[1].id == "background"
486        self.common_parameters[1].default = 0.0
487        self.defaults = self._get_defaults()
488
[95498a3]489    def check_angles(self):
490        """
491        Check that orientation angles are theta, phi and possibly psi.
492        """
[8698a0d]493        theta = phi = psi = -1
494        for k, p in enumerate(self.kernel_parameters):
495            if p.name == 'theta':
496                theta = k
497                if p.type != 'orientation':
498                    raise TypeError("theta must be an orientation parameter")
499            elif p.name == 'phi':
500                phi = k
501                if p.type != 'orientation':
502                    raise TypeError("phi must be an orientation parameter")
503            elif p.name == 'psi':
504                psi = k
505                if p.type != 'orientation':
506                    raise TypeError("psi must be an orientation parameter")
[108e70e]507            elif p.type == 'orientation':
508                raise TypeError("only theta, phi and psi can be orientation parameters")
[8698a0d]509        if theta >= 0 and phi >= 0:
[108e70e]510            last_par = len(self.kernel_parameters) - 1
[8698a0d]511            if phi != theta+1:
512                raise TypeError("phi must follow theta")
513            if psi >= 0 and psi != phi+1:
514                raise TypeError("psi must follow phi")
[95498a3]515            if (psi >= 0 and psi != last_par) or (psi < 0 and phi != last_par):
516                raise TypeError("orientation parameters must appear at the "
517                                "end of the parameter table")
[8698a0d]518        elif theta >= 0 or phi >= 0 or psi >= 0:
519            raise TypeError("oriented shapes must have both theta and phi and maybe psi")
520
[0bdddc2]521    def __getitem__(self, key):
522        # Find the parameter definition
523        for par in self.call_parameters:
524            if par.name == key:
[e65c3ba]525                return par
526        raise KeyError("unknown parameter %r"%key)
[0bdddc2]527
[e3571cb]528    def __contains__(self, key):
529        for par in self.call_parameters:
530            if par.name == key:
531                return True
[e65c3ba]532        return False
[e3571cb]533
[a86f6c0]534    def _set_vector_lengths(self):
[256dfe1]535        # type: () -> List[str]
[a86f6c0]536        """
537        Walk the list of kernel parameters, setting the length field of the
538        vector parameters from the upper limit of the reference parameter.
539
540        This needs to be done once the entire parameter table is available
541        since the reference may still be undefined when the parameter is
542        initially created.
543
[256dfe1]544        Returns the list of control parameter names.
545
[a86f6c0]546        Note: This modifies the underlying parameter object.
547        """
548        # Sort out the length of the vector parameters such as thickness[n]
549        for p in self.kernel_parameters:
550            if p.length_control:
[e65c3ba]551                ref = self._get_ref(p)
[a86f6c0]552                ref.is_control = True
[256dfe1]553                ref.polydisperse = False
[a86f6c0]554                low, high = ref.limits
555                if int(low) != low or int(high) != high or low < 0 or high > 20:
556                    raise ValueError("expected limits on %s to be within [0, 20]"
557                                     % ref.name)
[04045f4]558                p.length = int(high)
[a86f6c0]559
[e65c3ba]560    def _get_ref(self, p):
561        # type: (Parameter) -> Parameter
562        for ref in self.kernel_parameters:
563            if ref.id == p.length_control:
564                return ref
565        raise ValueError("no reference variable %r for %s"
566                         % (p.length_control, p.name))
567
[a86f6c0]568    def _get_defaults(self):
569        # type: () -> ParameterSet
570        """
571        Get a list of parameter defaults from the parameters.
572
573        Expands vector parameters into parameter id+number.
574        """
575        # Construct default values, including vector defaults
576        defaults = {}
577        for p in self.call_parameters:
578            if p.length == 1:
579                defaults[p.id] = p.default
580            else:
581                for k in range(1, p.length+1):
582                    defaults["%s%d"%(p.id, k)] = p.default
583        return defaults
584
585    def _get_call_parameters(self):
586        # type: () -> List[Parameter]
[b6dab59]587        full_list = self.common_parameters[:]
[a86f6c0]588        for p in self.kernel_parameters:
589            if p.length == 1:
590                full_list.append(p)
591            else:
592                for k in range(1, p.length+1):
593                    pk = Parameter(p.id+str(k), p.units, p.default,
594                                   p.limits, p.type, p.description)
595                    pk.polydisperse = p.polydisperse
596                    pk.relative_pd = p.relative_pd
[50ec515]597                    pk.choices = p.choices
[a86f6c0]598                    full_list.append(pk)
[32e3c9b]599
600        # Add the magnetic parameters to the end of the call parameter list.
[9eb3632]601        if self.nmagnetic > 0:
[32e3c9b]602            full_list.extend([
[610ef23]603                Parameter('up_frac_i', '', 0., [0., 1.],
[32e3c9b]604                          'magnetic', 'fraction of spin up incident'),
[610ef23]605                Parameter('up_frac_f', '', 0., [0., 1.],
[32e3c9b]606                          'magnetic', 'fraction of spin up final'),
[610ef23]607                Parameter('up_angle', 'degrees', 0., [0., 360.],
[32e3c9b]608                          'magnetic', 'spin up angle'),
609            ])
610            slds = [p for p in full_list if p.type == 'sld']
611            for p in slds:
612                full_list.extend([
[610ef23]613                    Parameter(p.id+'_M0', '1e-6/Ang^2', 0., [-np.inf, np.inf],
[32e3c9b]614                              'magnetic', 'magnetic amplitude for '+p.description),
[610ef23]615                    Parameter(p.id+'_mtheta', 'degrees', 0., [-90., 90.],
[bb4b509]616                              'magnetic', 'magnetic latitude for '+p.description),
[610ef23]617                    Parameter(p.id+'_mphi', 'degrees', 0., [-180., 180.],
[bb4b509]618                              'magnetic', 'magnetic longitude for '+p.description),
[32e3c9b]619                ])
620        #print("call parameters", full_list)
[a86f6c0]621        return full_list
622
[85fe7f8]623    def user_parameters(self, pars, is2d=True):
[0d99a6a]624        # type: (Dict[str, float], bool) -> List[Parameter]
[69aa451]625        """
[d19962c]626        Return the list of parameters for the given data type.
627
[c5ac2b2]628        Vector parameters are expanded in place.  If multiple parameters
[d19962c]629        share the same vector length, then the parameters will be interleaved
630        in the result.  The control parameters come first.  For example,
631        if the parameter table is ordered as::
632
633            sld_core
634            sld_shell[num_shells]
635            sld_solvent
636            thickness[num_shells]
637            num_shells
638
639        and *pars[num_shells]=2* then the returned list will be::
640
641            num_shells
642            scale
643            background
644            sld_core
645            sld_shell1
646            thickness1
647            sld_shell2
648            thickness2
649            sld_solvent
650
651        Note that shell/thickness pairs are grouped together in the result
652        even though they were not grouped in the incoming table.  The control
653        parameter is always returned first since the GUI will want to set it
654        early, and rerender the table when it is changed.
[32e3c9b]655
656        Parameters marked as sld will automatically have a set of associated
[610ef23]657        magnetic parameters (p_M0, p_mtheta, p_mphi), as well as polarization
658        information (up_theta, up_frac_i, up_frac_f).
[69aa451]659        """
[04045f4]660        # control parameters go first
[d19962c]661        control = [p for p in self.kernel_parameters if p.is_control]
662
663        # Gather entries such as name[n] into groups of the same n
[04045f4]664        dependent = {} # type: Dict[str, List[Parameter]]
665        dependent.update((p.id, []) for p in control)
[d19962c]666        for p in self.kernel_parameters:
667            if p.length_control is not None:
668                dependent[p.length_control].append(p)
669
670        # Gather entries such as name[4] into groups of the same length
[32e3c9b]671        fixed_length = {}  # type: Dict[int, List[Parameter]]
[d19962c]672        for p in self.kernel_parameters:
673            if p.length > 1 and p.length_control is None:
[32e3c9b]674                fixed_length.setdefault(p.length, []).append(p)
[d19962c]675
676        # Using the call_parameters table, we already have expanded forms
677        # for each of the vector parameters; put them in a lookup table
[f88e248]678        # Note: p.id and p.name are currently identical for the call parameters
679        expanded_pars = dict((p.id, p) for p in self.call_parameters)
[d19962c]680
[32e3c9b]681        def append_group(name):
682            """add the named parameter, and related magnetic parameters if any"""
683            result.append(expanded_pars[name])
684            if is2d:
[610ef23]685                for tag in '_M0', '_mtheta', '_mphi':
686                    if name+tag in expanded_pars:
687                        result.append(expanded_pars[name+tag])
[32e3c9b]688
[d19962c]689        # Gather the user parameters in order
[b6dab59]690        result = control + self.common_parameters
[d19962c]691        for p in self.kernel_parameters:
692            if not is2d and p.type in ('orientation', 'magnetic'):
693                pass
694            elif p.is_control:
695                pass # already added
696            elif p.length_control is not None:
697                table = dependent.get(p.length_control, [])
698                if table:
699                    # look up length from incoming parameters
[fb5914f]700                    table_length = int(pars.get(p.length_control, p.length))
[d19962c]701                    del dependent[p.length_control] # first entry seen
702                    for k in range(1, table_length+1):
703                        for entry in table:
[32e3c9b]704                            append_group(entry.id+str(k))
[d19962c]705                else:
706                    pass # already processed all entries
707            elif p.length > 1:
[32e3c9b]708                table = fixed_length.get(p.length, [])
[d19962c]709                if table:
710                    table_length = p.length
[32e3c9b]711                    del fixed_length[p.length]
[d19962c]712                    for k in range(1, table_length+1):
713                        for entry in table:
[32e3c9b]714                            append_group(entry.id+str(k))
[d19962c]715                else:
716                    pass # already processed all entries
717            else:
[32e3c9b]718                append_group(p.id)
719
[610ef23]720        if is2d and 'up_angle' in expanded_pars:
[32e3c9b]721            result.extend([
[610ef23]722                expanded_pars['up_frac_i'],
723                expanded_pars['up_frac_f'],
724                expanded_pars['up_angle'],
[32e3c9b]725            ])
[d19962c]726
727        return result
[69aa451]728
729def isstr(x):
[a86f6c0]730    # type: (Any) -> bool
[7edce22]731    """
732    Return True if the object is a string.
733    """
[69aa451]734    # TODO: 2-3 compatible tests for str, including unicode strings
735    return isinstance(x, str)
736
[da63656]737
[108e70e]738#: Set of variables defined in the model that might contain C code
[e44432d]739C_SYMBOLS = ['Imagnetic', 'Iq', 'Iqxy', 'Iqac', 'Iqabc', 'form_volume', 'shell_volume', 'c_code']
[108e70e]740
[da63656]741def _find_source_lines(model_info, kernel_module):
[e65c3ba]742    # type: (ModelInfo, ModuleType) -> None
[da63656]743    """
744    Identify the location of the C source inside the model definition file.
745
[108e70e]746    This code runs through the source of the kernel module looking for lines
747    that contain C code (because it is a c function definition).  Clearly
748    there are all sorts of reasons why this might not work (e.g., code
749    commented out in a triple-quoted line block, code built using string
750    concatenation, code defined in the branch of an 'if' block, code imported
751    from another file), but it should work properly in the 95% case, and for
752    the remainder, getting the incorrect line number will merely be
753    inconvenient.
[da63656]754    """
[108e70e]755    # Only need line numbers if we are creating a C module and the C symbols
756    # are defined.
757    if (callable(model_info.Iq)
758            or not any(hasattr(model_info, s) for s in C_SYMBOLS)):
[da63656]759        return
760
[108e70e]761    # load the module source if we can
[98ba1fc]762    try:
763        source = inspect.getsource(kernel_module)
764    except IOError:
765        return
[da63656]766
[108e70e]767    # look for symbol at the start of the line
768    for lineno, line in enumerate(source.split('\n')):
769        for name in C_SYMBOLS:
770            if line.startswith(name):
771                # Add 1 since some compilers complain about "#line 0"
772                model_info.lineno[name] = lineno + 1
773                break
[da63656]774
[6d6508e]775def make_model_info(kernel_module):
[7edce22]776    # type: (module) -> ModelInfo
777    """
778    Extract the model definition from the loaded kernel module.
779
780    Fill in default values for parts of the module that are not provided.
781
[108e70e]782    Note: vectorized Iq and Iqac/Iqabc functions will be created for python
[7edce22]783    models when the model is first called, not when the model is loaded.
784    """
[65314f7]785    if hasattr(kernel_module, "model_info"):
786        # Custom sum/multi models
787        return kernel_module.model_info
[01c8d9e]788
[6d6508e]789    info = ModelInfo()
[b6dab59]790
791    # Build the parameter table
[6d6508e]792    #print("make parameter table", kernel_module.parameters)
[c1a888b]793    parameters = make_parameter_table(getattr(kernel_module, 'parameters', []))
[b6dab59]794
795    # background defaults to zero for structure factor models
796    structure_factor = getattr(kernel_module, 'structure_factor', False)
797    if structure_factor:
[0535624]798        parameters.set_zero_background()
[b6dab59]799
800    # TODO: remove demo parameters
801    # The plots in the docs are generated from the model default values.
802    # Sascomp set parameters from the command line, and so doesn't need
803    # demo values for testing.
[4bfd277]804    demo = expand_pars(parameters, getattr(kernel_module, 'demo', None))
[b6dab59]805
[724257c]806    filename = abspath(kernel_module.__file__).replace('.pyc', '.py')
[6d6508e]807    kernel_id = splitext(basename(filename))[0]
808    name = getattr(kernel_module, 'name', None)
809    if name is None:
810        name = " ".join(w.capitalize() for w in kernel_id.split('_'))
811
812    info.id = kernel_id  # string used to load the kernel
[724257c]813    info.filename = filename
[6d6508e]814    info.name = name
815    info.title = getattr(kernel_module, 'title', name+" model")
816    info.description = getattr(kernel_module, 'description', 'no description')
817    info.parameters = parameters
818    info.demo = demo
819    info.composition = None
820    info.docs = kernel_module.__doc__
821    info.category = getattr(kernel_module, 'category', None)
822    info.structure_factor = getattr(kernel_module, 'structure_factor', False)
[c036ddb]823    # TODO: find Fq by inspection
[6e7ba14]824    info.effective_radius_type = getattr(kernel_module, 'effective_radius_type', None)
[c036ddb]825    info.have_Fq = getattr(kernel_module, 'have_Fq', False)
[04045f4]826    info.profile_axes = getattr(kernel_module, 'profile_axes', ['x', 'y'])
[d321747]827    # Note: custom.load_custom_kernel_module assumes the C sources are defined
828    # by this attribute.
[6d6508e]829    info.source = getattr(kernel_module, 'source', [])
[108e70e]830    info.c_code = getattr(kernel_module, 'c_code', None)
[5399809]831    info.effective_radius = getattr(kernel_module, 'effective_radius', None)
[e62a134]832    # TODO: check the structure of the tests
[6d6508e]833    info.tests = getattr(kernel_module, 'tests', [])
[9a943d0]834    info.form_volume = getattr(kernel_module, 'form_volume', None) # type: ignore
[e44432d]835    info.shell_volume = getattr(kernel_module, 'shell_volume', None) # type: ignore
[9a943d0]836    info.Iq = getattr(kernel_module, 'Iq', None) # type: ignore
837    info.Iqxy = getattr(kernel_module, 'Iqxy', None) # type: ignore
[108e70e]838    info.Iqac = getattr(kernel_module, 'Iqac', None) # type: ignore
839    info.Iqabc = getattr(kernel_module, 'Iqabc', None) # type: ignore
[32e3c9b]840    info.Imagnetic = getattr(kernel_module, 'Imagnetic', None) # type: ignore
[9a943d0]841    info.profile = getattr(kernel_module, 'profile', None) # type: ignore
842    info.sesans = getattr(kernel_module, 'sesans', None) # type: ignore
[5124c969]843    # Default single and opencl to True for C models.  Python models have callable Iq.
844    info.opencl = getattr(kernel_module, 'opencl', not callable(info.Iq))
845    info.single = getattr(kernel_module, 'single', not callable(info.Iq))
[0bdddc2]846    info.random = getattr(kernel_module, 'random', None)
[256dfe1]847
848    # multiplicity info
849    control_pars = [p.id for p in parameters.kernel_parameters if p.is_control]
850    default_control = control_pars[0] if control_pars else None
851    info.control = getattr(kernel_module, 'control', default_control)
[04045f4]852    info.hidden = getattr(kernel_module, 'hidden', None) # type: ignore
[6d6508e]853
[108e70e]854    if callable(info.Iq) and parameters.has_2d:
855        raise ValueError("oriented python models not supported")
856
857    info.lineno = {}
[da63656]858    _find_source_lines(info, kernel_module)
[6d6508e]859    return info
860
861class ModelInfo(object):
862    """
863    Interpret the model definition file, categorizing the parameters.
864
865    The module can be loaded with a normal python import statement if you
866    know which module you need, or with __import__('sasmodels.model.'+name)
867    if the name is in a string.
868
[7edce22]869    The structure should be mostly static, other than the delayed definition
[108e70e]870    of *Iq*, *Iqac* and *Iqabc* if they need to be defined.
[6d6508e]871    """
[a5b8477]872    #: Full path to the file defining the kernel, if any.
[340428e]873    filename = None         # type: Optional[str]
[a5b8477]874    #: Id of the kernel used to load it from the filesystem.
[a86f6c0]875    id = None               # type: str
[a5b8477]876    #: Display name of the model, which defaults to the model id but with
877    #: capitalization of the parts so for example core_shell defaults to
878    #: "Core Shell".
[a86f6c0]879    name = None             # type: str
[a5b8477]880    #: Short description of the model.
[a86f6c0]881    title = None            # type: str
[a5b8477]882    #: Long description of the model.
[a86f6c0]883    description = None      # type: str
[a5b8477]884    #: Model parameter table. Parameters are defined using a list of parameter
885    #: definitions, each of which is contains parameter name, units,
886    #: default value, limits, type and description.  See :class:`Parameter`
887    #: for details on the individual parameters.  The parameters are gathered
888    #: into a :class:`ParameterTable`, which provides various views into the
889    #: parameter list.
[a86f6c0]890    parameters = None       # type: ParameterTable
[a5b8477]891    #: Demo parameters as a *parameter:value* map used as the default values
892    #: for :mod:`compare`.  Any parameters not set in *demo* will use the
893    #: defaults from the parameter table.  That means no polydispersity, and
894    #: in the case of multiplicity models, a minimal model with no interesting
895    #: scattering.
[a86f6c0]896    demo = None             # type: Dict[str, float]
[a5b8477]897    #: Composition is None if this is an independent model, or it is a
898    #: tuple with comoposition type ('product' or 'misture') and a list of
899    #: :class:`ModelInfo` blocks for the composed objects.  This allows us
900    #: to rebuild a complete mixture or product model from the info block.
901    #: *composition* is not given in the model definition file, but instead
902    #: arises when the model is constructed using names such as
903    #: *sphere*hardsphere* or *cylinder+sphere*.
[a86f6c0]904    composition = None      # type: Optional[Tuple[str, List[ModelInfo]]]
[a5b8477]905    #: Name of the control parameter for a variant model such as :ref:`rpa`.
906    #: The *control* parameter should appear in the parameter table, with
907    #: limits defined as *[CASES]*, for case names such as
908    #: *CASES = ["diblock copolymer", "triblock copolymer", ...]*.
909    #: This should give *limits=[[case1, case2, ...]]*, but the
910    #: model loader translates this to *limits=[0, len(CASES)-1]*, and adds
911    #: *choices=CASES* to the :class:`Parameter` definition. Note that
912    #: models can use a list of cases as a parameter without it being a
913    #: control parameter.  Either way, the parameter is sent to the model
914    #: evaluator as *float(choice_num)*, where choices are numbered from 0.
915    #: See also :attr:`hidden`.
[04045f4]916    control = None          # type: str
[a5b8477]917    #: Different variants require different parameters.  In order to show
918    #: just the parameters needed for the variant selected by :attr:`control`,
919    #: you should provide a function *hidden(control) -> set(['a', 'b', ...])*
920    #: indicating which parameters need to be hidden.  For multiplicity
921    #: models, you need to use the complete name of the parameter, including
922    #: its number.  So for example, if variant "a" uses only *sld1* and *sld2*,
923    #: then *sld3*, *sld4* and *sld5* of multiplicity parameter *sld[5]*
924    #: should be in the hidden set.
925    hidden = None           # type: Optional[Callable[[int], Set[str]]]
926    #: Doc string from the top of the model file.  This should be formatted
927    #: using ReStructuredText format, with latex markup in ".. math"
928    #: environments, or in dollar signs.  This will be automatically
929    #: extracted to a .rst file by :func:`generate.make_docs`, then
930    #: converted to HTML or PDF by Sphinx.
[a86f6c0]931    docs = None             # type: str
[a5b8477]932    #: Location of the model description in the documentation.  This takes the
933    #: form of "section" or "section:subsection".  So for example,
934    #: :ref:`porod` uses *category="shape-independent"* so it is in the
[40a87fa]935    #: :ref:`shape-independent` section whereas
936    #: :ref:`capped-cylinder` uses: *category="shape:cylinder"*, which puts
[a5b8477]937    #: it in the :ref:`shape-cylinder` section.
[a86f6c0]938    category = None         # type: Optional[str]
[a5b8477]939    #: True if the model can be computed accurately with single precision.
[40a87fa]940    #: This is True by default, but models such as :ref:`bcc-paracrystal` set
[a5b8477]941    #: it to False because they require double precision calculations.
[a86f6c0]942    single = None           # type: bool
[bb4b509]943    #: True if the model can be run as an opencl model.  If for some reason
944    #: the model cannot be run in opencl (e.g., because the model passes
945    #: functions by reference), then set this to false.
946    opencl = None           # type: bool
[a5b8477]947    #: True if the model is a structure factor used to model the interaction
948    #: between form factor models.  This will default to False if it is not
949    #: provided in the file.
[a86f6c0]950    structure_factor = None # type: bool
[c036ddb]951    #: True if the model defines an Fq function with signature
952    #: void Fq(double q, double *F1, double *F2, ...)
953    have_Fq = False
[e44432d]954    #: List of options for computing the effective radius of the shape,
955    #: or None if the model is not usable as a form factor model.
956    effective_radius_type = None   # type: List[str]
[a5b8477]957    #: List of C source files used to define the model.  The source files
[108e70e]958    #: should define the *Iq* function, and possibly *Iqac* or *Iqabc* if the
959    #: model defines orientation parameters. Files containing the most basic
960    #: functions must appear first in the list, followed by the files that
[e44432d]961    #: use those functions.
[a86f6c0]962    source = None           # type: List[str]
[39a06c9]963    #: inline source code, added after all elements of source
964    c_code = None           # type: Optional[str]
[a5b8477]965    #: Returns the form volume for python-based models.  Form volume is needed
966    #: for volume normalization in the polydispersity integral.  If no
967    #: parameters are *volume* parameters, then form volume is not needed.
968    #: For C-based models, (with :attr:`sources` defined, or with :attr:`Iq`
969    #: defined using a string containing C code), form_volume must also be
970    #: C code, either defined as a string, or in the sources.
[f619de7]971    form_volume = None      # type: Union[None, str, Callable[[np.ndarray], float]]
[e44432d]972    #: Returns the shell volume for python-based models.  Form volume and
973    #: shell volume are needed for volume normalization in the polydispersity
974    #: integral and structure interactions for hollow shapes.  If no
975    #: parameters are *volume* parameters, then shell volume is not needed.
976    #: For C-based models, (with :attr:`sources` defined, or with :attr:`Iq`
977    #: defined using a string containing C code), shell_volume must also be
978    #: C code, either defined as a string, or in the sources.
979    shell_volume = None      # type: Union[None, str, Callable[[np.ndarray], float]]
[a5b8477]980    #: Returns *I(q, a, b, ...)* for parameters *a*, *b*, etc. defined
981    #: by the parameter table.  *Iq* can be defined as a python function, or
982    #: as a C function.  If it is defined in C, then set *Iq* to the body of
983    #: the C function, including the return statement.  This function takes
984    #: values for *q* and each of the parameters as separate *double* values
985    #: (which may be converted to float or long double by sasmodels).  All
986    #: source code files listed in :attr:`sources` will be loaded before the
987    #: *Iq* function is defined.  If *Iq* is not present, then sources should
988    #: define *static double Iq(double q, double a, double b, ...)* which
989    #: will return *I(q, a, b, ...)*.  Multiplicity parameters are sent as
990    #: pointers to doubles.  Constants in floating point expressions should
991    #: include the decimal point. See :mod:`generate` for more details.
[f619de7]992    Iq = None               # type: Union[None, str, Callable[[np.ndarray], np.ndarray]]
[108e70e]993    #: Returns *I(qab, qc, a, b, ...)*.  The interface follows :attr:`Iq`.
994    Iqac = None             # type: Union[None, str, Callable[[np.ndarray], np.ndarray]]
995    #: Returns *I(qa, qb, qc, a, b, ...)*.  The interface follows :attr:`Iq`.
996    Iqabc = None            # type: Union[None, str, Callable[[np.ndarray], np.ndarray]]
[32e3c9b]997    #: Returns *I(qx, qy, a, b, ...)*.  The interface follows :attr:`Iq`.
998    Imagnetic = None        # type: Union[None, str, Callable[[np.ndarray], np.ndarray]]
[a5b8477]999    #: Returns a model profile curve *x, y*.  If *profile* is defined, this
1000    #: curve will appear in response to the *Show* button in SasView.  Use
1001    #: :attr:`profile_axes` to set the axis labels.  Note that *y* values
1002    #: will be scaled by 1e6 before plotting.
[c1a888b]1003    profile = None          # type: Optional[Callable[[np.ndarray], None]]
[a5b8477]1004    #: Axis labels for the :attr:`profile` plot.  The default is *['x', 'y']*.
1005    #: Only the *x* component is used for now.
1006    profile_axes = None     # type: Tuple[str, str]
1007    #: Returns *sesans(z, a, b, ...)* for models which can directly compute
1008    #: the SESANS correlation function.  Note: not currently implemented.
[c1a888b]1009    sesans = None           # type: Optional[Callable[[np.ndarray], np.ndarray]]
[e65c3ba]1010    #: Returns a random parameter set for the model
1011    random = None           # type: Optional[Callable[[], Dict[str, float]]]
[108e70e]1012    #: Line numbers for symbols defining C code
1013    lineno = None           # type: Dict[str, int]
[39a06c9]1014    #: The set of tests that must pass.  The format of the tests is described
1015    #: in :mod:`model_test`.
1016    tests = None            # type: List[TestCondition]
[da63656]1017
[6d6508e]1018    def __init__(self):
[a86f6c0]1019        # type: () -> None
[6d6508e]1020        pass
1021
[04045f4]1022    def get_hidden_parameters(self, control):
[a5b8477]1023        """
1024        Returns the set of hidden parameters for the model.  *control* is the
1025        value of the control parameter.  Note that multiplicity models have
1026        an implicit control parameter, which is the parameter that controls
1027        the multiplicity.
1028        """
[04045f4]1029        if self.hidden is not None:
1030            hidden = self.hidden(control)
1031        else:
[ce176ca]1032            controls = [p for p in self.parameters.kernel_parameters
1033                        if p.is_control]
[04045f4]1034            if len(controls) != 1:
1035                raise ValueError("more than one control parameter")
1036            hidden = set(p.id+str(k)
1037                         for p in self.parameters.kernel_parameters
1038                         for k in range(control+1, p.length+1)
1039                         if p.length > 1)
[bd547d0]1040            for p in self.parameters.kernel_parameters:
1041                if p.length > 1 and p.type == "sld":
1042                    for k in range(control+1, p.length+1):
1043                        base = p.id+str(k)
1044                        hidden.update((base+"_M0", base+"_mtheta", base+"_mphi"))
[04045f4]1045        return hidden
Note: See TracBrowser for help on using the repository browser.