source: sasmodels/sasmodels/convert.py @ b3f4831

ticket_1156
Last change on this file since b3f4831 was b3f4831, checked in by Paul Kienzle <pkienzle@…>, 12 months ago

Merge branch 'master' into ticket_1156

  • Property mode set to 100644
File size: 24.9 KB
Line 
1"""
2Convert models to and from sasview.
3"""
4from __future__ import print_function, division
5
6import math
7import warnings
8
9import numpy as np
10
11from .conversion_table import CONVERSION_TABLE
12from .core import load_model_info
13
14# List of models which SasView versions don't contain the explicit 'scale' argument.
15# When converting such a model, please update this list.
16MODELS_WITHOUT_SCALE = [
17    'teubner_strey',
18    'broad_peak',
19    'two_lorentzian',
20    "two_power_law",
21    'gauss_lorentz_gel',
22    'be_polyelectrolyte',
23    'correlation_length',
24    'fractal_core_shell',
25    'binary_hard_sphere',
26    'raspberry'
27]
28
29# List of models which SasView versions don't contain the explicit 'background' argument.
30# When converting such a model, please update this list.
31MODELS_WITHOUT_BACKGROUND = [
32    'guinier',
33]
34
35MODELS_WITHOUT_VOLFRACTION = [
36    'fractal',
37    'vesicle',
38    'multilayer_vesicle',
39]
40
41MAGNETIC_SASVIEW_MODELS = [
42    'core_shell',
43    'core_multi_shell',
44    'cylinder',
45    'parallelepiped',
46    'sphere',
47]
48
49
50# Convert new style names for polydispersity info to old style names
51PD_DOT = [
52    ("_pd", ".width"),
53    ("_pd_n", ".npts"),
54    ("_pd_nsigma", ".nsigmas"),
55    ("_pd_type", ".type"),
56    (".lower", ".lower"),
57    (".upper", ".upper"),
58    (".fittable", ".fittable"),
59    (".std", ".std"),
60    (".units", ".units"),
61    ("", "")
62    ]
63
64def _rescale(par, scale):
65    return [pk*scale for pk in par] if isinstance(par, list) else par*scale
66
67def _is_sld(model_info, par):
68    """
69    Return True if parameter is a magnetic magnitude or SLD parameter.
70    """
71    if par.startswith('M0:'):
72        return True
73    if '_pd' in par or '.' in par:
74        return False
75    for p in model_info.parameters.call_parameters:
76        if p.id == par:
77            return p.type == 'sld'
78    # check through kernel parameters in case it is a named as a vector
79    for p in model_info.parameters.kernel_parameters:
80        if p.id == par:
81            return p.type == 'sld'
82    return False
83
84def _rescale_sld(model_info, pars, scale):
85    """
86    rescale all sld parameters in the new model definition by *scale* so the
87    numbers are nicer.  Relies on the fact that all sld parameters in the
88    new model definition end with sld.  For backward conversion use
89    *scale=1e-6*.  For forward conversion use *scale=1e6*.
90    """
91    return dict((par, (_rescale(v, scale) if _is_sld(model_info, par) else v))
92                for par, v in pars.items())
93
94
95def _get_translation_table(model_info, version=(3, 1, 2)):
96    conv_param = CONVERSION_TABLE.get(version, {}).get(model_info.id, [None, {}])
97    translation = conv_param[1].copy()
98    for p in model_info.parameters.kernel_parameters:
99        if p.length > 1:
100            newid = p.id
101            oldid = translation.get(p.id, p.id)
102            translation.pop(newid, None)
103            for k in range(1, p.length+1):
104                if newid+str(k) not in translation:
105                    translation[newid+str(k)] = oldid+str(k)
106    # Remove control parameter from the result
107    if model_info.control:
108        translation[model_info.control] = "CONTROL"
109    return translation
110
111# ========= FORWARD CONVERSION sasview 3.x => sasmodels ===========
112def _dot_pd_to_underscore_pd(par):
113    if par.endswith(".width"):
114        return par[:-6]+"_pd"
115    elif par.endswith(".type"):
116        return par[:-5]+"_pd_type"
117    elif par.endswith(".nsigmas"):
118        return par[:-8]+"_pd_nsigma"
119    elif par.endswith(".npts"):
120        return par[:-5]+"_pd_n"
121    else:
122        return par
123
124def _pd_to_underscores(pars):
125    return dict((_dot_pd_to_underscore_pd(k), v) for k, v in pars.items())
126
127def _convert_pars(pars, mapping):
128    """
129    Rename the parameters and any associated polydispersity attributes.
130    """
131    newpars = pars.copy()
132    for new, old in mapping.items():
133        if old == new:
134            continue
135        if old is None:
136            continue
137        for _, dot in PD_DOT:
138            source = old+dot
139            if source in newpars:
140                if new is not None:
141                    target = new+dot
142                else:
143                    target = None
144                if source != target:
145                    if target:
146                        newpars[target] = pars[old+dot]
147                    del newpars[source]
148    return newpars
149
150def _conversion_target(model_name, version=(3, 1, 2)):
151    """
152    Find the sasmodel name which translates into the sasview name.
153
154    Note: *CoreShellEllipsoidModel* translates into *core_shell_ellipsoid:1*.
155    This is necessary since there is only one variant in sasmodels for the
156    two variants in sasview.
157    """
158    for sasmodels_name, sasview_dict in \
159            CONVERSION_TABLE.get(version, {}).items():
160        if sasview_dict[0] == model_name:
161            return sasmodels_name
162    return None
163
164def _hand_convert(name, oldpars, version=(3, 1, 2)):
165    if version == (3, 1, 2):
166        oldpars = _hand_convert_3_1_2_to_4_1(name, oldpars)
167    if version < (4, 2, 0):
168        oldpars = _hand_convert_4_1_to_4_2(name, oldpars)
169        oldpars = _rename_magnetic_pars(oldpars)
170    return oldpars
171
172def _hand_convert_4_1_to_4_2(name, oldpars):
173    if name in ('bcc_paracrystal', 'fcc_paracrystal', 'sc_paracrystal'):
174        oldpars['lattice_spacing'] = oldpars.pop('dnn')
175        oldpars['lattice_distortion'] = oldpars.pop('d_factor')
176    return oldpars
177
178def _rename_magnetic_pars(pars):
179    """
180    Change from M0:par to par_M0, etc.
181    """
182    keys = list(pars.items())
183    for k in keys:
184        if k.startswith('M0:'):
185            pars[k[3:]+'_M0'] = pars.pop(k)
186        elif k.startswith('mtheta:'):
187            pars[k[7:]+'_mtheta'] = pars.pop(k)
188        elif k.startswith('mphi:'):
189            pars[k[5:]+'_mphi'] = pars.pop(k)
190        elif k.startswith('up:'):
191            pars['up_'+k[3:]] = pars.pop(k)
192    return pars
193
194def _hand_convert_3_1_2_to_4_1(name, oldpars):
195    if name == 'core_shell_parallelepiped':
196        # Make sure pd on rim parameters defaults to zero
197        # ... probably not necessary.
198        oldpars['rimA.width'] = 0.0
199        oldpars['rimB.width'] = 0.0
200        oldpars['rimC.width'] = 0.0
201    elif name == 'core_shell_ellipsoid:1':
202        # Reverse translation (from new to old), from core_shell_ellipsoid.c
203        #    equat_shell = equat_core + thick_shell
204        #    polar_core = equat_core * x_core
205        #    polar_shell = equat_core * x_core + thick_shell*x_polar_shell
206        # Forward translation (from old to new), inverting reverse translation:
207        #    thick_shell = equat_shell - equat_core
208        #    x_core = polar_core / equat_core
209        #    x_polar_shell = (polar_shell - polar_core)/(equat_shell - equat_core)
210        # Auto translation (old <=> new) happens after hand_convert
211        #    equat_shell <=> thick_shell
212        #    polar_core <=> x_core
213        #    polar_shell <=> x_polar_shell
214        # So...
215        equat_core, equat_shell = oldpars['equat_core'], oldpars['equat_shell']
216        polar_core, polar_shell = oldpars['polar_core'], oldpars['polar_shell']
217        oldpars['equat_shell'] = equat_shell - equat_core
218        oldpars['polar_core'] = polar_core / equat_core
219        oldpars['polar_shell'] = (polar_shell-polar_core)/(equat_shell-equat_core)
220    elif name == 'hollow_cylinder':
221        # now uses radius and thickness
222        thickness = oldpars['radius'] - oldpars['core_radius']
223        oldpars['radius'] = thickness
224        if 'radius.width' in oldpars:
225            pd = oldpars['radius.width']*oldpars['radius']/thickness
226            oldpars['radius.width'] = pd
227    elif name == 'multilayer_vesicle':
228        if 'scale' in oldpars:
229            oldpars['volfraction'] = oldpars['scale']
230            oldpars['scale'] = 1.0
231        if 'scale.lower' in oldpars:
232            oldpars['volfraction.lower'] = oldpars['scale.lower']
233        if 'scale.upper' in oldpars:
234            oldpars['volfraction.upper'] = oldpars['scale.upper']
235        if 'scale.fittable' in oldpars:
236            oldpars['volfraction.fittable'] = oldpars['scale.fittable']
237        if 'scale.std' in oldpars:
238            oldpars['volfraction.std'] = oldpars['scale.std']
239        if 'scale.units' in oldpars:
240            oldpars['volfraction.units'] = oldpars['scale.units']
241    elif name == 'pearl_necklace':
242        pass
243        #_remove_pd(oldpars, 'num_pearls', name)
244        #_remove_pd(oldpars, 'thick_string', name)
245    elif name == 'polymer_micelle':
246        if 'ndensity' in oldpars:
247            oldpars['ndensity'] /= 1e15
248        if 'ndensity.lower' in oldpars:
249            oldpars['ndensity.lower'] /= 1e15
250        if 'ndensity.upper' in oldpars:
251            oldpars['ndensity.upper'] /= 1e15
252    elif name == 'rpa':
253        # convert scattering lengths from femtometers to centimeters
254        for p in "L1", "L2", "L3", "L4":
255            if p in oldpars:
256                oldpars[p] /= 1e-13
257            if p + ".lower" in oldpars:
258                oldpars[p + ".lower"] /= 1e-13
259            if p + ".upper" in oldpars:
260                oldpars[p + ".upper"] /= 1e-13
261    elif name == 'spherical_sld':
262        j = 0
263        while "func_inter" + str(j) in oldpars:
264            name = "func_inter" + str(j)
265            new_name = "shape" + str(j + 1)
266            if oldpars[name] == 'Erf(|nu|*z)':
267                oldpars[new_name] = int(0)
268            elif oldpars[name] == 'RPower(z^|nu|)':
269                oldpars[new_name] = int(1)
270            elif oldpars[name] == 'LPower(z^|nu|)':
271                oldpars[new_name] = int(2)
272            elif oldpars[name] == 'RExp(-|nu|*z)':
273                oldpars[new_name] = int(3)
274            elif oldpars[name] == 'LExp(-|nu|*z)':
275                oldpars[new_name] = int(4)
276            else:
277                oldpars[new_name] = int(0)
278            oldpars.pop(name)
279            oldpars['n_shells'] = str(j + 1)
280            j += 1
281    elif name == 'teubner_strey':
282        # basically undoing the entire Teubner-Strey calculations here.
283        #    drho = (sld_a - sld_b)
284        #    k = 2.0*math.pi*xi/d
285        #    a2 = (1.0 + k**2)**2
286        #    c1 = 2.0 * xi**2 * (1.0 - k**2)
287        #    c2 = xi**4
288        #    prefactor = 8.0*math.pi*phi*(1.0-phi)*drho**2*c2/xi
289        #    scale = 1e-4*prefactor
290        #    oldpars['scale'] = a2/scale
291        #    oldpars['c1'] = c1/scale
292        #    oldpars['c2'] = c2/scale
293
294        # need xi, d, sld_a, sld_b, phi=volfraction_a
295        # assume contrast is 1.0e-6, scale=1, background=0
296        sld_a, sld_b = 1.0, 0.
297        drho = sld_a - sld_b
298
299        # find xi
300        p_scale = oldpars['scale']
301        p_c1 = oldpars['c1']
302        p_c2 = oldpars['c2']
303        i_1 = 0.5*p_c1/p_c2
304        i_2 = math.sqrt(math.fabs(p_scale/p_c2))
305        i_3 = 2/(i_1 + i_2)
306        xi = math.sqrt(math.fabs(i_3))
307
308        # find d from xi
309        k = math.sqrt(math.fabs(1 - 0.5*p_c1/p_c2*xi**2))
310        d = 2*math.pi*xi/k
311
312        # solve quadratic phi (1-phi) = xi/(1e-4 8 pi drho^2 c2)
313        # favour volume fraction in [0, 0.5]
314        c = xi / (1e-4 * 8.0 * math.pi * drho**2 * p_c2)
315        phi = 0.5 - math.sqrt(0.25 - c)
316
317        # scale sld_a by 1e-6 because the translator will scale it back
318        oldpars.update(volfraction_a=phi, xi=xi, d=d, sld_a=sld_a*1e-6,
319                       sld_b=sld_b, scale=1.0)
320        oldpars.pop('c1')
321        oldpars.pop('c2')
322
323    return oldpars
324
325def convert_model(name, pars, use_underscore=False, model_version=(3, 1, 2)):
326    """
327    Convert model from old style parameter names to new style.
328    """
329    newpars = pars
330    keys = sorted(CONVERSION_TABLE.keys())
331    for i, version in enumerate(keys):
332        # Don't allow indices outside list
333        next_i = i + 1
334        if next_i == len(keys):
335            next_i = i
336        # If the save state is from a later version, skip the check
337        if model_version <= keys[next_i]:
338            newname = _conversion_target(name, version)
339        else:
340            newname = None
341        # If no conversion is found, move on
342        if newname is None:
343            newname = name
344            continue
345        if ':' in newname:   # core_shell_ellipsoid:1
346            model_info = load_model_info(newname[:-2])
347            # Know the table exists and isn't multiplicity so grab it directly
348            # Can't use _get_translation_table since that will return the 'bare'
349            # version.
350            translation = CONVERSION_TABLE.get(version, {})[newname][1]
351        else:
352            model_info = load_model_info(newname)
353            translation = _get_translation_table(model_info, version)
354        newpars = _hand_convert(newname, newpars, version)
355        newpars = _convert_pars(newpars, translation)
356        # TODO: Still not convinced this is the best check
357        if not model_info.structure_factor and version == (3, 1, 2):
358            newpars = _rescale_sld(model_info, newpars, 1e6)
359        newpars.setdefault('scale', 1.0)
360        newpars.setdefault('background', 0.0)
361        if use_underscore:
362            newpars = _pd_to_underscores(newpars)
363        name = newname
364    return newname, newpars
365
366# ========= BACKWARD CONVERSION sasmodels => sasview 3.x ===========
367
368def _revert_pars(pars, mapping):
369    """
370    Rename the parameters and any associated polydispersity attributes.
371    """
372    newpars = pars.copy()
373
374    for new, old in mapping.items():
375        for underscore, dot in PD_DOT:
376            if old and old+underscore == new+dot:
377                continue
378            if new+underscore in newpars:
379                if old is not None:
380                    newpars[old+dot] = pars[new+underscore]
381                del newpars[new+underscore]
382    for k in list(newpars.keys()):
383        for underscore, dot in PD_DOT[1:]:  # skip "" => ""
384            if k.endswith(underscore):
385                newpars[k[:-len(underscore)]+dot] = newpars[k]
386                del newpars[k]
387    return newpars
388
389def revert_name(model_info):
390    oldname, _ = CONVERSION_TABLE.get(model_info.id, [None, {}])
391    return oldname
392
393def _remove_pd(pars, key, name):
394    """
395    Remove polydispersity from the parameter list.
396
397    Note: operates in place
398    """
399    # Bumps style parameter names
400    width = pars.pop(key+".width", 0.0)
401    n_points = pars.pop(key+".npts", 0)
402    if width != 0.0 and n_points != 0:
403        warnings.warn("parameter %s not polydisperse in sasview %s"%(key, name))
404    pars.pop(key+".nsigmas", None)
405    pars.pop(key+".type", None)
406    return pars
407
408def _trim_vectors(model_info, pars, oldpars):
409    _, translation = CONVERSION_TABLE.get(model_info.id, [None, {}])
410    for p in model_info.parameters.kernel_parameters:
411        if p.length_control is not None:
412            n = int(pars[p.length_control])
413            oldname = translation.get(p.id, p.id)
414            for k in range(n+1, p.length+1):
415                for _, old in PD_DOT:
416                    oldpars.pop(oldname+str(k)+old, None)
417    return oldpars
418
419def revert_pars(model_info, pars):
420    """
421    Convert model from new style parameter names to old style.
422    """
423    if model_info.composition is not None:
424        composition_type, parts = model_info.composition
425        if composition_type == 'product':
426            translation = _get_translation_table(parts[0])
427            # structure factor models include scale:scale_factor mapping
428            translation.update(_get_translation_table(parts[1]))
429        else:
430            raise NotImplementedError("cannot convert to sasview sum")
431    else:
432        translation = _get_translation_table(model_info)
433    oldpars = _revert_pars(_rescale_sld(model_info, pars, 1e-6), translation)
434    oldpars = _trim_vectors(model_info, pars, oldpars)
435
436    # Make sure the control parameter is an integer
437    if "CONTROL" in oldpars:
438        oldpars["CONTROL"] = int(oldpars["CONTROL"])
439
440    # Note: update compare.constrain_pars to match
441    name = model_info.id
442    if name in MODELS_WITHOUT_SCALE or model_info.structure_factor:
443        if oldpars.pop('scale', 1.0) != 1.0:
444            warnings.warn("parameter scale not used in sasview %s"%name)
445    if name in MODELS_WITHOUT_BACKGROUND or model_info.structure_factor:
446        if oldpars.pop('background', 0.0) != 0.0:
447            warnings.warn("parameter background not used in sasview %s"%name)
448
449    # Remove magnetic parameters from non-magnetic sasview models
450    if name not in MAGNETIC_SASVIEW_MODELS:
451        oldpars = dict((k, v) for k, v in oldpars.items() if ':' not in k)
452
453    # If it is a product model P*S, then check the individual forms for special
454    # cases.  Note: despite the structure factor alone not having scale or
455    # background, the product model does, so this is below the test for
456    # models without scale or background.
457    namelist = name.split('*') if '*' in name else [name]
458    for name in namelist:
459        if name in MODELS_WITHOUT_VOLFRACTION:
460            del oldpars['volfraction']
461        elif name == 'core_multi_shell':
462            # kill extra shells
463            for k in range(5, 11):
464                oldpars.pop('sld_shell'+str(k), 0)
465                oldpars.pop('thick_shell'+str(k), 0)
466                oldpars.pop('mtheta:sld'+str(k), 0)
467                oldpars.pop('mphi:sld'+str(k), 0)
468                oldpars.pop('M0:sld'+str(k), 0)
469                _remove_pd(oldpars, 'sld_shell'+str(k), 'sld')
470                _remove_pd(oldpars, 'thick_shell'+str(k), 'thickness')
471        elif name == 'core_shell_parallelepiped':
472            _remove_pd(oldpars, 'rimA', name)
473            _remove_pd(oldpars, 'rimB', name)
474            _remove_pd(oldpars, 'rimC', name)
475        elif name == 'hollow_cylinder':
476            # now uses radius and thickness
477            thickness = oldpars['core_radius']
478            oldpars['radius'] += thickness
479            oldpars['radius.width'] *= thickness/oldpars['radius']
480        #elif name in ['mono_gauss_coil', 'poly_gauss_coil']:
481        #    del oldpars['i_zero']
482        elif name == 'onion':
483            oldpars.pop('n_shells', None)
484        elif name == 'pearl_necklace':
485            _remove_pd(oldpars, 'num_pearls', name)
486            _remove_pd(oldpars, 'thick_string', name)
487        elif name == 'polymer_micelle':
488            if 'ndensity' in oldpars:
489                oldpars['ndensity'] *= 1e15
490        elif name == 'rpa':
491            # convert scattering lengths from femtometers to centimeters
492            for p in "L1", "L2", "L3", "L4":
493                if p in oldpars: oldpars[p] *= 1e-13
494            if pars['case_num'] < 2:
495                for k in ("a", "b"):
496                    for p in ("L", "N", "Phi", "b", "v"):
497                        oldpars.pop(p+k, None)
498                for k in "Kab,Kac,Kad,Kbc,Kbd".split(','):
499                    oldpars.pop(k, None)
500            elif pars['case_num'] < 5:
501                for k in ("a",):
502                    for p in ("L", "N", "Phi", "b", "v"):
503                        oldpars.pop(p+k, None)
504                for k in "Kab,Kac,Kad".split(','):
505                    oldpars.pop(k, None)
506        elif name == 'spherical_sld':
507            oldpars["CONTROL"] -= 1
508            # remove polydispersity from shells
509            for k in range(1, 11):
510                _remove_pd(oldpars, 'thick_flat'+str(k), 'thickness')
511                _remove_pd(oldpars, 'thick_inter'+str(k), 'interface')
512            # remove extra shells
513            for k in range(int(pars['n_shells']), 11):
514                oldpars.pop('sld_flat'+str(k), 0)
515                oldpars.pop('thick_flat'+str(k), 0)
516                oldpars.pop('thick_inter'+str(k), 0)
517                oldpars.pop('func_inter'+str(k), 0)
518                oldpars.pop('nu_inter'+str(k), 0)
519        elif name == 'stacked_disks':
520            _remove_pd(oldpars, 'n_stacking', name)
521        elif name == 'teubner_strey':
522            # basically redoing the entire Teubner-Strey calculations here.
523            volfraction = oldpars.pop('volfraction_a')
524            xi = oldpars.pop('xi')
525            d = oldpars.pop('d')
526            sld_a = oldpars.pop('sld_a')
527            sld_b = oldpars.pop('sld_b')
528            drho = 1e6*(sld_a - sld_b)  # conversion autoscaled these
529            k = 2.0*math.pi*xi/d
530            a2 = (1.0 + k**2)**2
531            c1 = 2.0 * xi**2 * (1.0 - k**2)
532            c2 = xi**4
533            prefactor = 8.0*math.pi*volfraction*(1.0-volfraction)*drho**2*c2/xi
534            scale = 1e-4*prefactor
535            oldpars['scale'] = a2/scale
536            oldpars['c1'] = c1/scale
537            oldpars['c2'] = c2/scale
538
539    #print("convert from",list(sorted(pars)))
540    #print("convert to",list(sorted(oldpars.items())))
541    return oldpars
542
543def constrain_new_to_old(model_info, pars):
544    """
545    Restrict parameter values to those that will match sasview.
546    """
547    name = model_info.id
548    # Note: update convert.revert_model to match
549    if name in MODELS_WITHOUT_SCALE or model_info.structure_factor:
550        pars['scale'] = 1
551    if name in MODELS_WITHOUT_BACKGROUND or model_info.structure_factor:
552        pars['background'] = 0
553    # sasview multiplies background by structure factor
554    if '*' in name:
555        pars['background'] = 0
556
557    # Shut off magnetism when comparing non-magnetic sasview models
558    if name not in MAGNETIC_SASVIEW_MODELS:
559        suppress_magnetism = False
560        for key in pars.keys():
561            if key.startswith("M0:"):
562                suppress_magnetism = suppress_magnetism or (pars[key] != 0)
563                pars[key] = 0
564        if suppress_magnetism:
565            warnings.warn("suppressing magnetism for comparison with sasview")
566
567    # Shut off theta polydispersity since algorithm has changed
568    if 'theta_pd_n' in pars:
569        if pars['theta_pd_n'] != 0:
570            warnings.warn("suppressing theta polydispersity for comparison with sasview")
571        pars['theta_pd_n'] = 0
572
573    # If it is a product model P*S, then check the individual forms for special
574    # cases.  Note: despite the structure factor alone not having scale or
575    # background, the product model does, so this is below the test for
576    # models without scale or background.
577    namelist = name.split('*') if '*' in name else [name]
578    for name in namelist:
579        if name in MODELS_WITHOUT_VOLFRACTION:
580            pars['volfraction'] = 1
581        if name == 'core_multi_shell':
582            pars['n'] = min(math.ceil(pars['n']), 4)
583        elif name == 'gel_fit':
584            pars['scale'] = 1
585        elif name == 'line':
586            pars['scale'] = 1
587            pars['background'] = 0
588        elif name == 'mono_gauss_coil':
589            pars['scale'] = 1
590        elif name == 'onion':
591            pars['n_shells'] = math.ceil(pars['n_shells'])
592        elif name == 'pearl_necklace':
593            pars['string_thickness_pd_n'] = 0
594            pars['number_of_pearls_pd_n'] = 0
595        elif name == 'poly_gauss_coil':
596            pars['scale'] = 1
597        elif name == 'rpa':
598            pars['case_num'] = int(pars['case_num'])
599        elif name == 'spherical_sld':
600            pars['n_shells'] = math.ceil(pars['n_shells'])
601            pars['n_steps'] = math.ceil(pars['n_steps'])
602            for k in range(1, 11):
603                pars['shape%d'%k] = math.trunc(pars['shape%d'%k]+0.5)
604            for k in range(2, 11):
605                pars['thickness%d_pd_n'%k] = 0
606                pars['interface%d_pd_n'%k] = 0
607        elif name == 'teubner_strey':
608            pars['scale'] = 1
609            if pars['volfraction_a'] > 0.5:
610                pars['volfraction_a'] = 1.0 - pars['volfraction_a']
611        elif name == 'unified_power_Rg':
612            pars['level'] = int(pars['level'])
613
614def _check_one(name, seed=None):
615    """
616    Generate a random set of parameters for *name*, and check that they can
617    be converted back to SasView 3.x and forward again to sasmodels.  Raises
618    an error if the parameters are changed.
619    """
620    from . import compare
621
622    model_info = load_model_info(name)
623
624    old_name = revert_name(model_info)
625    if old_name is None:
626        return
627
628    pars = compare.get_pars(model_info, use_demo=False)
629    if seed is not None:
630        np.random.seed(seed)
631    pars = compare.randomize_pars(model_info, pars)
632    if name == "teubner_strey":
633        # T-S model is underconstrained, so fix the assumptions.
634        pars['sld_a'], pars['sld_b'] = 1.0, 0.0
635    compare.constrain_pars(model_info, pars)
636    constrain_new_to_old(model_info, pars)
637    old_pars = revert_pars(model_info, pars)
638    new_name, new_pars = convert_model(old_name, old_pars, use_underscore=True)
639    if 1:
640        print("==== %s in ====="%name)
641        print(str(compare.parlist(model_info, pars, True)))
642        print("==== %s ====="%old_name)
643        for k, v in sorted(old_pars.items()):
644            print(k, v)
645        print("==== %s out ====="%new_name)
646        print(str(compare.parlist(model_info, new_pars, True)))
647    assert name == new_name, "%r != %r"%(name, new_name)
648    for k, v in new_pars.items():
649        assert k in pars, "%s: %r appeared from conversion"%(name, k)
650        if isinstance(v, float):
651            assert abs(v-pars[k]) <= abs(1e-12*v), \
652                "%s: %r  %s != %s"%(name, k, v, pars[k])
653        else:
654            assert v == pars[k], "%s: %r  %s != %s"%(name, k, v, pars[k])
655    for k, v in pars.items():
656        assert k in pars, "%s: %r not converted"%(name, k)
657
658def test_backward_forward():
659    from .core import list_models
660    L = lambda name: _check_one(name, seed=1)
661    for name in list_models('all'):
662        yield L, name
Note: See TracBrowser for help on using the repository browser.