source: sasmodels/sasmodels/modelinfo.py @ 6e7ff6d

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

reenable python models
  • Property mode set to 100644
File size: 20.2 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. The volume parameters are used for calls
192    to form_volume within the kernel (required for volume normalization)
193    and for calls to ER and VR for effective radius and volume ratio
194    respectively.
195
196    *description* is a short description of the parameter.  This will
197    be displayed in the parameter table and used as a tool tip for the
198    parameter value in the user interface.
199
200    Additional values can be set after the parameter is created:
201
202    * *length* is the length of the field if it is a vector field
203
204    * *length_control* is the parameter which sets the vector length
205
206    * *is_control* is True if the parameter is a control parameter for a vector
207
208    * *polydisperse* is true if the parameter accepts a polydispersity
209
210    * *relative_pd* is true if that polydispersity is a portion of the
211    value (so a 10% length dipsersity would use a polydispersity value of 0.1)
212    rather than absolute dispersisity (such as an angle plus or minus
213    15 degrees).
214
215    In the usual process these values are set by :func:`make_parameter_table`
216    and :func:`parse_parameter` therein.
217    """
218    def __init__(self, name, units='', default=None, limits=(-np.inf, np.inf),
219                 type='', description=''):
220        self.id = name.split('[')[0].strip()
221        self.name = name
222        self.units = units
223        self.default = default
224        self.limits = limits
225        self.type = type
226        self.description = description
227
228        # Length and length_control will be filled in once the complete
229        # parameter table is available.
230        self.length = 1
231        self.length_control = None
232        self.is_control = False
233
234        # TODO: need better control over whether a parameter is polydisperse
235        self.polydisperse = False
236        self.relative_pd = False
237
238        # choices are also set externally.
239        self.choices = []
240
241    def as_definition(self):
242        """
243        Declare space for the variable in a parameter structure.
244
245        For example, the parameter thickness with length 3 will
246        return "double thickness[3];", with no spaces before and
247        no new line character afterward.
248        """
249        if self.length == 1:
250            return "double %s;"%self.id
251        else:
252            return "double %s[%d];"%(self.id, self.length)
253
254    def as_function_argument(self):
255        """
256        Declare the variable as a function argument.
257
258        For example, the parameter thickness with length 3 will
259        return "double *thickness", with no spaces before and
260        no comma afterward.
261        """
262        if self.length == 1:
263            return "double %s"%self.id
264        else:
265            return "double *%s"%self.id
266
267    def as_call_reference(self, prefix=""):
268        # Note: if the parameter is a struct type, then we will need to use
269        # &prefix+id.  For scalars and vectors we can just use prefix+id.
270        return prefix + self.id
271
272    def __str__(self):
273        return "<%s>"%self.name
274
275    def __repr__(self):
276        return "P<%s>"%self.name
277
278
279class ParameterTable(object):
280    """
281    ParameterTable manages the list of available parameters.
282
283    There are a couple of complications which mean that the list of parameters
284    for the kernel differs from the list of parameters that the user sees.
285
286    (1) Common parameters.  Scale and background are implicit to every model,
287    but are not passed to the kernel.
288
289    (2) Vector parameters.  Vector parameters are passed to the kernel as a
290    pointer to an array, e.g., thick[], but they are seen by the user as n
291    separate parameters thick1, thick2, ...
292
293    Therefore, the parameter table is organized by how it is expected to be
294    used. The following information is needed to set up the kernel functions:
295
296    * *kernel_parameters* is the list of parameters in the kernel parameter
297    table, with vector parameter p declared as p[].
298
299    * *iq_parameters* is the list of parameters to the Iq(q, ...) function,
300    with vector parameter p sent as p[].
301
302    * *iqxy_parameters* is the list of parameters to the Iqxy(qx, qy, ...)
303    function, with vector parameter p sent as p[].
304
305    * *form_volume_parameters* is the list of parameters to the form_volume(...)
306    function, with vector parameter p sent as p[].
307
308    Problem details, which sets up the polydispersity loops, requires the
309    following:
310
311    * *theta_offset* is the offset of the theta parameter in the kernel parameter
312    table, with vector parameters counted as n individual parameters
313    p1, p2, ..., or offset is -1 if there is no theta parameter.
314
315    * *max_pd* is the maximum number of polydisperse parameters, with vector
316    parameters counted as n individual parameters p1, p2, ...  Note that
317    this number is limited to sasmodels.modelinfo.MAX_PD.
318
319    * *npars* is the total number of parameters to the kernel, with vector
320    parameters counted as n individual parameters p1, p2, ...
321
322    * *call_parameters* is the complete list of parameters to the kernel,
323    including scale and background, with vector parameters recorded as
324    individual parameters p1, p2, ...
325
326    * *active_1d* is the set of names that may be polydisperse for 1d data
327
328    * *active_2d* is the set of names that may be polydisperse for 2d data
329
330    User parameters are the set of parameters visible to the user, including
331    the scale and background parameters that the kernel does not see.  User
332    parameters don't use vector notation, and instead use p1, p2, ...
333
334    * *control_parameters* is the
335
336    """
337    # scale and background are implicit parameters
338    COMMON = [Parameter(*p) for p in COMMON_PARAMETERS]
339
340    def __init__(self, parameters):
341        self.kernel_parameters = parameters
342        self._set_vector_lengths()
343        self._make_call_parameter_list()
344        self._categorize_parameters()
345        self._set_defaults()
346        #self._name_table= dict((p.id, p) for p in parameters)
347
348    def _set_vector_lengths(self):
349        # Sort out the length of the vector parameters such as thickness[n]
350        for p in self.kernel_parameters:
351            if p.length_control:
352                for ref in self.kernel_parameters:
353                    if ref.id == p.length_control:
354                        break
355                else:
356                    raise ValueError("no reference variable %r for %s"
357                                     % (p.length_control, p.name))
358                ref.is_control = True
359                low, high = ref.limits
360                if int(low) != low or int(high) != high or low<0 or high>20:
361                    raise ValueError("expected limits on %s to be within [0, 20]"
362                                     % ref.name)
363                p.length = high
364
365    def _set_defaults(self):
366        # Construct default values, including vector defaults
367        defaults = {}
368        for p in self.call_parameters:
369            if p.length == 1:
370                defaults[p.id] = p.default
371            else:
372                for k in range(1, p.length+1):
373                    defaults["%s%d"%(p.id, k)] = p.default
374        self.defaults = defaults
375
376    def _make_call_parameter_list(self):
377        full_list = self.COMMON[:]
378        for p in self.kernel_parameters:
379            if p.length == 1:
380                full_list.append(p)
381            else:
382                for k in range(1, p.length+1):
383                    pk = Parameter(p.id+str(k), p.units, p.default,
384                                   p.limits, p.type, p.description)
385                    pk.polydisperse = p.polydisperse
386                    pk.relative_pd = p.relative_pd
387                    full_list.append(pk)
388        self.call_parameters = full_list
389
390    def _categorize_parameters(self):
391        # Set the kernel parameters.  Assumes background and scale are the
392        # first two parameters in the parameter list, but these are not sent
393        # to the underlying kernel functions.
394        self.iq_parameters = [p for p in self.kernel_parameters
395                              if p.type not in ('orientation', 'magnetic')]
396        self.iqxy_parameters = [p for p in self.kernel_parameters
397                                if p.type != 'magnetic']
398        self.form_volume_parameters = [p for p in self.kernel_parameters
399                                       if p.type == 'volume']
400
401        # Theta offset
402        offset = 0
403        for p in self.kernel_parameters:
404            if p.name == 'theta':
405                self.theta_offset = offset
406                break
407            offset += p.length
408        else:
409            self.theta_offset = -1
410
411        # number of polydisperse parameters
412        num_pd = sum(p.length for p in self.kernel_parameters if p.polydisperse)
413        # Don't use more polydisperse parameters than are available in the model
414        # Note: we can do polydispersity on arbitrary parameters, so it is not
415        # clear that this is a good idea; it does however make the poly_details
416        # code easier to write, so we will leave it in for now.
417        self.max_pd = min(num_pd, MAX_PD)
418
419        self.npars = sum(p.length for p in self.kernel_parameters)
420
421        # true if has 2D parameters
422        self.has_2d = any(p.type in ('orientation', 'magnetic')
423                          for p in self.kernel_parameters)
424
425        self.pd_1d = set(p.name for p in self.call_parameters
426                if p.polydisperse and p.type not in ('orientation', 'magnetic'))
427        self.pd_2d = set(p.name for p in self.call_parameters
428                         if p.polydisperse and p.type != 'magnetic')
429
430    def user_parameters(self, pars={}, is2d=True):
431        """
432        Return the list of parameters for the given data type.
433
434        Vector parameters are expanded as in place.  If multiple parameters
435        share the same vector length, then the parameters will be interleaved
436        in the result.  The control parameters come first.  For example,
437        if the parameter table is ordered as::
438
439            sld_core
440            sld_shell[num_shells]
441            sld_solvent
442            thickness[num_shells]
443            num_shells
444
445        and *pars[num_shells]=2* then the returned list will be::
446
447            num_shells
448            scale
449            background
450            sld_core
451            sld_shell1
452            thickness1
453            sld_shell2
454            thickness2
455            sld_solvent
456
457        Note that shell/thickness pairs are grouped together in the result
458        even though they were not grouped in the incoming table.  The control
459        parameter is always returned first since the GUI will want to set it
460        early, and rerender the table when it is changed.
461        """
462        control = [p for p in self.kernel_parameters if p.is_control]
463
464        # Gather entries such as name[n] into groups of the same n
465        dependent = dict((p.id, []) for p in control)
466        for p in self.kernel_parameters:
467            if p.length_control is not None:
468                dependent[p.length_control].append(p)
469
470        # Gather entries such as name[4] into groups of the same length
471        fixed = {}
472        for p in self.kernel_parameters:
473            if p.length > 1 and p.length_control is None:
474                fixed.setdefault(p.length, []).append(p)
475
476        # Using the call_parameters table, we already have expanded forms
477        # for each of the vector parameters; put them in a lookup table
478        expanded_pars = dict((p.name, p) for p in self.call_parameters)
479
480        # Gather the user parameters in order
481        result = control + self.COMMON
482        for p in self.kernel_parameters:
483            if not is2d and p.type in ('orientation', 'magnetic'):
484                pass
485            elif p.is_control:
486                pass # already added
487            elif p.length_control is not None:
488                table = dependent.get(p.length_control, [])
489                if table:
490                    # look up length from incoming parameters
491                    table_length = int(pars.get(p.length_control, p.length))
492                    del dependent[p.length_control] # first entry seen
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            elif p.length > 1:
499                table = fixed.get(p.length, [])
500                if table:
501                    table_length = p.length
502                    del fixed[p.length]
503                    for k in range(1, table_length+1):
504                        for entry in table:
505                            result.append(expanded_pars[entry.id+str(k)])
506                else:
507                    pass # already processed all entries
508            else:
509                result.append(p)
510
511        return result
512
513
514
515def isstr(x):
516    # TODO: 2-3 compatible tests for str, including unicode strings
517    return isinstance(x, str)
518
Note: See TracBrowser for help on using the repository browser.