source: sasmodels/sasmodels/generate.py @ ce896fd

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

improved handling of vector parameters; remove compile errors from onion.c

  • Property mode set to 100644
File size: 31.6 KB
Line 
1"""
2SAS model constructor.
3
4Small angle scattering models are defined by a set of kernel functions:
5
6    *Iq(q, p1, p2, ...)* returns the scattering at q for a form with
7    particular dimensions averaged over all orientations.
8
9    *Iqxy(qx, qy, p1, p2, ...)* returns the scattering at qx, qy for a form
10    with particular dimensions for a single orientation.
11
12    *Imagnetic(qx, qy, result[], p1, p2, ...)* returns the scattering for the
13    polarized neutron spin states (up-up, up-down, down-up, down-down) for
14    a form with particular dimensions for a single orientation.
15
16    *form_volume(p1, p2, ...)* returns the volume of the form with particular
17    dimension.
18
19    *ER(p1, p2, ...)* returns the effective radius of the form with
20    particular dimensions.
21
22    *VR(p1, p2, ...)* returns the volume ratio for core-shell style forms.
23
24    #define INVALID(v) (expr)  returns False if v.parameter is invalid
25    for some parameter or other (e.g., v.bell_radius < v.radius).  If
26    necessary, the expression can call a function.
27
28These functions are defined in a kernel module .py script and an associated
29set of .c files.  The model constructor will use them to create models with
30polydispersity across volume and orientation parameters, and provide
31scale and background parameters for each model.
32
33*Iq*, *Iqxy*, *Imagnetic* and *form_volume* should be stylized C-99
34functions written for OpenCL.  All functions need prototype declarations
35even if the are defined before they are used.  OpenCL does not support
36*#include* preprocessor directives, so instead the list of includes needs
37to be given as part of the metadata in the kernel module definition.
38The included files should be listed using a path relative to the kernel
39module, or if using "lib/file.c" if it is one of the standard includes
40provided with the sasmodels source.  The includes need to be listed in
41order so that functions are defined before they are used.
42
43Floating point values should be declared as *double*.  For single precision
44calculations, *double* will be replaced by *float*.  The single precision
45conversion will also tag floating point constants with "f" to make them
46single precision constants.  When using integral values in floating point
47expressions, they should be expressed as floating point values by including
48a decimal point.  This includes 0., 1. and 2.
49
50OpenCL has a *sincos* function which can improve performance when both
51the *sin* and *cos* values are needed for a particular argument.  Since
52this function does not exist in C99, all use of *sincos* should be
53replaced by the macro *SINCOS(value, sn, cn)* where *sn* and *cn* are
54previously declared *double* variables.  When compiled for systems without
55OpenCL, *SINCOS* will be replaced by *sin* and *cos* calls.   If *value* is
56an expression, it will appear twice in this case; whether or not it will be
57evaluated twice depends on the quality of the compiler.
58
59If the input parameters are invalid, the scattering calculator should
60return a negative number. Particularly with polydispersity, there are
61some sets of shape parameters which lead to nonsensical forms, such
62as a capped cylinder where the cap radius is smaller than the
63cylinder radius.  The polydispersity calculation will ignore these points,
64effectively chopping the parameter weight distributions at the boundary
65of the infeasible region.  The resulting scattering will be set to
66background.  This will work correctly even when polydispersity is off.
67
68*ER* and *VR* are python functions which operate on parameter vectors.
69The constructor code will generate the necessary vectors for computing
70them with the desired polydispersity.
71The kernel module must set variables defining the kernel meta data:
72
73    *id* is an implicit variable formed from the filename.  It will be
74    a valid python identifier, and will be used as the reference into
75    the html documentation, with '_' replaced by '-'.
76
77    *name* is the model name as displayed to the user.  If it is missing,
78    it will be constructed from the id.
79
80    *title* is a short description of the model, suitable for a tool tip,
81    or a one line model summary in a table of models.
82
83    *description* is an extended description of the model to be displayed
84    while the model parameters are being edited.
85
86    *parameters* is the list of parameters.  Parameters in the kernel
87    functions must appear in the same order as they appear in the
88    parameters list.  Two additional parameters, *scale* and *background*
89    are added to the beginning of the parameter list.  They will show up
90    in the documentation as model parameters, but they are never sent to
91    the kernel functions.  Note that *effect_radius* and *volfraction*
92    must occur first in structure factor calculations.
93
94    *category* is the default category for the model.  The category is
95    two level structure, with the form "group:section", indicating where
96    in the manual the model will be located.  Models are alphabetical
97    within their section.
98
99    *source* is the list of C-99 source files that must be joined to
100    create the OpenCL kernel functions.  The files defining the functions
101    need to be listed before the files which use the functions.
102
103    *ER* is a python function defining the effective radius.  If it is
104    not present, the effective radius is 0.
105
106    *VR* is a python function defining the volume ratio.  If it is not
107    present, the volume ratio is 1.
108
109    *form_volume*, *Iq*, *Iqxy*, *Imagnetic* are strings containing the
110    C source code for the body of the volume, Iq, and Iqxy functions
111    respectively.  These can also be defined in the last source file.
112
113    *Iq* and *Iqxy* also be instead be python functions defining the
114    kernel.  If they are marked as *Iq.vectorized = True* then the
115    kernel is passed the entire *q* vector at once, otherwise it is
116    passed values one *q* at a time.  The performance improvement of
117    this step is significant.
118
119    *demo* is a dictionary of parameter=value defining a set of
120    parameters to use by default when *compare* is called.  Any
121    parameter not set in *demo* gets the initial value from the
122    parameter list.  *demo* is mostly needed to set the default
123    polydispersity values for tests.
124
125    *oldname* is the name of the model in sasview before sasmodels
126    was split into its own package, and *oldpars* is a dictionary
127    of *parameter: old_parameter* pairs defining the new names for
128    the parameters.  This is used by *compare* to check the values
129    of the new model against the values of the old model before
130    you are ready to add the new model to sasmodels.
131
132
133An *model_info* dictionary is constructed from the kernel meta data and
134returned to the caller.
135
136The model evaluator, function call sequence consists of q inputs and the return vector,
137followed by the loop value/weight vector, followed by the values for
138the non-polydisperse parameters, followed by the lengths of the
139polydispersity loops.  To construct the call for 1D models, the
140categories *fixed-1d* and *pd-1d* list the names of the parameters
141of the non-polydisperse and the polydisperse parameters respectively.
142Similarly, *fixed-2d* and *pd-2d* provide parameter names for 2D models.
143The *pd-rel* category is a set of those parameters which give
144polydispersitiy as a portion of the value (so a 10% length dispersity
145would use a polydispersity value of 0.1) rather than absolute
146dispersity such as an angle plus or minus 15 degrees.
147
148The *volume* category lists the volume parameters in order for calls
149to volume within the kernel (used for volume normalization) and for
150calls to ER and VR for effective radius and volume ratio respectively.
151
152The *orientation* and *magnetic* categories list the orientation and
153magnetic parameters.  These are used by the sasview interface.  The
154blank category is for parameters such as scale which don't have any
155other marking.
156
157The doc string at the start of the kernel module will be used to
158construct the model documentation web pages.  Embedded figures should
159appear in the subdirectory "img" beside the model definition, and tagged
160with the kernel module name to avoid collision with other models.  Some
161file systems are case-sensitive, so only use lower case characters for
162file names and extensions.
163
164
165The function :func:`make` loads the metadata from the module and returns
166the kernel source.  The function :func:`make_doc` extracts the doc string
167and adds the parameter table to the top.  The function :func:`model_sources`
168returns a list of files required by the model.
169
170Code follows the C99 standard with the following extensions and conditions::
171
172    M_PI_180 = pi/180
173    M_4PI_3 = 4pi/3
174    square(x) = x*x
175    cube(x) = x*x*x
176    sinc(x) = sin(x)/x, with sin(0)/0 -> 1
177    all double precision constants must include the decimal point
178    all double declarations may be converted to half, float, or long double
179    FLOAT_SIZE is the number of bytes in the converted variables
180"""
181from __future__ import print_function
182
183#TODO: determine which functions are useful outside of generate
184#__all__ = ["model_info", "make_doc", "make_source", "convert_type"]
185
186from os.path import abspath, dirname, join as joinpath, exists, basename, \
187    splitext, getmtime
188import re
189import string
190import warnings
191
192import numpy as np
193
194from .modelinfo import ModelInfo, Parameter, make_parameter_table, set_demo
195
196# TODO: identify model files which have changed since loading and reload them.
197
198TEMPLATE_ROOT = dirname(__file__)
199
200F16 = np.dtype('float16')
201F32 = np.dtype('float32')
202F64 = np.dtype('float64')
203try:  # CRUFT: older numpy does not support float128
204    F128 = np.dtype('float128')
205except TypeError:
206    F128 = None
207
208# Conversion from units defined in the parameter table for each model
209# to units displayed in the sphinx documentation.
210RST_UNITS = {
211    "Ang": "|Ang|",
212    "1/Ang": "|Ang^-1|",
213    "1/Ang^2": "|Ang^-2|",
214    "1e-6/Ang^2": "|1e-6Ang^-2|",
215    "degrees": "degree",
216    "1/cm": "|cm^-1|",
217    "Ang/cm": "|Ang*cm^-1|",
218    "g/cm3": "|g/cm^3|",
219    "mg/m2": "|mg/m^2|",
220    "": "None",
221    }
222
223# Headers for the parameters tables in th sphinx documentation
224PARTABLE_HEADERS = [
225    "Parameter",
226    "Description",
227    "Units",
228    "Default value",
229    ]
230
231# Minimum width for a default value (this is shorter than the column header
232# width, so will be ignored).
233PARTABLE_VALUE_WIDTH = 10
234
235# Documentation header for the module, giving the model name, its short
236# description and its parameter table.  The remainder of the doc comes
237# from the module docstring.
238DOC_HEADER = """.. _%(id)s:
239
240%(name)s
241=======================================================
242
243%(title)s
244
245%(parameters)s
246
247%(returns)s
248
249%(docs)s
250"""
251
252def format_units(units):
253    """
254    Convert units into ReStructured Text format.
255    """
256    return "string" if isinstance(units, list) else RST_UNITS.get(units, units)
257
258def make_partable(pars):
259    """
260    Generate the parameter table to include in the sphinx documentation.
261    """
262    column_widths = [
263        max(len(p.name) for p in pars),
264        max(len(p.description) for p in pars),
265        max(len(format_units(p.units)) for p in pars),
266        PARTABLE_VALUE_WIDTH,
267        ]
268    column_widths = [max(w, len(h))
269                     for w, h in zip(column_widths, PARTABLE_HEADERS)]
270
271    sep = " ".join("="*w for w in column_widths)
272    lines = [
273        sep,
274        " ".join("%-*s" % (w, h)
275                 for w, h in zip(column_widths, PARTABLE_HEADERS)),
276        sep,
277        ]
278    for p in pars:
279        lines.append(" ".join([
280            "%-*s" % (column_widths[0], p.name),
281            "%-*s" % (column_widths[1], p.description),
282            "%-*s" % (column_widths[2], format_units(p.units)),
283            "%*g" % (column_widths[3], p.default),
284            ]))
285    lines.append(sep)
286    return "\n".join(lines)
287
288def _search(search_path, filename):
289    """
290    Find *filename* in *search_path*.
291
292    Raises ValueError if file does not exist.
293    """
294    for path in search_path:
295        target = joinpath(path, filename)
296        if exists(target):
297            return target
298    raise ValueError("%r not found in %s" % (filename, search_path))
299
300
301def model_sources(model_info):
302    """
303    Return a list of the sources file paths for the module.
304    """
305    search_path = [dirname(model_info['filename']),
306                   abspath(joinpath(dirname(__file__), 'models'))]
307    return [_search(search_path, f) for f in model_info['source']]
308
309def timestamp(model_info):
310    """
311    Return a timestamp for the model corresponding to the most recently
312    changed file or dependency.
313    """
314    source_files = (model_sources(model_info)
315                    + model_templates()
316                    + [model_info['filename']])
317    newest = max(getmtime(f) for f in source_files)
318    return newest
319
320def convert_type(source, dtype):
321    """
322    Convert code from double precision to the desired type.
323
324    Floating point constants are tagged with 'f' for single precision or 'L'
325    for long double precision.
326    """
327    if dtype == F16:
328        fbytes = 2
329        source = _convert_type(source, "float", "f")
330    elif dtype == F32:
331        fbytes = 4
332        source = _convert_type(source, "float", "f")
333    elif dtype == F64:
334        fbytes = 8
335        # no need to convert if it is already double
336    elif dtype == F128:
337        fbytes = 16
338        source = _convert_type(source, "long double", "L")
339    else:
340        raise ValueError("Unexpected dtype in source conversion: %s"%dtype)
341    return ("#define FLOAT_SIZE %d\n"%fbytes)+source
342
343
344def _convert_type(source, type_name, constant_flag):
345    """
346    Replace 'double' with *type_name* in *source*, tagging floating point
347    constants with *constant_flag*.
348    """
349    # Convert double keyword to float/long double/half.
350    # Accept an 'n' # parameter for vector # values, where n is 2, 4, 8 or 16.
351    # Assume complex numbers are represented as cdouble which is typedef'd
352    # to double2.
353    source = re.sub(r'(^|[^a-zA-Z0-9_]c?)double(([248]|16)?($|[^a-zA-Z0-9_]))',
354                    r'\1%s\2'%type_name, source)
355    # Convert floating point constants to single by adding 'f' to the end,
356    # or long double with an 'L' suffix.  OS/X complains if you don't do this.
357    source = re.sub(r'[^a-zA-Z_](\d*[.]\d+|\d+[.]\d*)([eE][+-]?\d+)?',
358                    r'\g<0>%s'%constant_flag, source)
359    return source
360
361
362def kernel_name(model_info, is_2d):
363    """
364    Name of the exported kernel symbol.
365    """
366    return model_info['name'] + "_" + ("Iqxy" if is_2d else "Iq")
367
368
369def indent(s, depth):
370    """
371    Indent a string of text with *depth* additional spaces on each line.
372    """
373    spaces = " "*depth
374    sep = "\n" + spaces
375    return spaces + sep.join(s.split("\n"))
376
377
378_template_cache = {}
379def load_template(filename):
380    path = joinpath(TEMPLATE_ROOT, filename)
381    mtime = getmtime(path)
382    if filename not in _template_cache or mtime > _template_cache[filename][0]:
383        with open(path) as fid:
384            _template_cache[filename] = (mtime, fid.read(), path)
385    return _template_cache[filename][1]
386
387def model_templates():
388    # TODO: fails DRY; templates are listed in two places.
389    # should instead have model_info contain a list of paths
390    return [joinpath(TEMPLATE_ROOT, filename)
391            for filename in ('kernel_header.c', 'kernel_iq.c')]
392
393
394_FN_TEMPLATE = """\
395double %(name)s(%(pars)s);
396double %(name)s(%(pars)s) {
397    %(body)s
398}
399
400
401"""
402
403def _gen_fn(name, pars, body):
404    """
405    Generate a function given pars and body.
406
407    Returns the following string::
408
409         double fn(double a, double b, ...);
410         double fn(double a, double b, ...) {
411             ....
412         }
413    """
414    par_decl = ', '.join(p.as_function_argument() for p in pars) if pars else 'void'
415    return _FN_TEMPLATE % {'name': name, 'body': body, 'pars': par_decl}
416
417def _call_pars(prefix, pars):
418    """
419    Return a list of *prefix.parameter* from parameter items.
420    """
421    return [p.as_call_reference(prefix) for p in pars]
422
423_IQXY_PATTERN = re.compile("^((inline|static) )? *(double )? *Iqxy *([(]|$)",
424                           flags=re.MULTILINE)
425def _have_Iqxy(sources):
426    """
427    Return true if any file defines Iqxy.
428
429    Note this is not a C parser, and so can be easily confused by
430    non-standard syntax.  Also, it will incorrectly identify the following
431    as having Iqxy::
432
433        /*
434        double Iqxy(qx, qy, ...) { ... fill this in later ... }
435        */
436
437    If you want to comment out an Iqxy function, use // on the front of the
438    line instead.
439    """
440    for code in sources:
441        if _IQXY_PATTERN.search(code):
442            return True
443    else:
444        return False
445
446def make_source(model_info):
447    """
448    Generate the OpenCL/ctypes kernel from the module info.
449
450    Uses source files found in the given search path.
451    """
452    if callable(model_info['Iq']):
453        return None
454
455    # TODO: need something other than volume to indicate dispersion parameters
456    # No volume normalization despite having a volume parameter.
457    # Thickness is labelled a volume in order to trigger polydispersity.
458    # May want a separate dispersion flag, or perhaps a separate category for
459    # disperse, but not volume.  Volume parameters also use relative values
460    # for the distribution rather than the absolute values used by angular
461    # dispersion.  Need to be careful that necessary parameters are available
462    # for computing volume even if we allow non-disperse volume parameters.
463
464    partable = model_info['parameters']
465
466    # Identify parameters for Iq, Iqxy, Iq_magnetic and form_volume.
467    # Note that scale and volume are not possible types.
468
469    # Load templates and user code
470    kernel_header = load_template('kernel_header.c')
471    kernel_code = load_template('kernel_iq.c')
472    user_code = [open(f).read() for f in model_sources(model_info)]
473
474    # Build initial sources
475    source = [kernel_header] + user_code
476
477    # Make parameters for q, qx, qy so that we can use them in declarations
478    q, qx, qy = [Parameter(name=v) for v in ('q', 'qx', 'qy')]
479    # Generate form_volume function, etc. from body only
480    if model_info['form_volume'] is not None:
481        pars = partable.form_volume_parameters
482        source.append(_gen_fn('form_volume', pars, model_info['form_volume']))
483    if model_info['Iq'] is not None:
484        pars = [q] + partable.iq_parameters
485        source.append(_gen_fn('Iq', pars, model_info['Iq']))
486    if model_info['Iqxy'] is not None:
487        pars = [qx, qy] + partable.iqxy_parameters
488        source.append(_gen_fn('Iqxy', pars, model_info['Iqxy']))
489
490    # Define the parameter table
491    source.append("#define PARAMETER_TABLE \\")
492    source.append("\\\n".join(p.as_definition()
493                              for p in partable.kernel_parameters))
494
495    # Define the function calls
496    if partable.form_volume_parameters:
497        refs = _call_pars("v.", partable.form_volume_parameters)
498        call_volume = "#define CALL_VOLUME(v) form_volume(%s)" % (",".join(refs))
499    else:
500        # Model doesn't have volume.  We could make the kernel run a little
501        # faster by not using/transferring the volume normalizations, but
502        # the ifdef's reduce readability more than is worthwhile.
503        call_volume = "#define CALL_VOLUME(v) 1.0"
504    source.append(call_volume)
505
506    refs = ["q[i]"] + _call_pars("v.", partable.iq_parameters)
507    call_iq = "#define CALL_IQ(q,i,v) Iq(%s)" % (",".join(refs))
508    if _have_Iqxy(user_code):
509        # Call 2D model
510        refs = ["q[2*i]", "q[2*i+1]"] + _call_pars("v.", partable.iqxy_parameters)
511        call_iqxy = "#define CALL_IQ(q,i,v) Iqxy(%s)" % (",".join(refs))
512    else:
513        # Call 1D model with sqrt(qx^2 + qy^2)
514        warnings.warn("Creating Iqxy = Iq(sqrt(qx^2 + qy^2))")
515        # still defined:: refs = ["q[i]"] + _call_pars("v", iq_parameters)
516        pars_sqrt = ["sqrt(q[2*i]*q[2*i]+q[2*i+1]*q[2*i+1])"] + refs[1:]
517        call_iqxy = "#define CALL_IQ(q,i,v) Iq(%s)" % (",".join(pars_sqrt))
518
519    # Fill in definitions for numbers of parameters
520    source.append("#define MAX_PD %s"%partable.max_pd)
521    source.append("#define NPARS %d"%partable.npars)
522
523    # TODO: allow mixed python/opencl kernels?
524
525    # define the Iq kernel
526    source.append("#define KERNEL_NAME %s_Iq"%model_info['name'])
527    source.append(call_iq)
528    source.append(kernel_code)
529    source.append("#undef CALL_IQ")
530    source.append("#undef KERNEL_NAME")
531
532    # define the Iqxy kernel from the same source with different #defines
533    source.append("#define KERNEL_NAME %s_Iqxy"%model_info['name'])
534    source.append(call_iqxy)
535    source.append(kernel_code)
536    source.append("#undef CALL_IQ")
537    source.append("#undef KERNEL_NAME")
538
539    return '\n'.join(source)
540
541class CoordinationDetails(object):
542    def __init__(self, model_info):
543        parameters = model_info['parameters']
544        max_pd = parameters.max_pd
545        npars = parameters.npars
546        par_offset = 4*max_pd
547        self.details = np.zeros(par_offset + 3*npars + 4, 'i4')
548
549        # generate views on different parts of the array
550        self._pd_par     = self.details[0*max_pd:1*max_pd]
551        self._pd_length  = self.details[1*max_pd:2*max_pd]
552        self._pd_offset  = self.details[2*max_pd:3*max_pd]
553        self._pd_stride  = self.details[3*max_pd:4*max_pd]
554        self._par_offset = self.details[par_offset+0*npars:par_offset+1*npars]
555        self._par_coord  = self.details[par_offset+1*npars:par_offset+2*npars]
556        self._pd_coord   = self.details[par_offset+2*npars:par_offset+3*npars]
557
558        # theta_par is fixed
559        self.details[-1] = parameters.theta_offset
560
561    @property
562    def ctypes(self): return self.details.ctypes
563    @property
564    def pd_par(self): return self._pd_par
565    @property
566    def pd_length(self): return self._pd_length
567    @property
568    def pd_offset(self): return self._pd_offset
569    @property
570    def pd_stride(self): return self._pd_stride
571    @property
572    def pd_coord(self): return self._pd_coord
573    @property
574    def par_coord(self): return self._par_coord
575    @property
576    def par_offset(self): return self._par_offset
577    @property
578    def num_coord(self): return self.details[-2]
579    @num_coord.setter
580    def num_coord(self, v): self.details[-2] = v
581    @property
582    def total_pd(self): return self.details[-3]
583    @total_pd.setter
584    def total_pd(self, v): self.details[-3] = v
585    @property
586    def num_active(self): return self.details[-4]
587    @num_active.setter
588    def num_active(self, v): self.details[-4] = v
589
590    def show(self):
591        print("total_pd", self.total_pd)
592        print("num_active", self.num_active)
593        print("pd_par", self.pd_par)
594        print("pd_length", self.pd_length)
595        print("pd_offset", self.pd_offset)
596        print("pd_stride", self.pd_stride)
597        print("par_offsets", self.par_offset)
598        print("num_coord", self.num_coord)
599        print("par_coord", self.par_coord)
600        print("pd_coord", self.pd_coord)
601        print("theta par", self.details[-1])
602
603def mono_details(model_info):
604    details = CoordinationDetails(model_info)
605    # The zero defaults for monodisperse systems are mostly fine
606    details.par_offset[:] = np.arange(2, len(details.par_offset)+2)
607    return details
608
609def poly_details(model_info, weights):
610    #print("weights",weights)
611    weights = weights[2:] # Skip scale and background
612
613    # Decreasing list of polydispersity lengths
614    # Note: the reversing view, x[::-1], does not require a copy
615    pd_length = np.array([len(w) for w in weights])
616    num_active = np.sum(pd_length>1)
617    if num_active > model_info['parameters'].max_pd:
618        raise ValueError("Too many polydisperse parameters")
619
620    pd_offset = np.cumsum(np.hstack((0, pd_length)))
621    idx = np.argsort(pd_length)[::-1][:num_active]
622    par_length = np.array([max(len(w),1) for w in weights])
623    pd_stride = np.cumprod(np.hstack((1, par_length[idx])))
624    par_offsets = np.cumsum(np.hstack((2, par_length)))
625
626    details = CoordinationDetails(model_info)
627    details.pd_par[:num_active] = idx
628    details.pd_length[:num_active] = pd_length[idx]
629    details.pd_offset[:num_active] = pd_offset[idx]
630    details.pd_stride[:num_active] = pd_stride[:-1]
631    details.par_offset[:] = par_offsets[:-1]
632    details.total_pd = pd_stride[-1]
633    details.num_active = num_active
634    # Without constraints coordinated parameters are just the pd parameters
635    details.par_coord[:num_active] = idx
636    details.pd_coord[:num_active] = 2**np.arange(num_active)
637    details.num_coord = num_active
638    #details.show()
639    return details
640
641def constrained_poly_details(model_info, weights, constraints):
642    # Need to find the independently varying pars and sort them
643    # Need to build a coordination list for the dependent variables
644    # Need to generate a constraints function which takes values
645    # and weights, returning par blocks
646    raise NotImplementedError("Can't handle constraints yet")
647
648
649def create_default_functions(model_info):
650    """
651    Autogenerate missing functions, such as Iqxy from Iq.
652
653    This only works for Iqxy when Iq is written in python. :func:`make_source`
654    performs a similar role for Iq written in C.
655    """
656    if callable(model_info['Iq']) and model_info['Iqxy'] is None:
657        partable = model_info['parameters']
658        if partable.has_2d:
659            raise ValueError("Iqxy model is missing")
660        Iq = model_info['Iq']
661        def Iqxy(qx, qy, **kw):
662            return Iq(np.sqrt(qx**2 + qy**2), **kw)
663        model_info['Iqxy'] = Iqxy
664
665
666def make_model_info(kernel_module):
667    """
668    Interpret the model definition file, categorizing the parameters.
669
670    The module can be loaded with a normal python import statement if you
671    know which module you need, or with __import__('sasmodels.model.'+name)
672    if the name is in a string.
673
674    The *model_info* structure contains the following fields:
675
676    * *id* is the id of the kernel
677    * *name* is the display name of the kernel
678    * *filename* is the full path to the module defining the file (if any)
679    * *title* is a short description of the kernel
680    * *description* is a long description of the kernel (this doesn't seem
681      very useful since the Help button on the model page brings you directly
682      to the documentation page)
683    * *docs* is the docstring from the module.  Use :func:`make_doc` to
684    * *category* specifies the model location in the docs
685    * *parameters* is the model parameter table
686    * *single* is True if the model allows single precision
687    * *structure_factor* is True if the model is useable in a product
688    * *variant_info* contains the information required to select between
689      model variants (e.g., the list of cases) or is None if there are no
690      model variants
691    * *par_type* categorizes the model parameters. See
692      :func:`categorize_parameters` for details.
693    * *demo* contains the *{parameter: value}* map used in compare (and maybe
694      for the demo plot, if plots aren't set up to use the default values).
695      If *demo* is not given in the file, then the default values will be used.
696    * *tests* is a set of tests that must pass
697    * *source* is the list of library files to include in the C model build
698    * *Iq*, *Iqxy*, *form_volume*, *ER*, *VR* and *sesans* are python functions
699      implementing the kernel for the module, or None if they are not
700      defined in python
701    * *oldname* is the model name in pre-4.0 Sasview
702    * *oldpars* is the *{new: old}* parameter translation table
703      from pre-4.0 Sasview
704    * *composition* is None if the model is independent, otherwise it is a
705      tuple with composition type ('product' or 'mixture') and a list of
706      *model_info* blocks for the composition objects.  This allows us to
707      build complete product and mixture models from just the info.
708    """
709    # TODO: maybe turn model_info into a class ModelDefinition
710    #print("make parameter table", kernel_module.parameters)
711    parameters = make_parameter_table(kernel_module.parameters)
712    filename = abspath(kernel_module.__file__)
713    kernel_id = splitext(basename(filename))[0]
714    name = getattr(kernel_module, 'name', None)
715    if name is None:
716        name = " ".join(w.capitalize() for w in kernel_id.split('_'))
717    model_info = dict(
718        id=kernel_id,  # string used to load the kernel
719        filename=abspath(kernel_module.__file__),
720        name=name,
721        title=getattr(kernel_module, 'title', name+" model"),
722        description=getattr(kernel_module, 'description', 'no description'),
723        parameters=parameters,
724        composition=None,
725        docs=kernel_module.__doc__,
726        category=getattr(kernel_module, 'category', None),
727        single=getattr(kernel_module, 'single', True),
728        structure_factor=getattr(kernel_module, 'structure_factor', False),
729        profile_axes=getattr(kernel_module, 'profile_axes', ['x','y']),
730        variant_info=getattr(kernel_module, 'invariant_info', None),
731        demo=getattr(kernel_module, 'demo', None),
732        source=getattr(kernel_module, 'source', []),
733        oldname=getattr(kernel_module, 'oldname', None),
734        oldpars=getattr(kernel_module, 'oldpars', {}),
735        tests=getattr(kernel_module, 'tests', []),
736        )
737    set_demo(model_info, getattr(kernel_module, 'demo', None))
738    # Check for optional functions
739    functions = "ER VR form_volume Iq Iqxy profile sesans".split()
740    model_info.update((k, getattr(kernel_module, k, None)) for k in functions)
741    create_default_functions(model_info)
742    # Precalculate the monodisperse parameters
743    # TODO: make this a lazy evaluator
744    # make_model_info is called for every model on sasview startup
745    model_info['mono_details'] = mono_details(model_info)
746    return model_info
747
748section_marker = re.compile(r'\A(?P<first>[%s])(?P=first)*\Z'
749                            %re.escape(string.punctuation))
750def _convert_section_titles_to_boldface(lines):
751    """
752    Do the actual work of identifying and converting section headings.
753    """
754    prior = None
755    for line in lines:
756        if prior is None:
757            prior = line
758        elif section_marker.match(line):
759            if len(line) >= len(prior):
760                yield "".join(("**", prior, "**"))
761                prior = None
762            else:
763                yield prior
764                prior = line
765        else:
766            yield prior
767            prior = line
768    if prior is not None:
769        yield prior
770
771def convert_section_titles_to_boldface(s):
772    """
773    Use explicit bold-face rather than section headings so that the table of
774    contents is not polluted with section names from the model documentation.
775
776    Sections are identified as the title line followed by a line of punctuation
777    at least as long as the title line.
778    """
779    return "\n".join(_convert_section_titles_to_boldface(s.split('\n')))
780
781def make_doc(model_info):
782    """
783    Return the documentation for the model.
784    """
785    Iq_units = "The returned value is scaled to units of |cm^-1| |sr^-1|, absolute scale."
786    Sq_units = "The returned value is a dimensionless structure factor, $S(q)$."
787    docs = convert_section_titles_to_boldface(model_info['docs'])
788    subst = dict(id=model_info['id'].replace('_', '-'),
789                 name=model_info['name'],
790                 title=model_info['title'],
791                 parameters=make_partable(model_info['parameters']),
792                 returns=Sq_units if model_info['structure_factor'] else Iq_units,
793                 docs=docs)
794    return DOC_HEADER % subst
795
796
797def demo_time():
798    """
799    Show how long it takes to process a model.
800    """
801    from .models import cylinder
802    import datetime
803    tic = datetime.datetime.now()
804    make_source(make_model_info(cylinder))
805    toc = (datetime.datetime.now() - tic).total_seconds()
806    print("time: %g"%toc)
807
808def main():
809    """
810    Program which prints the source produced by the model.
811    """
812    import sys
813    from sasmodels.core import make_model_by_name
814    if len(sys.argv) <= 1:
815        print("usage: python -m sasmodels.generate modelname")
816    else:
817        name = sys.argv[1]
818        model_info = make_model_by_name(name)
819        source = make_source(model_info)
820        print(source)
821
822if __name__ == "__main__":
823    main()
Note: See TracBrowser for help on using the repository browser.