""" Sasview model constructor. Given a module defining an OpenCL kernel such as sasmodels.models.cylinder, create a sasview model class to run that kernel as follows:: from sasmodels.sasview_model import make_class from sasmodels.models import cylinder CylinderModel = make_class(cylinder, dtype='single') The model parameters for sasmodels are different from those in sasview. When reloading previously saved models, the parameters should be converted using :func:`sasmodels.convert.convert`. """ import math from copy import deepcopy import warnings import numpy as np from . import core def make_class(model_info, dtype='single', namestyle='name'): """ Load the sasview model defined in *kernel_module*. Returns a class that can be used directly as a sasview model. Defaults to using the new name for a model. Setting *namestyle='oldname'* will produce a class with a name compatible with SasView. """ model = core.build_model(model_info, dtype=dtype) def __init__(self, multfactor=1): SasviewModel.__init__(self, model) attrs = dict(__init__=__init__) ConstructedModel = type(model.info[namestyle], (SasviewModel,), attrs) return ConstructedModel class SasviewModel(object): """ Sasview wrapper for opencl/ctypes model. """ def __init__(self, model): """Initialization""" self._model = model self.name = model.info['name'] self.oldname = model.info['oldname'] self.description = model.info['description'] self.category = None self.multiplicity_info = None self.is_multifunc = False ## interpret the parameters ## TODO: reorganize parameter handling self.details = dict() self.params = dict() self.dispersion = dict() partype = model.info['partype'] for name, units, default, limits, _, _ in model.info['parameters']: self.params[name] = default self.details[name] = [units] + limits for name in partype['pd-2d']: self.dispersion[name] = { 'width': 0, 'npts': 35, 'nsigmas': 3, 'type': 'gaussian', } self.orientation_params = ( partype['orientation'] + [n + '.width' for n in partype['orientation']] + partype['magnetic']) self.magnetic_params = partype['magnetic'] self.fixed = [n + '.width' for n in partype['pd-2d']] self.non_fittable = [] ## independent parameter name and unit [string] self.input_name = model.info.get("input_name", "Q") self.input_unit = model.info.get("input_unit", "A^{-1}") self.output_name = model.info.get("output_name", "Intensity") self.output_unit = model.info.get("output_unit", "cm^{-1}") ## _persistency_dict is used by sas.perspectives.fitting.basepage ## to store dispersity reference. ## TODO: _persistency_dict to persistency_dict throughout sasview self._persistency_dict = {} ## New fields introduced for opencl rewrite self.cutoff = 1e-5 def __str__(self): """ :return: string representation """ return self.name def is_fittable(self, par_name): """ Check if a given parameter is fittable or not :param par_name: the parameter name to check """ return par_name.lower() in self.fixed #For the future #return self.params[str(par_name)].is_fittable() # pylint: disable=no-self-use def getProfile(self): """ Get SLD profile : return: (z, beta) where z is a list of depth of the transition points beta is a list of the corresponding SLD values """ return None, None def setParam(self, name, value): """ Set the value of a model parameter :param name: name of the parameter :param value: value of the parameter """ # Look for dispersion parameters toks = name.split('.') if len(toks) == 2: for item in self.dispersion.keys(): if item.lower() == toks[0].lower(): for par in self.dispersion[item]: if par.lower() == toks[1].lower(): self.dispersion[item][par] = value return else: # Look for standard parameter for item in self.params.keys(): if item.lower() == name.lower(): self.params[item] = value return raise ValueError("Model does not contain parameter %s" % name) def getParam(self, name): """ Set the value of a model parameter :param name: name of the parameter """ # Look for dispersion parameters toks = name.split('.') if len(toks) == 2: for item in self.dispersion.keys(): if item.lower() == toks[0].lower(): for par in self.dispersion[item]: if par.lower() == toks[1].lower(): return self.dispersion[item][par] else: # Look for standard parameter for item in self.params.keys(): if item.lower() == name.lower(): return self.params[item] raise ValueError("Model does not contain parameter %s" % name) def getParamList(self): """ Return a list of all available parameters for the model """ param_list = self.params.keys() # WARNING: Extending the list with the dispersion parameters param_list.extend(self.getDispParamList()) return param_list def getDispParamList(self): """ Return a list of all available parameters for the model """ # TODO: fix test so that parameter order doesn't matter ret = ['%s.%s' % (d.lower(), p) for d in self._model.info['partype']['pd-2d'] for p in ('npts', 'nsigmas', 'width')] #print(ret) return ret def clone(self): """ Return a identical copy of self """ return deepcopy(self) def run(self, x=0.0): """ Evaluate the model :param x: input q, or [q,phi] :return: scattering function P(q) **DEPRECATED**: use calculate_Iq instead """ if isinstance(x, (list, tuple)): # pylint: disable=unpacking-non-sequence q, phi = x return self.calculate_Iq([q * math.cos(phi)], [q * math.sin(phi)])[0] else: return self.calculate_Iq([float(x)])[0] def runXY(self, x=0.0): """ Evaluate the model in cartesian coordinates :param x: input q, or [qx, qy] :return: scattering function P(q) **DEPRECATED**: use calculate_Iq instead """ if isinstance(x, (list, tuple)): return self.calculate_Iq([float(x[0])], [float(x[1])])[0] else: return self.calculate_Iq([float(x)])[0] def evalDistribution(self, qdist): r""" Evaluate a distribution of q-values. :param qdist: array of q or a list of arrays [qx,qy] * For 1D, a numpy array is expected as input :: evalDistribution(q) where *q* is a numpy array. * For 2D, a list of *[qx,qy]* is expected with 1D arrays as input :: qx = [ qx[0], qx[1], qx[2], ....] qy = [ qy[0], qy[1], qy[2], ....] If the model is 1D only, then .. math:: q = \sqrt{q_x^2+q_y^2} """ if isinstance(qdist, (list, tuple)): # Check whether we have a list of ndarrays [qx,qy] qx, qy = qdist partype = self._model.info['partype'] if not partype['orientation'] and not partype['magnetic']: return self.calculate_Iq(np.sqrt(qx ** 2 + qy ** 2)) else: return self.calculate_Iq(qx, qy) elif isinstance(qdist, np.ndarray): # We have a simple 1D distribution of q-values return self.calculate_Iq(qdist) else: raise TypeError("evalDistribution expects q or [qx, qy], not %r" % type(qdist)) def calculate_Iq(self, *args): """ Calculate Iq for one set of q with the current parameters. If the model is 1D, use *q*. If 2D, use *qx*, *qy*. This should NOT be used for fitting since it copies the *q* vectors to the card for each evaluation. """ q_vectors = [np.asarray(q) for q in args] fn = self._model(q_vectors) pars = [self.params[v] for v in fn.fixed_pars] pd_pars = [self._get_weights(p) for p in fn.pd_pars] result = fn(pars, pd_pars, self.cutoff) fn.input.release() fn.release() return result def calculate_ER(self): """ Calculate the effective radius for P(q)*S(q) :return: the value of the effective radius """ ER = self._model.info.get('ER', None) if ER is None: return 1.0 else: values, weights = self._dispersion_mesh() fv = ER(*values) #print(values[0].shape, weights.shape, fv.shape) return np.sum(weights * fv) / np.sum(weights) def calculate_VR(self): """ Calculate the volf ratio for P(q)*S(q) :return: the value of the volf ratio """ VR = self._model.info.get('VR', None) if VR is None: return 1.0 else: values, weights = self._dispersion_mesh() whole, part = VR(*values) return np.sum(weights * part) / np.sum(weights * whole) def set_dispersion(self, parameter, dispersion): """ Set the dispersion object for a model parameter :param parameter: name of the parameter [string] :param dispersion: dispersion object of type Dispersion """ if parameter.lower() in (s.lower() for s in self.params.keys()): # TODO: Store the disperser object directly in the model. # The current method of creating one on the fly whenever it is # needed is kind of funky. # Note: can't seem to get disperser parameters from sasview # (1) Could create a sasview model that has not yet # been # converted, assign the disperser to one of its polydisperse # parameters, then retrieve the disperser parameters from the # sasview model. (2) Could write a disperser parameter retriever # in sasview. (3) Could modify sasview to use sasmodels.weights # dispersers. # For now, rely on the fact that the sasview only ever uses # new dispersers in the set_dispersion call and create a new # one instead of trying to assign parameters. from . import weights disperser = weights.dispersers[dispersion.__class__.__name__] dispersion = weights.models[disperser]() self.dispersion[parameter] = dispersion.get_pars() else: raise ValueError("%r is not a dispersity or orientation parameter") def _dispersion_mesh(self): """ Create a mesh grid of dispersion parameters and weights. Returns [p1,p2,...],w where pj is a vector of values for parameter j and w is a vector containing the products for weights for each parameter set in the vector. """ pars = self._model.info['partype']['volume'] return core.dispersion_mesh([self._get_weights(p) for p in pars]) def _get_weights(self, par): """ Return dispersion weights :param par parameter name """ from . import weights relative = self._model.info['partype']['pd-rel'] limits = self._model.info['limits'] dis = self.dispersion[par] value, weight = weights.get_weights( dis['type'], dis['npts'], dis['width'], dis['nsigmas'], self.params[par], limits[par], par in relative) return value, weight / np.sum(weight)