source: sasmodels/sasmodels/convert.py @ 6592f56

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

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