source: sasmodels/sasmodels/modelinfo.py @ 225bf94

Last change on this file since 225bf94 was 225bf94, checked in by GitHub <noreply@…>, 7 months ago

Merge branch 'beta_approx' into master

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