source: sasmodels/sasmodels/modelinfo.py @ ee8f734

core_shell_microgelscostrafo411magnetic_modelrelease_v0.94release_v0.95ticket-1257-vesicle-productticket_1156ticket_1265_superballticket_822_more_unit_tests
Last change on this file since ee8f734 was ee8f734, checked in by Paul Kienzle <pkienzle@…>, 8 years ago

only use bare except when annotating an exception
  • Property mode set to 100644
File size: 19.8 KB
Line 
1
2import numpy as np
3
4# TODO: turn ModelInfo into a proper class
5ModelInfo = dict
6
7MAX_PD = 4
8
9COMMON_PARAMETERS = [
10    ["scale", "", 1, [0, np.inf], "", "Source intensity"],
11    ["background", "1/cm", 1e-3, [0, np.inf], "", "Source background"],
12]
13assert (len(COMMON_PARAMETERS) == 2
14        and COMMON_PARAMETERS[0][0]=="scale"
15        and COMMON_PARAMETERS[1][0]=="background"), "don't change common parameters"
16# assumptions about common parameters exist throughout the code, such as:
17# (1) kernel functions Iq, Iqxy, form_volume, ... don't see them
18# (2) kernel drivers assume scale is par[0] and background is par[1]
19# (3) mixture models drop the background on components and replace the scale
20#     with a scale that varies from [-inf, inf]
21# (4) product models drop the background and reassign scale
22# and maybe other places.
23# Note that scale and background cannot be coordinated parameters whose value
24# depends on the some polydisperse parameter with the current implementation
25
26def make_parameter_table(pars):
27    processed = []
28    for p in pars:
29        if not isinstance(p, list) or len(p) != 6:
30            raise ValueError("Parameter should be [name, units, default, limits, type, desc], but got %r"
31                             %str(p))
32        processed.append(parse_parameter(*p))
33    partable = ParameterTable(processed)
34    return partable
35
36def parse_parameter(name, units='', default=None,
37                    limits=(-np.inf, np.inf), type='', description=''):
38    # Parameter is a user facing class.  Do robust type checking.
39    if not isstr(name):
40        raise ValueError("expected string for parameter name %r"%name)
41    if not isstr(units):
42        raise ValueError("expected units to be a string for %s"%name)
43    # if limits is a list of strings, then this is a choice list
44    # field, and limits are 1 to length of string list
45    if isinstance(limits, list) and all(isstr(k) for k in limits):
46        choices = limits
47        limits = [0, len(choices)-1]
48    else:
49        choices = []
50    # TODO: maybe allow limits of None for (-inf, inf)
51    try:
52        low, high = limits
53        if not isinstance(low, (int, float)):
54            raise TypeError("low is not numeric")
55        if not isinstance(high, (int, float)):
56            raise TypeError("high is not numeric")
57        if low >= high:
58            raise ValueError("require low < high")
59    except Exception:
60        raise ValueError("invalid limits %s for %s"%(limits, name))
61
62    if not isinstance(default, (int, float)):
63        raise ValueError("expected default %r to be a number for %s"
64                         % (default, name))
65    if default < low or default > high:
66        raise ValueError("default value %r not in range for %s"
67                         % (default, name))
68
69    if type not in ("volume", "orientation", "sld", "magnetic", ""):
70        raise ValueError("unexpected type %r for %s" % (type, name))
71
72    if not isstr(description):
73        raise ValueError("expected description to be a string")
74
75
76    # Parameter id for name[n] does not include [n]
77    if "[" in name:
78        if not name.endswith(']'):
79            raise ValueError("Expected name[len] for vector parameter %s"%name)
80        pid, ref = name[:-1].split('[', 1)
81        ref = ref.strip()
82    else:
83        pid, ref = name, None
84
85
86    # automatically identify sld types
87    if type=='' and (pid.startswith('sld') or pid.endswith('sld')):
88        type = 'sld'
89
90    # Check if using a vector definition, name[k], as the parameter name
91    if ref:
92        if ref == '':
93            raise ValueError("Need to specify vector length for %s"%name)
94        try:
95            length = int(ref)
96            control = None
97        except Exception:
98            length = None
99            control = ref
100    else:
101        length = 1
102        control = None
103
104    # Build the parameter
105    parameter = Parameter(name=name, units=units, default=default,
106                          limits=limits, type=type, description=description)
107
108    # TODO: need better control over whether a parameter is polydisperse
109    parameter.polydisperse = type in ('orientation', 'volume')
110    parameter.relative_pd = type in ('volume')
111    parameter.choices = choices
112    parameter.length = length
113    parameter.length_control = control
114
115    return parameter
116
117
118def set_demo(model_info, demo):
119    """
120    Assign demo parameters to model_info['demo']
121
122    If demo is not defined, then use defaults.
123
124    If demo is defined, override defaults with value from demo.
125
126    If demo references vector fields, such as thickness[n], then support
127    different ways of assigning the demo values, including assigning a
128    specific value (e.g., thickness3=50.0), assigning a new value to all
129    (e.g., thickness=50.0) or assigning values using list notation.
130    """
131    partable = model_info['parameters']
132    if demo is None:
133        result = partable.defaults
134    else:
135        pars = dict((p.id, p) for p in partable.kernel_parameters)
136        result = partable.defaults.copy()
137        vectors = dict((name,value) for name,value in demo.items()
138                       if name in pars and pars[name].length > 1)
139        if vectors:
140            scalars = dict((name, value) for name, value in demo.items()
141                           if name not in pars or pars[name].length == 1)
142            for name, value in vectors.items():
143                if np.isscalar(value):
144                    # support for the form
145                    #    demo(thickness=0, thickness2=50)
146                    for k in range(1, pars[name].length+1):
147                        key = name+str(k)
148                        if key not in scalars:
149                            scalars[key] = vectors
150                else:
151                    # supoprt for the form
152                    #    demo(thickness=[20,10,3])
153                    for (k,v) in enumerate(value):
154                        scalars[name+str(k)] = v
155            result.update(scalars)
156        else:
157            result.update(demo)
158
159    model_info['demo'] = result
160
161class Parameter(object):
162    """
163    The available kernel parameters are defined as a list, with each parameter
164    defined as a sublist with the following elements:
165
166    *name* is the name that will be used in the call to the kernel
167    function and the name that will be displayed to the user.  Names
168    should be lower case, with words separated by underscore.  If
169    acronyms are used, the whole acronym should be upper case.
170
171    *units* should be one of *degrees* for angles, *Ang* for lengths,
172    *1e-6/Ang^2* for SLDs.
173
174    *default value* will be the initial value for  the model when it
175    is selected, or when an initial value is not otherwise specified.
176
177    *limits = [lb, ub]* are the hard limits on the parameter value, used to
178    limit the polydispersity density function.  In the fit, the parameter limits
179    given to the fit are the limits  on the central value of the parameter.
180    If there is polydispersity, it will evaluate parameter values outside
181    the fit limits, but not outside the hard limits specified in the model.
182    If there are no limits, use +/-inf imported from numpy.
183
184    *type* indicates how the parameter will be used.  "volume" parameters
185    will be used in all functions.  "orientation" parameters will be used
186    in *Iqxy* and *Imagnetic*.  "magnetic* parameters will be used in
187    *Imagnetic* only.  If *type* is the empty string, the parameter will
188    be used in all of *Iq*, *Iqxy* and *Imagnetic*.  "sld" parameters
189    can automatically be promoted to magnetic parameters, each of which
190    will have a magnitude and a direction, which may be different from
191    other sld parameters.
192
193    *description* is a short description of the parameter.  This will
194    be displayed in the parameter table and used as a tool tip for the
195    parameter value in the user interface.
196
197    Additional values can be set after the parameter is created:
198
199    *length* is the length of the field if it is a vector field
200    *length_control* is the parameter which sets the vector length
201    *is_control* is True if the parameter is a control parameter for a vector
202    *polydisperse* is true if the parameter accepts a polydispersity
203    *relative_pd* is true if that polydispersity is relative
204
205    In the usual process these values are set by :func:`make_parameter_table`
206    and :func:`parse_parameter` therein.
207    """
208    def __init__(self, name, units='', default=None, limits=(-np.inf, np.inf),
209                 type='', description=''):
210        self.id = name.split('[')[0].strip()
211        self.name = name
212        self.units = units
213        self.default = default
214        self.limits = limits
215        self.type = type
216        self.description = description
217
218        # Length and length_control will be filled in once the complete
219        # parameter table is available.
220        self.length = 1
221        self.length_control = None
222        self.is_control = False
223
224        # TODO: need better control over whether a parameter is polydisperse
225        self.polydisperse = False
226        self.relative_pd = False
227
228        # choices are also set externally.
229        self.choices = []
230
231    def as_definition(self):
232        """
233        Declare space for the variable in a parameter structure.
234
235        For example, the parameter thickness with length 3 will
236        return "double thickness[3];", with no spaces before and
237        no new line character afterward.
238        """
239        if self.length == 1:
240            return "double %s;"%self.id
241        else:
242            return "double %s[%d];"%(self.id, self.length)
243
244    def as_function_argument(self):
245        """
246        Declare the variable as a function argument.
247
248        For example, the parameter thickness with length 3 will
249        return "double *thickness", with no spaces before and
250        no comma afterward.
251        """
252        if self.length == 1:
253            return "double %s"%self.id
254        else:
255            return "double *%s"%self.id
256
257    def as_call_reference(self, prefix=""):
258        # Note: if the parameter is a struct type, then we will need to use
259        # &prefix+id.  For scalars and vectors we can just use prefix+id.
260        return prefix + self.id
261
262    def __str__(self):
263        return "<%s>"%self.name
264
265    def __repr__(self):
266        return "P<%s>"%self.name
267
268
269class ParameterTable(object):
270    """
271    ParameterTable manages the list of available parameters.
272
273    There are a couple of complications which mean that the list of parameters
274    for the kernel differs from the list of parameters that the user sees.
275
276    (1) Common parameters.  Scale and background are implicit to every model,
277    but are not passed to the kernel.
278
279    (2) Vector parameters.  Vector parameters are passed to the kernel as a
280    pointer to an array, e.g., thick[], but they are seen by the user as n
281    separate parameters thick1, thick2, ...
282
283    Therefore, the parameter table is organized by how it is expected to be
284    used. The following information is needed to set up the kernel functions:
285
286    * *kernel_parameters* is the list of parameters in the kernel parameter
287    table, with vector parameter p declared as p[].
288
289    * *iq_parameters* is the list of parameters to the Iq(q, ...) function,
290    with vector parameter p sent as p[].
291
292    * *iqxy_parameters* is the list of parameters to the Iqxy(qx, qy, ...)
293    function, with vector parameter p sent as p[].
294
295    * *form_volume_parameters* is the list of parameters to the form_volume(...)
296    function, with vector parameter p sent as p[].
297
298    Problem details, which sets up the polydispersity loops, requires the
299    following:
300
301    * *theta_offset* is the offset of the theta parameter in the kernel parameter
302    table, with vector parameters counted as n individual parameters
303    p1, p2, ..., or offset is -1 if there is no theta parameter.
304
305    * *max_pd* is the maximum number of polydisperse parameters, with vector
306    parameters counted as n individual parameters p1, p2, ...  Note that
307    this number is limited to sasmodels.modelinfo.MAX_PD.
308
309    * *npars* is the total number of parameters to the kernel, with vector
310    parameters counted as n individual parameters p1, p2, ...
311
312    * *call_parameters* is the complete list of parameters to the kernel,
313    including scale and background, with vector parameters recorded as
314    individual parameters p1, p2, ...
315
316    * *active_1d* is the set of names that may be polydisperse for 1d data
317
318    * *active_2d* is the set of names that may be polydisperse for 2d data
319
320    User parameters are the set of parameters visible to the user, including
321    the scale and background parameters that the kernel does not see.  User
322    parameters don't use vector notation, and instead use p1, p2, ...
323
324    * *control_parameters* is the
325
326    """
327    # scale and background are implicit parameters
328    COMMON = [Parameter(*p) for p in COMMON_PARAMETERS]
329
330    def __init__(self, parameters):
331        self.kernel_parameters = parameters
332        self._set_vector_lengths()
333        self._make_call_parameter_list()
334        self._categorize_parameters()
335        self._set_defaults()
336        #self._name_table= dict((p.id, p) for p in parameters)
337
338    def _set_vector_lengths(self):
339        # Sort out the length of the vector parameters such as thickness[n]
340        for p in self.kernel_parameters:
341            if p.length_control:
342                for ref in self.kernel_parameters:
343                    if ref.id == p.length_control:
344                        break
345                else:
346                    raise ValueError("no reference variable %r for %s"
347                                     % (p.length_control, p.name))
348                ref.is_control = True
349                low, high = ref.limits
350                if int(low) != low or int(high) != high or low<0 or high>20:
351                    raise ValueError("expected limits on %s to be within [0, 20]"
352                                     % ref.name)
353                p.length = high
354
355    def _set_defaults(self):
356        # Construct default values, including vector defaults
357        defaults = {}
358        for p in self.call_parameters:
359            if p.length == 1:
360                defaults[p.id] = p.default
361            else:
362                for k in range(1, p.length+1):
363                    defaults["%s%d"%(p.id, k)] = p.default
364        self.defaults = defaults
365
366    def _make_call_parameter_list(self):
367        full_list = self.COMMON[:]
368        for p in self.kernel_parameters:
369            if p.length == 1:
370                full_list.append(p)
371            else:
372                for k in range(1, p.length+1):
373                    pk = Parameter(p.id+str(k), p.units, p.default,
374                                   p.limits, p.type, p.description)
375                    pk.polydisperse = p.polydisperse
376                    pk.relative_pd = p.relative_pd
377                    full_list.append(pk)
378        self.call_parameters = full_list
379
380    def _categorize_parameters(self):
381        # Set the kernel parameters.  Assumes background and scale are the
382        # first two parameters in the parameter list, but these are not sent
383        # to the underlying kernel functions.
384        self.iq_parameters = [p for p in self.kernel_parameters
385                              if p.type not in ('orientation', 'magnetic')]
386        self.iqxy_parameters = [p for p in self.kernel_parameters
387                                if p.type != 'magnetic']
388        self.form_volume_parameters = [p for p in self.kernel_parameters
389                                       if p.type == 'volume']
390
391        # Theta offset
392        offset = 0
393        for p in self.kernel_parameters:
394            if p.name == 'theta':
395                self.theta_offset = offset
396                break
397            offset += p.length
398        else:
399            self.theta_offset = -1
400
401        # number of polydisperse parameters
402        num_pd = sum(p.length for p in self.kernel_parameters if p.polydisperse)
403        # Don't use more polydisperse parameters than are available in the model
404        # Note: we can do polydispersity on arbitrary parameters, so it is not
405        # clear that this is a good idea; it does however make the poly_details
406        # code easier to write, so we will leave it in for now.
407        self.max_pd = min(num_pd, MAX_PD)
408
409        self.npars = sum(p.length for p in self.kernel_parameters)
410
411        # true if has 2D parameters
412        self.has_2d = any(p.type in ('orientation', 'magnetic')
413                          for p in self.kernel_parameters)
414
415        self.pd_1d = set(p.name for p in self.call_parameters
416                if p.polydisperse and p.type not in ('orientation', 'magnetic'))
417        self.pd_2d = set(p.name for p in self.call_parameters
418                         if p.polydisperse and p.type != 'magnetic')
419
420    def user_parameters(self, pars={}, is2d=True):
421        """
422        Return the list of parameters for the given data type.
423
424        Vector parameters are expanded as in place.  If multiple parameters
425        share the same vector length, then the parameters will be interleaved
426        in the result.  The control parameters come first.  For example,
427        if the parameter table is ordered as::
428
429            sld_core
430            sld_shell[num_shells]
431            sld_solvent
432            thickness[num_shells]
433            num_shells
434
435        and *pars[num_shells]=2* then the returned list will be::
436
437            num_shells
438            scale
439            background
440            sld_core
441            sld_shell1
442            thickness1
443            sld_shell2
444            thickness2
445            sld_solvent
446
447        Note that shell/thickness pairs are grouped together in the result
448        even though they were not grouped in the incoming table.  The control
449        parameter is always returned first since the GUI will want to set it
450        early, and rerender the table when it is changed.
451        """
452        control = [p for p in self.kernel_parameters if p.is_control]
453
454        # Gather entries such as name[n] into groups of the same n
455        dependent = dict((p.id, []) for p in control)
456        for p in self.kernel_parameters:
457            if p.length_control is not None:
458                dependent[p.length_control].append(p)
459
460        # Gather entries such as name[4] into groups of the same length
461        fixed = {}
462        for p in self.kernel_parameters:
463            if p.length > 1 and p.length_control is None:
464                fixed.setdefault(p.length, []).append(p)
465
466        # Using the call_parameters table, we already have expanded forms
467        # for each of the vector parameters; put them in a lookup table
468        expanded_pars = dict((p.name, p) for p in self.call_parameters)
469
470        # Gather the user parameters in order
471        result = control + self.COMMON
472        for p in self.kernel_parameters:
473            if not is2d and p.type in ('orientation', 'magnetic'):
474                pass
475            elif p.is_control:
476                pass # already added
477            elif p.length_control is not None:
478                table = dependent.get(p.length_control, [])
479                if table:
480                    # look up length from incoming parameters
481                    table_length = int(pars.get(p.length_control, p.length))
482                    del dependent[p.length_control] # first entry seen
483                    for k in range(1, table_length+1):
484                        for entry in table:
485                            result.append(expanded_pars[entry.id+str(k)])
486                else:
487                    pass # already processed all entries
488            elif p.length > 1:
489                table = fixed.get(p.length, [])
490                if table:
491                    table_length = p.length
492                    del fixed[p.length]
493                    for k in range(1, table_length+1):
494                        for entry in table:
495                            result.append(expanded_pars[entry.id+str(k)])
496                else:
497                    pass # already processed all entries
498            else:
499                result.append(p)
500
501        return result
502
503
504
505def isstr(x):
506    # TODO: 2-3 compatible tests for str, including unicode strings
507    return isinstance(x, str)
508
Note: See TracBrowser for help on using the repository browser.