source: sasmodels/sasmodels/convert.py @ 1b49bf8

core_shell_microgelscostrafo411magnetic_modelticket-1257-vesicle-productticket_1156ticket_1265_superballticket_822_more_unit_tests
Last change on this file since 1b49bf8 was 1b49bf8, checked in by Paul Kienzle <pkienzle@…>, 8 years ago

update conversion tables corresponding to the great parameter shift
  • Property mode set to 100644
File size: 13.4 KB
Line 
1"""
2Convert models to and from sasview.
3"""
4from __future__ import print_function, division
5
6import math
7import warnings
8
9from .conversion_table import CONVERSION_TABLE
10
11# List of models which SasView versions don't contain the explicit 'scale' argument.
12# When converting such a model, please update this list.
13MODELS_WITHOUT_SCALE = [
14    'teubner_strey',
15    'broad_peak',
16    'two_lorentzian',
17    "two_power_law",
18    'gauss_lorentz_gel',
19    'be_polyelectrolyte',
20    'correlation_length',
21    'fractal_core_shell',
22    'binary_hard_sphere',
23    'raspberry'
24]
25
26# List of models which SasView versions don't contain the explicit 'background' argument.
27# When converting such a model, please update this list.
28MODELS_WITHOUT_BACKGROUND = [
29    'guinier',
30]
31
32MODELS_WITHOUT_VOLFRACTION = [
33    'fractal',
34    'vesicle',
35    'multilayer_vesicle',
36]
37
38MAGNETIC_SASVIEW_MODELS = [
39    'core_shell',
40    'core_multi_shell',
41    'cylinder',
42    'parallelepiped',
43    'sphere',
44]
45
46
47# Convert new style names for polydispersity info to old style names
48PD_DOT = [
49    ("", ""),
50    ("_pd", ".width"),
51    ("_pd_n", ".npts"),
52    ("_pd_nsigma", ".nsigmas"),
53    ("_pd_type", ".type"),
54    ]
55
56def _convert_pars(pars, mapping):
57    """
58    Rename the parameters and any associated polydispersity attributes.
59    """
60    newpars = pars.copy()
61    for new, old in mapping.items():
62        if old == new: continue
63        for underscore, dot in PD_DOT:
64            if old+dot in newpars:
65                if new is not None:
66                    newpars[new+underscore] = pars[old+dot]
67                del newpars[old+dot]
68    return newpars
69
70def convert_model(name, pars):
71    """
72    Convert model from old style parameter names to new style.
73    """
74    _, _ = name, pars # lint
75    raise NotImplementedError
76    # need to load all new models in order to determine old=>new
77    # model name mapping
78
79def _unscale(par, scale):
80    return [pk*scale for pk in par] if isinstance(par, list) else par*scale
81
82def _is_sld(modelinfo, id):
83    if id.startswith('M0:'):
84        return True
85    if (id.endswith('_pd') or id.endswith('_pd_n') or id.endswith('_pd_nsigma')
86            or id.endswith('_pd_width') or id.endswith('_pd_type')):
87        return False
88    for p in modelinfo.parameters.call_parameters:
89        if p.id == id:
90            return p.type == 'sld'
91    # check through kernel parameters in case it is a named as a vector
92    for p in modelinfo.parameters.kernel_parameters:
93        if p.id == id:
94            return p.type == 'sld'
95    raise ValueError("unknown parameter %r in conversion"%id)
96
97def _unscale_sld(modelinfo, pars):
98    """
99    rescale all sld parameters in the new model definition by 1e6 so the
100    numbers are nicer.  Relies on the fact that all sld parameters in the
101    new model definition end with sld.
102    """
103    return dict((id, (_unscale(v, 1e-6) if _is_sld(modelinfo, id) else v))
104                for id, v in pars.items())
105
106def _remove_pd(pars, key, name):
107    """
108    Remove polydispersity from the parameter list.
109
110    Note: operates in place
111    """
112    # Bumps style parameter names
113    width = pars.pop(key+".width", 0.0)
114    n_points = pars.pop(key+".npts", 0)
115    if width != 0.0 and n_points != 0:
116        warnings.warn("parameter %s not polydisperse in sasview %s"%(key, name))
117    pars.pop(key+".nsigmas", None)
118    pars.pop(key+".type", None)
119    return pars
120
121def _revert_pars(pars, mapping):
122    """
123    Rename the parameters and any associated polydispersity attributes.
124    """
125    newpars = pars.copy()
126
127    for new, old in mapping.items():
128        for underscore, dot in PD_DOT:
129            if old and old+underscore == new+dot:
130                continue
131            if new+underscore in newpars:
132                if old is not None:
133                    newpars[old+dot] = pars[new+underscore]
134                del newpars[new+underscore]
135    for k in list(newpars.keys()):
136        for underscore, dot in PD_DOT[1:]:  # skip "" => ""
137            if k.endswith(underscore):
138                newpars[k[:-len(underscore)]+dot] = newpars[k]
139                del newpars[k]
140    return newpars
141
142def revert_name(model_info):
143    oldname, _ = CONVERSION_TABLE.get(model_info.id, [None, {}])
144    return oldname
145
146def _get_translation_table(model_info):
147    _, translation = CONVERSION_TABLE.get(model_info.id, [None, {}])
148    translation = translation.copy()
149    for p in model_info.parameters.kernel_parameters:
150        if p.length > 1:
151            newid = p.id
152            oldid = translation.get(p.id, p.id)
153            translation.pop(newid, None)
154            for k in range(1, p.length+1):
155                if newid+str(k) not in translation:
156                    translation[newid+str(k)] = oldid+str(k)
157    # Remove control parameter from the result
158    if model_info.control:
159        translation[model_info.control] = "CONTROL"
160    return translation
161
162def _trim_vectors(model_info, pars, oldpars):
163    _, translation = CONVERSION_TABLE.get(model_info.id, [None, {}])
164    for p in model_info.parameters.kernel_parameters:
165        if p.length_control is not None:
166            n = int(pars[p.length_control])
167            oldname = translation.get(p.id, p.id)
168            for k in range(n+1, p.length+1):
169                for _, old in PD_DOT:
170                    oldpars.pop(oldname+str(k)+old, None)
171    return oldpars
172
173def revert_pars(model_info, pars):
174    """
175    Convert model from new style parameter names to old style.
176    """
177    if model_info.composition is not None:
178        composition_type, parts = model_info.composition
179        if composition_type == 'product':
180            translation = {'scale':'scale_factor'}
181            translation.update(_get_translation_table(parts[0]))
182            translation.update(_get_translation_table(parts[1]))
183        else:
184            raise NotImplementedError("cannot convert to sasview sum")
185    else:
186        translation = _get_translation_table(model_info)
187    oldpars = _revert_pars(_unscale_sld(model_info, pars), translation)
188    oldpars = _trim_vectors(model_info, pars, oldpars)
189
190    # Make sure the control parameter is an integer
191    if "CONTROL" in oldpars:
192        oldpars["CONTROL"] = int(oldpars["CONTROL"])
193
194    # Note: update compare.constrain_pars to match
195    name = model_info.id
196    if name in MODELS_WITHOUT_SCALE or model_info.structure_factor:
197        if oldpars.pop('scale', 1.0) != 1.0:
198            warnings.warn("parameter scale not used in sasview %s"%name)
199    if name in MODELS_WITHOUT_BACKGROUND or model_info.structure_factor:
200        if oldpars.pop('background', 0.0) != 0.0:
201            warnings.warn("parameter background not used in sasview %s"%name)
202
203    # Remove magnetic parameters from non-magnetic sasview models
204    if name not in MAGNETIC_SASVIEW_MODELS:
205        oldpars = dict((k, v) for k, v in oldpars.items() if ':' not in k)
206
207    # If it is a product model P*S, then check the individual forms for special
208    # cases.  Note: despite the structure factor alone not having scale or
209    # background, the product model does, so this is below the test for
210    # models without scale or background.
211    namelist = name.split('*') if '*' in name else [name]
212    for name in namelist:
213        if name in MODELS_WITHOUT_VOLFRACTION:
214            del oldpars['volfraction']
215        elif name == 'core_multi_shell':
216            # kill extra shells
217            for k in range(5, 11):
218                oldpars.pop('sld_shell'+str(k), 0)
219                oldpars.pop('thick_shell'+str(k), 0)
220                oldpars.pop('mtheta:sld'+str(k), 0)
221                oldpars.pop('mphi:sld'+str(k), 0)
222                oldpars.pop('M0:sld'+str(k), 0)
223                _remove_pd(oldpars, 'sld_shell'+str(k), 'sld')
224                _remove_pd(oldpars, 'thick_shell'+str(k), 'thickness')
225        elif name == 'core_shell_parallelepiped':
226            _remove_pd(oldpars, 'rimA', name)
227            _remove_pd(oldpars, 'rimB', name)
228            _remove_pd(oldpars, 'rimC', name)
229        elif name == 'hollow_cylinder':
230            # now uses radius and thickness
231            oldpars['radius'] += oldpars['core_radius']
232        elif name in ['mono_gauss_coil', 'poly_gauss_coil']:
233            del oldpars['i_zero']
234        elif name == 'onion':
235            oldpars.pop('n_shells', None)
236        elif name == 'pearl_necklace':
237            _remove_pd(oldpars, 'num_pearls', name)
238            _remove_pd(oldpars, 'thick_string', name)
239        elif name == 'polymer_micelle':
240            if 'ndensity' in oldpars:
241                oldpars['ndensity'] *= 1e15
242        elif name == 'rpa':
243            # convert scattering lengths from femtometers to centimeters
244            for p in "L1", "L2", "L3", "L4":
245                if p in oldpars: oldpars[p] *= 1e-13
246            if pars['case_num'] < 2:
247                for k in ("a", "b"):
248                    for p in ("L", "N", "Phi", "b", "v"):
249                        oldpars.pop(p+k, None)
250                for k in "Kab,Kac,Kad,Kbc,Kbd".split(','):
251                    oldpars.pop(k, None)
252            elif pars['case_num'] < 5:
253                for k in ("a",):
254                    for p in ("L", "N", "Phi", "b", "v"):
255                        oldpars.pop(p+k, None)
256                for k in "Kab,Kac,Kad".split(','):
257                    oldpars.pop(k, None)
258        elif name == 'spherical_sld':
259            oldpars["CONTROL"] -= 1
260            # remove polydispersity from shells
261            for k in range(1, 11):
262                _remove_pd(oldpars, 'thick_flat'+str(k), 'thickness')
263                _remove_pd(oldpars, 'thick_inter'+str(k), 'interface')
264            # remove extra shells
265            for k in range(int(pars['n_shells']), 11):
266                oldpars.pop('sld_flat'+str(k), 0)
267                oldpars.pop('thick_flat'+str(k), 0)
268                oldpars.pop('thick_inter'+str(k), 0)
269                oldpars.pop('func_inter'+str(k), 0)
270                oldpars.pop('nu_inter'+str(k), 0)
271        elif name == 'stacked_disks':
272            _remove_pd(oldpars, 'n_stacking', name)
273        elif name == 'teubner_strey':
274            # basically redoing the entire Teubner-Strey calculations here.
275            volfraction = oldpars.pop('volfraction_a')
276            xi = oldpars.pop('xi')
277            d = oldpars.pop('d')
278            sld_a = oldpars.pop('sld_a')
279            sld_b = oldpars.pop('sld_b')
280            drho = 1e6*(sld_a - sld_b)  # conversion autoscaled these
281            k = 2.0*math.pi*xi/d
282            a2 = (1.0 + k**2)**2
283            c1 = 2.0 * xi**2 * (1.0 - k**2)
284            c2 = xi**4
285            prefactor = 8.0*math.pi*volfraction*(1.0-volfraction)*drho**2*c2/xi
286            scale = 1e-4*prefactor
287            oldpars['scale'] = a2/scale
288            oldpars['c1'] = c1/scale
289            oldpars['c2'] = c2/scale
290
291    #print("convert from",list(sorted(pars)))
292    #print("convert to",list(sorted(oldpars.items())))
293    return oldpars
294
295def constrain_new_to_old(model_info, pars):
296    """
297    Restrict parameter values to those that will match sasview.
298    """
299    name = model_info.id
300    # Note: update convert.revert_model to match
301    if name in MODELS_WITHOUT_SCALE or model_info.structure_factor:
302        pars['scale'] = 1
303    if name in MODELS_WITHOUT_BACKGROUND or model_info.structure_factor:
304        pars['background'] = 0
305    # sasview multiplies background by structure factor
306    if '*' in name:
307        pars['background'] = 0
308
309    # Shut off magnetism when comparing non-magnetic sasview models
310    if name not in MAGNETIC_SASVIEW_MODELS:
311        suppress_magnetism = False
312        for key in pars.keys():
313            if key.startswith("M0:"):
314                suppress_magnetism = suppress_magnetism or (pars[key] != 0)
315                pars[key] = 0
316        if suppress_magnetism:
317            warnings.warn("suppressing magnetism for comparison with sasview")
318
319    # Shut off theta polydispersity since algorithm has changed
320    if 'theta_pd_n' in pars:
321        if pars['theta_pd_n'] != 0:
322            warnings.warn("suppressing theta polydispersity for comparison with sasview")
323        pars['theta_pd_n'] = 0
324
325    # If it is a product model P*S, then check the individual forms for special
326    # cases.  Note: despite the structure factor alone not having scale or
327    # background, the product model does, so this is below the test for
328    # models without scale or background.
329    namelist = name.split('*') if '*' in name else [name]
330    for name in namelist:
331        if name in MODELS_WITHOUT_VOLFRACTION:
332            pars['volfraction'] = 1
333        if name == 'core_multi_shell':
334            pars['n'] = min(math.ceil(pars['n']), 4)
335        elif name == 'gel_fit':
336            pars['scale'] = 1
337        elif name == 'line':
338            pars['scale'] = 1
339            pars['background'] = 0
340        elif name == 'mono_gauss_coil':
341            pars['i_zero'] = 1
342        elif name == 'onion':
343            pars['n_shells'] = math.ceil(pars['n_shells'])
344        elif name == 'pearl_necklace':
345            pars['string_thickness_pd_n'] = 0
346            pars['number_of_pearls_pd_n'] = 0
347        elif name == 'poly_gauss_coil':
348            pars['i_zero'] = 1
349        elif name == 'rpa':
350            pars['case_num'] = int(pars['case_num'])
351        elif name == 'spherical_sld':
352            pars['n_shells'] = math.ceil(pars['n_shells'])
353            pars['n_steps'] = math.ceil(pars['n_steps'])
354            for k in range(1, 11):
355                pars['shape%d'%k] = math.trunc(pars['shape%d'%k]+0.5)
356            for k in range(2, 11):
357                pars['thickness%d_pd_n'%k] = 0
358                pars['interface%d_pd_n'%k] = 0
359        elif name == 'teubner_strey':
360            pars['scale'] = 1
361
Note: See TracBrowser for help on using the repository browser.