source: sasmodels/sasmodels/generate.py @ 7bb290c

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

Convert subsection headings to emphasized text in model docs

  • Property mode set to 100644
File size: 25.7 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
24These functions are defined in a kernel module .py script and an associated
25set of .c files.  The model constructor will use them to create models with
26polydispersity across volume and orientation parameters, and provide
27scale and background parameters for each model.
28
29*Iq*, *Iqxy*, *Imagnetic* and *form_volume* should be stylized C-99
30functions written for OpenCL.  All functions need prototype declarations
31even if the are defined before they are used.  OpenCL does not support
32*#include* preprocessor directives, so instead the list of includes needs
33to be given as part of the metadata in the kernel module definition.
34The included files should be listed using a path relative to the kernel
35module, or if using "lib/file.c" if it is one of the standard includes
36provided with the sasmodels source.  The includes need to be listed in
37order so that functions are defined before they are used.
38
39Floating point values should be declared as *double*.  For single precision
40calculations, *double* will be replaced by *float*.  The single precision
41conversion will also tag floating point constants with "f" to make them
42single precision constants.  When using integral values in floating point
43expressions, they should be expressed as floating point values by including
44a decimal point.  This includes 0., 1. and 2.
45
46OpenCL has a *sincos* function which can improve performance when both
47the *sin* and *cos* values are needed for a particular argument.  Since
48this function does not exist in C99, all use of *sincos* should be
49replaced by the macro *SINCOS(value,sn,cn)* where *sn* and *cn* are
50previously declared *double* variables.  When compiled for systems without
51OpenCL, *SINCOS* will be replaced by *sin* and *cos* calls.   If *value* is
52an expression, it will appear twice in this case; whether or not it will be
53evaluated twice depends on the quality of the compiler.
54
55If the input parameters are invalid, the scattering calculator should
56return a negative number. Particularly with polydispersity, there are
57some sets of shape parameters which lead to nonsensical forms, such
58as a capped cylinder where the cap radius is smaller than the
59cylinder radius.  The polydispersity calculation will ignore these points,
60effectively chopping the parameter weight distributions at the boundary
61of the infeasible region.  The resulting scattering will be set to
62background.  This will work correctly even when polydispersity is off.
63
64*ER* and *VR* are python functions which operate on parameter vectors.
65The constructor code will generate the necessary vectors for computing
66them with the desired polydispersity.
67
68The available kernel parameters are defined as a list, with each parameter
69defined as a sublist with the following elements:
70
71    *name* is the name that will be used in the call to the kernel
72    function and the name that will be displayed to the user.  Names
73    should be lower case, with words separated by underscore.  If
74    acronyms are used, the whole acronym should be upper case.
75
76    *units* should be one of *degrees* for angles, *Ang* for lengths,
77    *1e-6/Ang^2* for SLDs.
78
79    *default value* will be the initial value for  the model when it
80    is selected, or when an initial value is not otherwise specified.
81
82    [*lb*, *ub*] are the hard limits on the parameter value, used to limit
83    the polydispersity density function.  In the fit, the parameter limits
84    given to the fit are the limits  on the central value of the parameter.
85    If there is polydispersity, it will evaluate parameter values outside
86    the fit limits, but not outside the hard limits specified in the model.
87    If there are no limits, use +/-inf imported from numpy.
88
89    *type* indicates how the parameter will be used.  "volume" parameters
90    will be used in all functions.  "orientation" parameters will be used
91    in *Iqxy* and *Imagnetic*.  "magnetic* parameters will be used in
92    *Imagnetic* only.  If *type* is the empty string, the parameter will
93    be used in all of *Iq*, *Iqxy* and *Imagnetic*.
94
95    *description* is a short description of the parameter.  This will
96    be displayed in the parameter table and used as a tool tip for the
97    parameter value in the user interface.
98
99The kernel module must set variables defining the kernel meta data:
100
101    *id* is an implicit variable formed from the filename.  It will be
102    a valid python identifier, and will be used as the reference into
103    the html documentation, with '_' replaced by '-'.
104
105    *name* is the model name as displayed to the user.  If it is missing,
106    it will be constructed from the id.
107
108    *title* is a short description of the model, suitable for a tool tip,
109    or a one line model summary in a table of models.
110
111    *description* is an extended description of the model to be displayed
112    while the model parameters are being edited.
113
114    *parameters* is the list of parameters.  Parameters in the kernel
115    functions must appear in the same order as they appear in the
116    parameters list.  Two additional parameters, *scale* and *background*
117    are added to the beginning of the parameter list.  They will show up
118    in the documentation as model parameters, but they are never sent to
119    the kernel functions.
120
121    *category* is the default category for the model.  Models in the
122    *structure-factor* category do not have *scale* and *background*
123    added.
124
125    *source* is the list of C-99 source files that must be joined to
126    create the OpenCL kernel functions.  The files defining the functions
127    need to be listed before the files which use the functions.
128
129    *ER* is a python function defining the effective radius.  If it is
130    not present, the effective radius is 0.
131
132    *VR* is a python function defining the volume ratio.  If it is not
133    present, the volume ratio is 1.
134
135    *form_volume*, *Iq*, *Iqxy*, *Imagnetic* are strings containing the
136    C source code for the body of the volume, Iq, and Iqxy functions
137    respectively.  These can also be defined in the last source file.
138
139    *Iq* and *Iqxy* also be instead be python functions defining the
140    kernel.  If they are marked as *Iq.vectorized = True* then the
141    kernel is passed the entire *q* vector at once, otherwise it is
142    passed values one *q* at a time.  The performance improvement of
143    this step is significant.
144
145    *demo* is a dictionary of parameter=value defining a set of
146    parameters to use by default when *compare* is called.  Any
147    parameter not set in *demo* gets the initial value from the
148    parameter list.  *demo* is mostly needed to set the default
149    polydispersity values for tests.
150
151    *oldname* is the name of the model in sasview before sasmodels
152    was split into its own package, and *oldpars* is a dictionary
153    of *parameter: old_parameter* pairs defining the new names for
154    the parameters.  This is used by *compare* to check the values
155    of the new model against the values of the old model before
156    you are ready to add the new model to sasmodels.
157
158
159An *info* dictionary is constructed from the kernel meta data and
160returned to the caller.
161
162The model evaluator, function call sequence consists of q inputs and the return vector,
163followed by the loop value/weight vector, followed by the values for
164the non-polydisperse parameters, followed by the lengths of the
165polydispersity loops.  To construct the call for 1D models, the
166categories *fixed-1d* and *pd-1d* list the names of the parameters
167of the non-polydisperse and the polydisperse parameters respectively.
168Similarly, *fixed-2d* and *pd-2d* provide parameter names for 2D models.
169The *pd-rel* category is a set of those parameters which give
170polydispersitiy as a portion of the value (so a 10% length dispersity
171would use a polydispersity value of 0.1) rather than absolute
172dispersity such as an angle plus or minus 15 degrees.
173
174The *volume* category lists the volume parameters in order for calls
175to volume within the kernel (used for volume normalization) and for
176calls to ER and VR for effective radius and volume ratio respectively.
177
178The *orientation* and *magnetic* categories list the orientation and
179magnetic parameters.  These are used by the sasview interface.  The
180blank category is for parameters such as scale which don't have any
181other marking.
182
183The doc string at the start of the kernel module will be used to
184construct the model documentation web pages.  Embedded figures should
185appear in the subdirectory "img" beside the model definition, and tagged
186with the kernel module name to avoid collision with other models.  Some
187file systems are case-sensitive, so only use lower case characters for
188file names and extensions.
189
190
191The function :func:`make` loads the metadata from the module and returns
192the kernel source.  The function :func:`doc` extracts the doc string
193and adds the parameter table to the top.  The function :func:`sources`
194returns a list of files required by the model.
195"""
196
197# TODO: identify model files which have changed since loading and reload them.
198
199__all__ = ["make", "doc", "sources", "use_single", "use_long_double"]
200
201import sys
202from os.path import abspath, dirname, join as joinpath, exists, basename, \
203    splitext
204import re
205import string
206
207import numpy as np
208C_KERNEL_TEMPLATE_PATH = joinpath(dirname(__file__), 'kernel_template.c')
209
210F32 = np.dtype('float32')
211F64 = np.dtype('float64')
212try:  # CRUFT: older numpy does not support float128
213    F128 = np.dtype('float128')
214except TypeError:
215    F128 = None
216
217
218# Scale and background, which are parameters common to every form factor
219COMMON_PARAMETERS = [
220    ["scale", "", 1, [0, np.inf], "", "Source intensity"],
221    ["background", "1/cm", 0, [0, np.inf], "", "Source background"],
222    ]
223
224
225# Conversion from units defined in the parameter table for each model
226# to units displayed in the sphinx documentation.
227RST_UNITS = {
228    "Ang": "|Ang|",
229    "1/Ang": "|Ang^-1|",
230    "1/Ang^2": "|Ang^-2|",
231    "1e-6/Ang^2": "|1e-6Ang^-2|",
232    "degrees": "degree",
233    "1/cm": "|cm^-1|",
234    "": "None",
235    }
236
237# Headers for the parameters tables in th sphinx documentation
238PARTABLE_HEADERS = [
239    "Parameter",
240    "Description",
241    "Units",
242    "Default value",
243    ]
244
245# Minimum width for a default value (this is shorter than the column header
246# width, so will be ignored).
247PARTABLE_VALUE_WIDTH = 10
248
249# Documentation header for the module, giving the model name, its short
250# description and its parameter table.  The remainder of the doc comes
251# from the module docstring.
252DOC_HEADER = """.. _%(id)s:
253
254%(name)s
255=======================================================
256
257%(title)s
258
259%(parameters)s
260
261%(returns)s
262
263%(docs)s
264"""
265
266def format_units(par):
267    return RST_UNITS.get(par, par)
268
269def make_partable(pars):
270    """
271    Generate the parameter table to include in the sphinx documentation.
272    """
273    column_widths = [
274        max(len(p[0]) for p in pars),
275        max(len(p[-1]) for p in pars),
276        max(len(format_units(p[1])) for p in pars),
277        PARTABLE_VALUE_WIDTH,
278        ]
279    column_widths = [max(w, len(h))
280                     for w, h in zip(column_widths, PARTABLE_HEADERS)]
281
282    sep = " ".join("="*w for w in column_widths)
283    lines = [
284        sep,
285        " ".join("%-*s" % (w, h) for w, h in zip(column_widths, PARTABLE_HEADERS)),
286        sep,
287        ]
288    for p in pars:
289        lines.append(" ".join([
290            "%-*s" % (column_widths[0], p[0]),
291            "%-*s" % (column_widths[1], p[-1]),
292            "%-*s" % (column_widths[2], format_units(p[1])),
293            "%*g" % (column_widths[3], p[2]),
294            ]))
295    lines.append(sep)
296    return "\n".join(lines)
297
298def _search(search_path, filename):
299    """
300    Find *filename* in *search_path*.
301
302    Raises ValueError if file does not exist.
303    """
304    for path in search_path:
305        target = joinpath(path, filename)
306        if exists(target):
307            return target
308    raise ValueError("%r not found in %s" % (filename, search_path))
309
310def sources(info):
311    """
312    Return a list of the sources file paths for the module.
313    """
314    search_path = [dirname(info['filename']),
315                   abspath(joinpath(dirname(__file__), 'models'))]
316    return [_search(search_path, f) for f in info['source']]
317
318def use_single(source):
319    """
320    Convert code from double precision to single precision.
321    """
322    # Convert double keyword to float.  Accept an 'n' parameter for vector
323    # values, where n is 2, 4, 8 or 16. Assume complex numbers are represented
324    # as cdouble which is typedef'd to double2.
325    source = re.sub(r'(^|[^a-zA-Z0-9_]c?)double(([248]|16)?($|[^a-zA-Z0-9_]))',
326                    r'\1float\2', source)
327    # Convert floating point constants to single by adding 'f' to the end.
328    # OS/X driver complains if you don't do this.
329    source = re.sub(r'[^a-zA-Z_](\d*[.]\d+|\d+[.]\d*)([eE][+-]?\d+)?',
330                    r'\g<0>f', source)
331    return source
332
333def use_long_double(source):
334    """
335    Convert code from double precision to long double precision.
336    """
337    # Convert double keyword to float.  Accept an 'n' parameter for vector
338    # values, where n is 2, 4, 8 or 16. Assume complex numbers are represented
339    # as cdouble which is typedef'd to double2.
340    source = re.sub(r'(^|[^a-zA-Z0-9_]c?)double(([248]|16)?($|[^a-zA-Z0-9_]))',
341                    r'\1long double\2', source)
342    # Convert floating point constants to single by adding 'f' to the end.
343    # OS/X driver complains if you don't do this.
344    source = re.sub(r'[^a-zA-Z_](\d*[.]\d+|\d+[.]\d*)([eE][+-]?\d+)?',
345                    r'\g<0>L', source)
346    return source
347
348
349def kernel_name(info, is_2D):
350    """
351    Name of the exported kernel symbol.
352    """
353    return info['name'] + "_" + ("Iqxy" if is_2D else "Iq")
354
355
356def categorize_parameters(pars):
357    """
358    Build parameter categories out of the the parameter definitions.
359
360    Returns a dictionary of categories.
361    """
362    partype = {
363        'volume': [], 'orientation': [], 'magnetic': [], '': [],
364        'fixed-1d': [], 'fixed-2d': [], 'pd-1d': [], 'pd-2d': [],
365        'pd-rel': set(),
366    }
367
368    for p in pars:
369        name, ptype = p[0], p[4]
370        if ptype == 'volume':
371            partype['pd-1d'].append(name)
372            partype['pd-2d'].append(name)
373            partype['pd-rel'].add(name)
374        elif ptype == 'magnetic':
375            partype['fixed-2d'].append(name)
376        elif ptype == 'orientation':
377            partype['pd-2d'].append(name)
378        elif ptype == '':
379            partype['fixed-1d'].append(name)
380            partype['fixed-2d'].append(name)
381        else:
382            raise ValueError("unknown parameter type %r" % ptype)
383        partype[ptype].append(name)
384
385    return partype
386
387def indent(s, depth):
388    """
389    Indent a string of text with *depth* additional spaces on each line.
390    """
391    spaces = " "*depth
392    sep = "\n" + spaces
393    return spaces + sep.join(s.split("\n"))
394
395
396def build_polydispersity_loops(pd_pars):
397    """
398    Build polydispersity loops
399
400    Returns loop opening and loop closing
401    """
402    LOOP_OPEN = """\
403for (int %(name)s_i=0; %(name)s_i < N%(name)s; %(name)s_i++) {
404  const double %(name)s = loops[2*(%(name)s_i%(offset)s)];
405  const double %(name)s_w = loops[2*(%(name)s_i%(offset)s)+1];\
406"""
407    depth = 4
408    offset = ""
409    loop_head = []
410    loop_end = []
411    for name in pd_pars:
412        subst = {'name': name, 'offset': offset}
413        loop_head.append(indent(LOOP_OPEN % subst, depth))
414        loop_end.insert(0, (" "*depth) + "}")
415        offset += '+N' + name
416        depth += 2
417    return "\n".join(loop_head), "\n".join(loop_end)
418
419C_KERNEL_TEMPLATE = None
420def make_model(info):
421    """
422    Generate the code for the kernel defined by info, using source files
423    found in the given search path.
424    """
425    # TODO: need something other than volume to indicate dispersion parameters
426    # No volume normalization despite having a volume parameter.
427    # Thickness is labelled a volume in order to trigger polydispersity.
428    # May want a separate dispersion flag, or perhaps a separate category for
429    # disperse, but not volume.  Volume parameters also use relative values
430    # for the distribution rather than the absolute values used by angular
431    # dispersion.  Need to be careful that necessary parameters are available
432    # for computing volume even if we allow non-disperse volume parameters.
433
434    # Load template
435    global C_KERNEL_TEMPLATE
436    if C_KERNEL_TEMPLATE is None:
437        with open(C_KERNEL_TEMPLATE_PATH) as fid:
438            C_KERNEL_TEMPLATE = fid.read()
439
440    # Load additional sources
441    source = [open(f).read() for f in sources(info)]
442
443    # Prepare defines
444    defines = []
445    partype = info['partype']
446    pd_1d = partype['pd-1d']
447    pd_2d = partype['pd-2d']
448    fixed_1d = partype['fixed-1d']
449    fixed_2d = partype['fixed-1d']
450
451    iq_parameters = [p[0]
452                     for p in info['parameters'][2:] # skip scale, background
453                     if p[0] in set(fixed_1d + pd_1d)]
454    iqxy_parameters = [p[0]
455                       for p in info['parameters'][2:] # skip scale, background
456                       if p[0] in set(fixed_2d + pd_2d)]
457    volume_parameters = [p[0]
458                         for p in info['parameters']
459                         if p[4] == 'volume']
460
461    # Fill in defintions for volume parameters
462    if volume_parameters:
463        defines.append(('VOLUME_PARAMETERS',
464                        ','.join(volume_parameters)))
465        defines.append(('VOLUME_WEIGHT_PRODUCT',
466                        '*'.join(p + '_w' for p in volume_parameters)))
467
468    # Generate form_volume function from body only
469    if info['form_volume'] is not None:
470        if volume_parameters:
471            vol_par_decl = ', '.join('double ' + p for p in volume_parameters)
472        else:
473            vol_par_decl = 'void'
474        defines.append(('VOLUME_PARAMETER_DECLARATIONS',
475                        vol_par_decl))
476        fn = """\
477double form_volume(VOLUME_PARAMETER_DECLARATIONS);
478double form_volume(VOLUME_PARAMETER_DECLARATIONS) {
479    %(body)s
480}
481""" % {'body':info['form_volume']}
482        source.append(fn)
483
484    # Fill in definitions for Iq parameters
485    defines.append(('IQ_KERNEL_NAME', info['name'] + '_Iq'))
486    defines.append(('IQ_PARAMETERS', ', '.join(iq_parameters)))
487    if fixed_1d:
488        defines.append(('IQ_FIXED_PARAMETER_DECLARATIONS',
489                        ', \\\n    '.join('const double %s' % p for p in fixed_1d)))
490    if pd_1d:
491        defines.append(('IQ_WEIGHT_PRODUCT',
492                        '*'.join(p + '_w' for p in pd_1d)))
493        defines.append(('IQ_DISPERSION_LENGTH_DECLARATIONS',
494                        ', \\\n    '.join('const int N%s' % p for p in pd_1d)))
495        defines.append(('IQ_DISPERSION_LENGTH_SUM',
496                        '+'.join('N' + p for p in pd_1d)))
497        open_loops, close_loops = build_polydispersity_loops(pd_1d)
498        defines.append(('IQ_OPEN_LOOPS',
499                        open_loops.replace('\n', ' \\\n')))
500        defines.append(('IQ_CLOSE_LOOPS',
501                        close_loops.replace('\n', ' \\\n')))
502    if info['Iq'] is not None:
503        defines.append(('IQ_PARAMETER_DECLARATIONS',
504                        ', '.join('double ' + p for p in iq_parameters)))
505        fn = """\
506double Iq(double q, IQ_PARAMETER_DECLARATIONS);
507double Iq(double q, IQ_PARAMETER_DECLARATIONS) {
508    %(body)s
509}
510""" % {'body':info['Iq']}
511        source.append(fn)
512
513    # Fill in definitions for Iqxy parameters
514    defines.append(('IQXY_KERNEL_NAME', info['name'] + '_Iqxy'))
515    defines.append(('IQXY_PARAMETERS', ', '.join(iqxy_parameters)))
516    if fixed_2d:
517        defines.append(('IQXY_FIXED_PARAMETER_DECLARATIONS',
518                        ', \\\n    '.join('const double %s' % p for p in fixed_2d)))
519    if pd_2d:
520        defines.append(('IQXY_WEIGHT_PRODUCT',
521                        '*'.join(p + '_w' for p in pd_2d)))
522        defines.append(('IQXY_DISPERSION_LENGTH_DECLARATIONS',
523                        ', \\\n    '.join('const int N%s' % p for p in pd_2d)))
524        defines.append(('IQXY_DISPERSION_LENGTH_SUM',
525                        '+'.join('N' + p for p in pd_2d)))
526        open_loops, close_loops = build_polydispersity_loops(pd_2d)
527        defines.append(('IQXY_OPEN_LOOPS',
528                        open_loops.replace('\n', ' \\\n')))
529        defines.append(('IQXY_CLOSE_LOOPS',
530                        close_loops.replace('\n', ' \\\n')))
531    if info['Iqxy'] is not None:
532        defines.append(('IQXY_PARAMETER_DECLARATIONS',
533                        ', '.join('double ' + p for p in iqxy_parameters)))
534        fn = """\
535double Iqxy(double qx, double qy, IQXY_PARAMETER_DECLARATIONS);
536double Iqxy(double qx, double qy, IQXY_PARAMETER_DECLARATIONS) {
537    %(body)s
538}
539""" % {'body':info['Iqxy']}
540        source.append(fn)
541
542    # Need to know if we have a theta parameter for Iqxy; it is not there
543    # for the magnetic sphere model, for example, which has a magnetic
544    # orientation but no shape orientation.
545    if 'theta' in pd_2d:
546        defines.append(('IQXY_HAS_THETA', '1'))
547
548    #for d in defines: print d
549    DEFINES = '\n'.join('#define %s %s' % (k, v) for k, v in defines)
550    SOURCES = '\n\n'.join(source)
551    return C_KERNEL_TEMPLATE % {
552        'DEFINES':DEFINES,
553        'SOURCES':SOURCES,
554        }
555
556def make_info(kernel_module):
557    """
558    Interpret the model definition file, categorizing the parameters.
559    """
560    #print kernelfile
561    category = getattr(kernel_module, 'category', None)
562    parameters = COMMON_PARAMETERS + kernel_module.parameters
563    # Default the demo parameters to the starting values for the individual
564    # parameters if an explicit demo parameter set has not been specified.
565    demo_parameters = getattr(kernel_module, 'demo', None)
566    if demo_parameters is None:
567        demo_parameters = dict((p[0],p[2]) for p in parameters)
568    filename = abspath(kernel_module.__file__)
569    kernel_id = splitext(basename(filename))[0]
570    name = getattr(kernel_module, 'name', None)
571    if name is None:
572        name = " ".join(w.capitalize() for w in kernel_id.split('_'))
573    info = dict(
574        id = kernel_id,  # string used to load the kernel
575        filename=abspath(kernel_module.__file__),
576        name=name,
577        title=kernel_module.title,
578        description=kernel_module.description,
579        category=category,
580        parameters=parameters,
581        demo=demo_parameters,
582        source=getattr(kernel_module, 'source', []),
583        oldname=kernel_module.oldname,
584        oldpars=kernel_module.oldpars,
585        )
586    # Fill in attributes which default to None
587    info.update((k, getattr(kernel_module, k, None))
588                for k in ('ER', 'VR', 'form_volume', 'Iq', 'Iqxy'))
589    # Fill in the derived attributes
590    info['limits'] = dict((p[0], p[3]) for p in info['parameters'])
591    info['partype'] = categorize_parameters(info['parameters'])
592    info['defaults'] = dict((p[0], p[2]) for p in info['parameters'])
593    return info
594
595def make(kernel_module):
596    """
597    Build an OpenCL/ctypes function from the definition in *kernel_module*.
598
599    The module can be loaded with a normal python import statement if you
600    know which module you need, or with __import__('sasmodels.model.'+name)
601    if the name is in a string.
602    """
603    info = make_info(kernel_module)
604    # Assume if one part of the kernel is python then all parts are.
605    source = make_model(info) if not callable(info['Iq']) else None
606    return source, info
607
608section_marker = re.compile(r'\A(?P<first>[%s])(?P=first)*\Z'
609                            %re.escape(string.punctuation))
610def _convert_section_titles_to_boldface(lines):
611    prior = None
612    for line in lines:
613        if prior is None:
614            prior = line
615        elif section_marker.match(line):
616            if len(line) >= len(prior):
617                yield "".join( ("**",prior,"**") )
618                prior = None
619            else:
620                yield prior
621                prior = line
622        else:
623            yield prior
624            prior = line
625    if prior is not None:
626        yield prior
627
628def convert_section_titles_to_boldface(string):
629    return "\n".join(_convert_section_titles_to_boldface(string.split('\n')))
630
631def doc(kernel_module):
632    """
633    Return the documentation for the model.
634    """
635    Iq_units = "The returned value is scaled to units of |cm^-1| |sr^-1|, absolute scale."
636    Sq_units = "The returned value is a dimensionless structure factor, $S(q)$."
637    info = make_info(kernel_module)
638    is_Sq = ("structure-factor" in info['category'])
639    #docs = kernel_module.__doc__
640    docs = convert_section_titles_to_boldface(kernel_module.__doc__)
641    subst = dict(id=info['id'].replace('_', '-'),
642                 name=info['name'],
643                 title=info['title'],
644                 parameters=make_partable(info['parameters']),
645                 returns=Sq_units if is_Sq else Iq_units,
646                 docs=docs)
647    return DOC_HEADER % subst
648
649
650
651def demo_time():
652    from .models import cylinder
653    import datetime
654    tic = datetime.datetime.now()
655    make(cylinder)
656    toc = (datetime.datetime.now() - tic).total_seconds()
657    print "time:", toc
658
659def main():
660    if len(sys.argv) <= 1:
661        print "usage: python -m sasmodels.generate modelname"
662    else:
663        name = sys.argv[1]
664        import sasmodels.models
665        __import__('sasmodels.models.' + name)
666        model = getattr(sasmodels.models, name)
667        source, _ = make(model)
668        print source
669
670if __name__ == "__main__":
671    main()
Note: See TracBrowser for help on using the repository browser.