Changes in / [8696a87:504abee] in sasmodels


Ignore:
Files:
1 deleted
18 edited

Legend:

Unmodified
Added
Removed

  • extra/pylint.rc

    r823e620 rd4666ca  
    1010 
    1111# Profiled execution. 
    12 #profile=no 
     12profile=no 
    1313 
    1414# Add files or directories to the blacklist. They should be base names, not 
    … …  
    9292# Add a comment according to your evaluation note. This is used by the global 
    9393# evaluation report (RP0004). 
    94 #comment=no 
     94comment=no 
    9595 
    9696# Template used to display messages. This is a python new-style format string 
    … …  
    101101 
    102102# Required attributes for module, separated by a comma 
    103 #required-attributes= 
     103required-attributes= 
    104104 
    105105# List of builtins function names that should not be used, separated by a comma 
    … …  
    108108# Good variable names which should always be accepted, separated by a comma 
    109109good-names=_, 
    110     input, id, 
    111     h, i, j, k, n, p, v, w, x, y, z, 
    112     q, qx, qy, qz, Q, Qx, Qy, Qz, 
    113     dt, dx, dy, dz, dq, 
    114     Iq, dIq, Iqxy, Iq_calc, 
     110    input, 
     111    h,i,j,k,n,w,x,y,z, 
     112    q,qx,qy,qz, 
     113    dt,dx,dy,dz,id, 
     114    Iq,dIq,Qx,Qy,Qz,Iqxy, 
     115    p, 
    115116    ER, call_ER, VR, call_VR, 
    116117    Rmax, SElength, 
    … …  
    280281# (useful for modules/projects where namespaces are manipulated during runtime 
    281282# and thus existing member attributes cannot be deduced by static analysis 
    282 ignored-modules=numpy,np,numpy.random, 
    283     bumps,sas, 
     283ignored-modules=numpy,np 
    284284 
    285285# List of classes names for which member attributes should not be checked 
    … …  
    289289# When zope mode is activated, add a predefined set of Zope acquired attributes 
    290290# to generated-members. 
    291 #zope=no 
     291zope=no 
    292292 
    293293# List of members which are set dynamically and missed by pylint inference 
    … …  
    319319# List of interface methods to ignore, separated by a comma. This is used for 
    320320# instance to not check methods defines in Zope's Interface base class. 
    321 #ignore-iface-methods=isImplementedBy,deferred,extends,names,namesAndDescriptions,queryDescriptionFor,getBases,getDescriptionFor,getDoc,getName,getTaggedValue,getTaggedValueTags,isEqualOrExtendedBy,setTaggedValue,isImplementedByInstancesOf,adaptWith,is_implemented_by 
     321ignore-iface-methods=isImplementedBy,deferred,extends,names,namesAndDescriptions,queryDescriptionFor,getBases,getDescriptionFor,getDoc,getName,getTaggedValue,getTaggedValueTags,isEqualOrExtendedBy,setTaggedValue,isImplementedByInstancesOf,adaptWith,is_implemented_by 
    322322 
    323323# List of method names used to declare (i.e. assign) instance attributes. 

  • extra/pylint_numpy.py

    r823e620 r63b32bb  
    88    #print("processing",module.name) 
    99    if module.name.startswith('numpy'): 
    10         if module.name == 'numpy': 
    11             import numpy 
    12         elif module.name == 'numpy.random': 
    13             import numpy.random 
    14             from numpy.random import seed, get_state, set_state 
     10        if module.name == 'numpy':  import numpy 
     11        elif module.name == 'numpy.random': import numpy.random 
    1512 

  • extra/pylint_pyopencl.py

    r823e620 r3c56da87  
    99    if module.name == 'pyopencl': 
    1010        import pyopencl 
    11         import pyopencl as cl 

  • extra/run-pylint.py

    r823e620 rf01e53d  
    66 
    77def main(): 
    8     extra = abspath(dirname(__file__)) 
    9     root = abspath(joinpath(extra, '..')) 
    10  
    118    envpath = os.environ.get('PYTHONPATH',None) 
    129    path = [envpath] if envpath else [] 
    13     path.append(extra) 
    14  
    15     #bumps = abspath(joinpath(root, "..", "bumps")) 
    16     #periodictable = abspath(joinpath(root, "..", "periodictable")) 
    17     #sasview = abspath(joinpath(root, "..", "sasview", "src")) 
    18     #path.extend((bumps, periodictable, sasview)) 
    19  
     10    path.append(abspath(dirname(__file__))) # so we can find the plugins 
    2011    os.environ['PYTHONPATH'] = ':'.join(path) 
    21  
    22     # Run the lint command 
     12    root = abspath(joinpath(dirname(__file__), '..')) 
     13    os.chdir(root) 
    2314    cmd = "pylint --rcfile extra/pylint.rc -f parseable sasmodels" 
    24     os.chdir(root) 
    2515    status = os.system(cmd) 
    2616    sys.exit(status) 

  • sasmodels/compare.py

    r190fc2b r190fc2b  
    604604def columnize(L, indent="", width=79): 
    605605    """ 
    606     Format a list of strings into columns. 
    607  
    608     Returns a string with carriage returns ready for printing. 
     606    Format a list of strings into columns for printing. 
    609607    """ 
    610608    column_width = max(len(w) for w in L) + 1 

  • sasmodels/core.py

    reafc9fa r190fc2b  
    3838    """ 
    3939    Load a model definition given the model name. 
    40  
    41     This returns a handle to the module defining the model.  This can be 
    42     used with functions in generate to build the docs or extract model info. 
    4340    """ 
    4441    __import__('sasmodels.models.'+model_name) 
    … …  
    118115    Return a computation kernel from the model definition and the q input. 
    119116    """ 
    120     return model(q_vectors) 
     117    model_input = model.make_input(q_vectors) 
     118    return model(model_input) 
    121119 
    122120def get_weights(info, pars, name): 

  • sasmodels/data.py

    r5c962df r69ec80f  
    8787 
    8888class Data1D(object): 
    89     """ 
    90     1D data object. 
    91  
    92     Note that this definition matches the attributes from sasview, with 
    93     some generic 1D data vectors and some SAS specific definitions.  Some 
    94     refactoring to allow consistent naming conventions between 1D, 2D and 
    95     SESANS data would be helpful. 
    96  
    97     **Attributes** 
    98  
    99     *x*, *dx*: $q$ vector and gaussian resolution 
    100  
    101     *y*, *dy*: $I(q)$ vector and measurement uncertainty 
    102  
    103     *mask*: values to include in plotting/analysis 
    104  
    105     *dxl*: slit widths for slit smeared data, with *dx* ignored 
    106  
    107     *qmin*, *qmax*: range of $q$ values in *x* 
    108  
    109     *filename*: label for the data line 
    110  
    111     *_xaxis*, *_xunit*: label and units for the *x* axis 
    112  
    113     *_yaxis*, *_yunit*: label and units for the *y* axis 
    114     """ 
    11589    def __init__(self, x=None, y=None, dx=None, dy=None): 
    11690        self.x, self.y, self.dx, self.dy = x, y, dx, dy 
    … …  
    11993        self.qmin = x.min() if x is not None else np.NaN 
    12094        self.qmax = x.max() if x is not None else np.NaN 
    121         # TODO: why is 1D mask False and 2D mask True? 
    122         self.mask = (np.isnan(y) if y is not None 
    123                      else np.zeros_like(x, 'b') if x is not None 
    124                      else None) 
     95        self.mask = np.isnan(y) if y is not None else None 
    12596        self._xaxis, self._xunit = "x", "" 
    12697        self._yaxis, self._yunit = "y", "" 
    … …  
    143114 
    144115class Data2D(object): 
    145     """ 
    146     2D data object. 
    147  
    148     Note that this definition matches the attributes from sasview. Some 
    149     refactoring to allow consistent naming conventions between 1D, 2D and 
    150     SESANS data would be helpful. 
    151  
    152     **Attributes** 
    153  
    154     *qx_data*, *dqx_data*: $q_x$ matrix and gaussian resolution 
    155  
    156     *qy_data*, *dqy_data*: $q_y$ matrix and gaussian resolution 
    157  
    158     *data*, *err_data*: $I(q)$ matrix and measurement uncertainty 
    159  
    160     *mask*: values to exclude from plotting/analysis 
    161  
    162     *qmin*, *qmax*: range of $q$ values in *x* 
    163  
    164     *filename*: label for the data line 
    165  
    166     *_xaxis*, *_xunit*: label and units for the *x* axis 
    167  
    168     *_yaxis*, *_yunit*: label and units for the *y* axis 
    169  
    170     *_zaxis*, *_zunit*: label and units for the *y* axis 
    171  
    172     *Q_unit*, *I_unit*: units for Q and intensity 
    173  
    174     *x_bins*, *y_bins*: grid steps in *x* and *y* directions 
    175     """ 
    176116    def __init__(self, x=None, y=None, z=None, dx=None, dy=None, dz=None): 
    177117        self.qx_data, self.dqx_data = x, dx 
    178118        self.qy_data, self.dqy_data = y, dy 
    179119        self.data, self.err_data = z, dz 
    180         self.mask = (~np.isnan(z) if z is not None 
    181                      else np.ones_like(x) if x is not None 
    182                      else None) 
     120        self.mask = ~np.isnan(z) if z is not None else None 
    183121        self.q_data = np.sqrt(x**2 + y**2) 
    184122        self.qmin = 1e-16 
    … …  
    188126        self.Q_unit = "1/A" 
    189127        self.I_unit = "1/cm" 
    190         self.xaxis("Q_x", "1/A") 
    191         self.yaxis("Q_y", "1/A") 
    192         self.zaxis("Intensity", "1/cm") 
     128        self.xaxis("Q_x", "A^{-1}") 
     129        self.yaxis("Q_y", "A^{-1}") 
     130        self.zaxis("Intensity", r"\text{cm}^{-1}") 
    193131        self._xaxis, self._xunit = "x", "" 
    194132        self._yaxis, self._yunit = "y", "" 
    … …  
    219157 
    220158class Vector(object): 
    221     """ 
    222     3-space vector of *x*, *y*, *z* 
    223     """ 
    224159    def __init__(self, x=None, y=None, z=None): 
    225160        self.x, self.y, self.z = x, y, z 
    … …  
    300235    """ 
    301236    Plot data loaded by the sasview loader. 
    302  
    303     *data* is a sasview data object, either 1D, 2D or SESANS. 
    304  
    305     *view* is log or linear. 
    306  
    307     *limits* sets the intensity limits on the plot; if None then the limits 
    308     are inferred from the data. 
    309237    """ 
    310238    # Note: kind of weird using the plot result functions to plot just the 
    … …  
    321249def plot_theory(data, theory, resid=None, view='log', 
    322250                use_data=True, limits=None): 
    323     """ 
    324     Plot theory calculation. 
    325  
    326     *data* is needed to define the graph properties such as labels and 
    327     units, and to define the data mask. 
    328  
    329     *theory* is a matrix of the same shape as the data. 
    330  
    331     *view* is log or linear 
    332  
    333     *use_data* is True if the data should be plotted as well as the theory. 
    334  
    335     *limits* sets the intensity limits on the plot; if None then the limits 
    336     are inferred from the data. 
    337     """ 
    338251    if hasattr(data, 'lam'): 
    339252        _plot_result_sesans(data, theory, resid, use_data=True, limits=limits) 
    … …  
    345258 
    346259def protect(fn): 
    347     """ 
    348     Decorator to wrap calls in an exception trapper which prints the 
    349     exception and continues.  Keyboard interrupts are ignored. 
    350     """ 
    351260    def wrapper(*args, **kw): 
    352         """ 
    353         Trap and print errors from function. 
    354         """ 
    355261        try: 
    356262            return fn(*args, **kw) 
    357         except KeyboardInterrupt: 
    358             raise 
    359263        except: 
    360264            traceback.print_exc() 
    … …  
    428332@protect 
    429333def _plot_result_sesans(data, theory, resid, use_data, limits=None): 
    430     """ 
    431     Plot SESANS results. 
    432     """ 
    433334    import matplotlib.pyplot as plt 
    434335    use_data = use_data and data.y is not None 
    … …  
    558459 
    559460def demo(): 
    560     """ 
    561     Load and plot a SAS dataset. 
    562     """ 
    563461    data = load_data('DEC07086.DAT') 
    564462    set_beam_stop(data, 0.004) 

  • sasmodels/direct_model.py

    reafc9fa r9404dd3  
    1 """ 
    2 Class interface to the model calculator. 
    3  
    4 Calling a model is somewhat non-trivial since the functions called depend 
    5 on the data type.  For 1D data the *Iq* kernel needs to be called, for 
    6 2D data the *Iqxy* kernel needs to be called, and for SESANS data the 
    7 *Iq* kernel needs to be called followed by a Hankel transform.  Before 
    8 the kernel is called an appropriate *q* calculation vector needs to be 
    9 constructed.  This is not the simple *q* vector where you have measured 
    10 the data since the resolution calculation will require values beyond the 
    11 range of the measured data.  After the calculation the resolution calculator 
    12 must be called to return the predicted value for each measured data point. 
    13  
    14 :class:`DirectModel` is a callable object that takes *parameter=value* 
    15 keyword arguments and returns the appropriate theory values for the data. 
    16  
    17 :class:`DataMixin` does the real work of interpreting the data and calling 
    18 the model calculator.  This is used by :class:`DirectModel`, which uses 
    19 direct parameter values and by :class:`bumps_model.Experiment` which wraps 
    20 the parameter values in boxes so that the user can set fitting ranges, etc. 
    21 on the individual parameters and send the model to the Bumps optimizers. 
    22 """ 
    23 from __future__ import print_function 
     1import warnings 
    242 
    253import numpy as np 
    … …  
    4119    currently used by :class:`sasview_model.SasviewModel` since this will 
    4220    require a number of changes to SasView before we can do it. 
    43  
    44     :meth:`_interpret_data` initializes the data structures necessary 
    45     to manage the calculations.  This sets attributes in the child class 
    46     such as *data_type* and *resolution*. 
    47  
    48     :meth:`_calc_theory` evaluates the model at the given control values. 
    49  
    50     :meth:`_set_data` sets the intensity data in the data object, 
    51     possibly with random noise added.  This is useful for simulating a 
    52     dataset with the results from :meth:`_calc_theory`. 
    5321    """ 
    5422    def _interpret_data(self, data, model): 
    55         # pylint: disable=attribute-defined-outside-init 
    56  
    5723        self._data = data 
    5824        self._model = model 
    … …  
    7036        if self.data_type == 'sesans': 
    7137            q = sesans.make_q(data.sample.zacceptance, data.Rmax) 
    72             index = slice(None, None) 
    73             res = None 
     38            self.index = slice(None, None) 
    7439            if data.y is not None: 
    75                 Iq, dIq = data.y, data.dy 
    76             else: 
    77                 Iq, dIq = None, None 
     40                self.Iq = data.y 
     41                self.dIq = data.dy 
    7842            #self._theory = np.zeros_like(q) 
    7943            q_vectors = [q] 
    … …  
    8549            qmax = getattr(data, 'qmax', np.inf) 
    8650            accuracy = getattr(data, 'accuracy', 'Low') 
    87             index = ~data.mask & (q >= qmin) & (q <= qmax) 
     51            self.index = ~data.mask & (q >= qmin) & (q <= qmax) 
    8852            if data.data is not None: 
    89                 index &= ~np.isnan(data.data) 
    90                 Iq = data.data[index] 
    91                 dIq = data.err_data[index] 
     53                self.index &= ~np.isnan(data.data) 
     54                self.Iq = data.data[self.index] 
     55                self.dIq = data.err_data[self.index] 
     56            self.resolution = resolution2d.Pinhole2D(data=data, index=self.index, 
     57                                                     nsigma=3.0, accuracy=accuracy) 
     58            #self._theory = np.zeros_like(self.Iq) 
     59            q_vectors = self.resolution.q_calc 
     60        elif self.data_type == 'Iq': 
     61            self.index = (data.x >= data.qmin) & (data.x <= data.qmax) 
     62            if data.y is not None: 
     63                self.index &= ~np.isnan(data.y) 
     64                self.Iq = data.y[self.index] 
     65                self.dIq = data.dy[self.index] 
     66            if getattr(data, 'dx', None) is not None: 
     67                q, dq = data.x[self.index], data.dx[self.index] 
     68                if (dq>0).any(): 
     69                    self.resolution = resolution.Pinhole1D(q, dq) 
     70                else: 
     71                    self.resolution = resolution.Perfect1D(q) 
     72            elif (getattr(data, 'dxl', None) is not None and 
     73                          getattr(data, 'dxw', None) is not None): 
     74                self.resolution = resolution.Slit1D(data.x[self.index], 
     75                                                    width=data.dxh[self.index], 
     76                                                    height=data.dxw[self.index]) 
    9277            else: 
    93                 Iq, dIq = None, None 
    94             res = resolution2d.Pinhole2D(data=data, index=index, 
    95                                          nsigma=3.0, accuracy=accuracy) 
    96             #self._theory = np.zeros_like(self.Iq) 
    97             q_vectors = res.q_calc 
    98         elif self.data_type == 'Iq': 
    99             index = (data.x >= data.qmin) & (data.x <= data.qmax) 
    100             if data.y is not None: 
    101                 index &= ~np.isnan(data.y) 
    102                 Iq = data.y[index] 
    103                 dIq = data.dy[index] 
    104             else: 
    105                 Iq, dIq = None, None 
    106             if getattr(data, 'dx', None) is not None: 
    107                 q, dq = data.x[index], data.dx[index] 
    108                 if (dq > 0).any(): 
    109                     res = resolution.Pinhole1D(q, dq) 
    110                 else: 
    111                     res = resolution.Perfect1D(q) 
    112             elif (getattr(data, 'dxl', None) is not None 
    113                   and getattr(data, 'dxw', None) is not None): 
    114                 res = resolution.Slit1D(data.x[index], 
    115                                         width=data.dxh[index], 
    116                                         height=data.dxw[index]) 
    117             else: 
    118                 res = resolution.Perfect1D(data.x[index]) 
     78                self.resolution = resolution.Perfect1D(data.x[self.index]) 
    11979 
    12080            #self._theory = np.zeros_like(self.Iq) 
    121             q_vectors = [res.q_calc] 
     81            q_vectors = [self.resolution.q_calc] 
    12282        else: 
    12383            raise ValueError("Unknown data type") # never gets here 
    … …  
    12787        self._kernel_inputs = [v for v in q_vectors] 
    12888        self._kernel = None 
    129         self.Iq, self.dIq, self.index = Iq, dIq, index 
    130         self.resolution = res 
    13189 
    13290    def _set_data(self, Iq, noise=None): 
    133         # pylint: disable=attribute-defined-outside-init 
    13491        if noise is not None: 
    13592            self.dIq = Iq*noise*0.01 
    … …  
    149106    def _calc_theory(self, pars, cutoff=0.0): 
    150107        if self._kernel is None: 
    151             self._kernel = make_kernel(self._model, self._kernel_inputs)  # pylint: disable=attribute-defined-outside-init 
     108            q_input = self._model.make_input(self._kernel_inputs) 
     109            self._kernel = self._model(q_input) 
    152110 
    153111        Iq_calc = call_kernel(self._kernel, pars, cutoff=cutoff) 
    … …  
    162120 
    163121class DirectModel(DataMixin): 
    164     """ 
    165     Create a calculator object for a model. 
    166  
    167     *data* is 1D SAS, 2D SAS or SESANS data 
    168  
    169     *model* is a model calculator return from :func:`generate.load_model` 
    170  
    171     *cutoff* is the polydispersity weight cutoff. 
    172     """ 
    173122    def __init__(self, data, model, cutoff=1e-5): 
    174123        self.model = model 
    175124        self.cutoff = cutoff 
    176         # Note: _interpret_data defines the model attributes 
    177125        self._interpret_data(data, model) 
    178  
     126        self.kernel = make_kernel(self.model, self._kernel_inputs) 
    179127    def __call__(self, **pars): 
    180128        return self._calc_theory(pars, cutoff=self.cutoff) 
    181  
    182129    def ER(self, **pars): 
    183         """ 
    184         Compute the equivalent radius for the model. 
    185  
    186         Return 0. if not defined. 
    187         """ 
    188130        return call_ER(self.model.info, pars) 
    189  
    190131    def VR(self, **pars): 
    191         """ 
    192         Compute the equivalent volume for the model, including polydispersity 
    193         effects. 
    194  
    195         Return 1. if not defined. 
    196         """ 
    197132        return call_VR(self.model.info, pars) 
    198  
    199133    def simulate_data(self, noise=None, **pars): 
    200         """ 
    201         Generate simulated data for the model. 
    202         """ 
    203134        Iq = self.__call__(**pars) 
    204135        self._set_data(Iq, noise=noise) 
    205136 
    206 def main(): 
    207     """ 
    208     Program to evaluate a particular model at a set of q values. 
    209     """ 
     137def demo(): 
    210138    import sys 
    211139    from .data import empty_data1D, empty_data2D 
    … …  
    216144    model_name = sys.argv[1] 
    217145    call = sys.argv[2].upper() 
    218     if call not in ("ER", "VR"): 
     146    if call not in ("ER","VR"): 
    219147        try: 
    220148            values = [float(v) for v in call.split(',')] 
    … …  
    225153            data = empty_data1D([q]) 
    226154        elif len(values) == 2: 
    227             qx, qy = values 
    228             data = empty_data2D([qx], [qy]) 
     155            qx,qy = values 
     156            data = empty_data2D([qx],[qy]) 
    229157        else: 
    230158            print("use q or qx,qy or ER or VR") 
    … …  
    236164    model = load_model(model_definition, dtype='single') 
    237165    calculator = DirectModel(data, model) 
    238     pars = dict((k, float(v)) 
     166    pars = dict((k,float(v)) 
    239167                for pair in sys.argv[3:] 
    240                 for k, v in [pair.split('=')]) 
     168                for k,v in [pair.split('=')]) 
    241169    if call == "ER": 
    242170        print(calculator.ER(**pars)) 
    … …  
    248176 
    249177if __name__ == "__main__": 
    250     main() 
     178    demo() 

  • sasmodels/exception.py

    r823e620 reaca9eb  
    88except NameError: 
    99    class WindowsError(Exception): 
    10         """ 
    11         Fake WindowsException when not on Windows. 
    12         """ 
    1310        pass 
    1411 
    … …  
    3431    # TODO: try to incorporate current stack trace in the raised exception 
    3532    if isinstance(exc, WindowsError): 
    36         raise OSError(str(exc) + " " + msg) 
     33        raise OSError(str(exc)+" "+msg) 
    3734 
    3835    args = exc.args 
    … …  
    4138    else: 
    4239        try: 
    43             arg0 = " ".join((args[0], msg)) 
     40            arg0 = " ".join((args[0],msg)) 
    4441            exc.args = tuple([arg0] + list(args[1:])) 
    4542        except: 
    46             exc.args = (" ".join((str(exc), msg)),) 
     43            exc.args = (" ".join((str(exc),msg)),) 

  • sasmodels/generate.py

    reafc9fa r190fc2b  
    77    particular dimensions averaged over all orientations. 
    88 
    9     *Iqxy(qx, qy, p1, p2, ...)* returns the scattering at qx, qy for a form 
     9    *Iqxy(qx, qy, p1, p2, ...)* returns the scattering at qx,qy for a form 
    1010    with particular dimensions for a single orientation. 
    1111 
    … …  
    4747the *sin* and *cos* values are needed for a particular argument.  Since 
    4848this function does not exist in C99, all use of *sincos* should be 
    49 replaced by the macro *SINCOS(value, sn, cn)* where *sn* and *cn* are 
     49replaced by the macro *SINCOS(value,sn,cn)* where *sn* and *cn* are 
    5050previously declared *double* variables.  When compiled for systems without 
    5151OpenCL, *SINCOS* will be replaced by *sin* and *cos* calls.   If *value* is 
    … …  
    194194returns a list of files required by the model. 
    195195""" 
    196 from __future__ import print_function 
    197196 
    198197# TODO: identify model files which have changed since loading and reload them. 
    … …  
    313312    raise ValueError("%r not found in %s" % (filename, search_path)) 
    314313 
    315 def model_sources(info): 
     314def sources(info): 
    316315    """ 
    317316    Return a list of the sources file paths for the module. 
    … …  
    373372 
    374373 
    375 def kernel_name(info, is_2d): 
     374def kernel_name(info, is_2D): 
    376375    """ 
    377376    Name of the exported kernel symbol. 
    378377    """ 
    379     return info['name'] + "_" + ("Iqxy" if is_2d else "Iq") 
     378    return info['name'] + "_" + ("Iqxy" if is_2D else "Iq") 
    380379 
    381380 
    … …  
    420419 
    421420 
    422 LOOP_OPEN = """\ 
     421def build_polydispersity_loops(pd_pars): 
     422    """ 
     423    Build polydispersity loops 
     424 
     425    Returns loop opening and loop closing 
     426    """ 
     427    LOOP_OPEN = """\ 
    423428for (int %(name)s_i=0; %(name)s_i < N%(name)s; %(name)s_i++) { 
    424429  const double %(name)s = loops[2*(%(name)s_i%(offset)s)]; 
    425430  const double %(name)s_w = loops[2*(%(name)s_i%(offset)s)+1];\ 
    426431""" 
    427 def build_polydispersity_loops(pd_pars): 
    428     """ 
    429     Build polydispersity loops 
    430  
    431     Returns loop opening and loop closing 
    432     """ 
    433432    depth = 4 
    434433    offset = "" 
    … …  
    465464 
    466465    # Load additional sources 
    467     source = [open(f).read() for f in model_sources(info)] 
     466    source = [open(f).read() for f in sources(info)] 
    468467 
    469468    # Prepare defines 
    … …  
    573572 
    574573    #for d in defines: print(d) 
    575     defines = '\n'.join('#define %s %s' % (k, v) for k, v in defines) 
    576     sources = '\n\n'.join(source) 
     574    DEFINES = '\n'.join('#define %s %s' % (k, v) for k, v in defines) 
     575    SOURCES = '\n\n'.join(source) 
    577576    return C_KERNEL_TEMPLATE % { 
    578         'DEFINES': defines, 
    579         'SOURCES': sources, 
     577        'DEFINES':DEFINES, 
     578        'SOURCES':SOURCES, 
    580579        } 
    581580 
    … …  
    591590    demo_parameters = getattr(kernel_module, 'demo', None) 
    592591    if demo_parameters is None: 
    593         demo_parameters = dict((p[0], p[2]) for p in parameters) 
     592        demo_parameters = dict((p[0],p[2]) for p in parameters) 
    594593    filename = abspath(kernel_module.__file__) 
    595594    kernel_id = splitext(basename(filename))[0] 
    … …  
    598597        name = " ".join(w.capitalize() for w in kernel_id.split('_')) 
    599598    info = dict( 
    600         id=kernel_id,  # string used to load the kernel 
     599        id = kernel_id,  # string used to load the kernel 
    601600        filename=abspath(kernel_module.__file__), 
    602601        name=name, 
    … …  
    644643        elif section_marker.match(line): 
    645644            if len(line) >= len(prior): 
    646                 yield "".join(("**", prior, "**")) 
     645                yield "".join( ("**",prior,"**") ) 
    647646                prior = None 
    648647            else: 

  • sasmodels/kernelcl.py

    reafc9fa rcde11f0f  
    11""" 
    2 GPU driver for C kernels 
     2GPU support through OpenCL 
    33 
    44There should be a single GPU environment running on the system.  This 
    … …  
    152152 
    153153 
     154def make_result(self, size): 
     155    self.res = np.empty(size, dtype=self.dtype) 
     156    self.res_b = cl.Buffer(self.program.context, mf.READ_WRITE, self.res.nbytes) 
     157    return self.res, self.res_b 
     158 
     159 
    154160# for now, this returns one device in the context 
    155161# TODO: create a context that contains all devices on all platforms 
    … …  
    177183 
    178184    def has_type(self, dtype): 
    179         """ 
    180         Return True if all devices support a given type. 
    181         """ 
    182185        dtype = generate.F32 if dtype == 'fast' else np.dtype(dtype) 
    183186        return all(has_type(d, dtype) for d in self.context.devices) 
    184187 
    185188    def _create_some_context(self): 
    186         """ 
    187         Protected call to cl.create_some_context without interactivity.  Use 
    188         this if PYOPENCL_CTX is set in the environment.  Sets the *context* 
    189         attribute. 
    190         """ 
    191189        try: 
    192190            self.context = cl.create_some_context(interactive=False) 
    … …  
    197195 
    198196    def compile_program(self, name, source, dtype, fast=False): 
    199         """ 
    200         Compile the program for the device in the given context. 
    201         """ 
    202197        key = "%s-%s-%s"%(name, dtype, fast) 
    203198        if key not in self.compiled: 
    … …  
    209204 
    210205    def release_program(self, name): 
    211         """ 
    212         Free memory associated with the program on the device. 
    213         """ 
    214206        if name in self.compiled: 
    215207            self.compiled[name].release() 
    … …  
    217209 
    218210def _get_default_context(): 
    219     """ 
    220     Get an OpenCL context, preferring GPU over CPU. 
    221     """ 
    222211    default = None 
    223212    for platform in cl.get_platforms(): 
    … …  
    252241        self.info = info 
    253242        self.source = source 
    254         self.dtype = generate.F32 if dtype == 'fast' else np.dtype(dtype) 
     243        self.dtype = generate.F32 if dtype=='fast' else np.dtype(dtype) 
    255244        self.fast = (dtype == 'fast') 
    256245        self.program = None # delay program creation 
    257246 
    258247    def __getstate__(self): 
    259         return self.info, self.source, self.dtype, self.fast 
     248        state = self.__dict__.copy() 
     249        state['program'] = None 
     250        return state 
    260251 
    261252    def __setstate__(self, state): 
    262         self.info, self.source, self.dtype, self.fast = state 
    263         self.program = None 
    264  
    265     def __call__(self, q_vectors): 
     253        self.__dict__ = state.copy() 
     254 
     255    def __call__(self, q_input): 
     256        if self.dtype != q_input.dtype: 
     257            raise TypeError("data is %s kernel is %s" 
     258                            % (q_input.dtype, self.dtype)) 
    266259        if self.program is None: 
    267260            compiler = environment().compile_program 
    268261            self.program = compiler(self.info['name'], self.source, self.dtype, 
    269262                                    self.fast) 
    270         is_2d = len(q_vectors) == 2 
    271         kernel_name = generate.kernel_name(self.info, is_2d) 
     263        kernel_name = generate.kernel_name(self.info, q_input.is_2D) 
    272264        kernel = getattr(self.program, kernel_name) 
    273         return GpuKernel(kernel, self.info, q_vectors, self.dtype) 
     265        return GpuKernel(kernel, self.info, q_input) 
    274266 
    275267    def release(self): 
    276         """ 
    277         Free the resources associated with the model. 
    278         """ 
    279268        if self.program is not None: 
    280269            environment().release_program(self.info['name']) 
    281270            self.program = None 
    282271 
    283     def __del__(self): 
    284         self.release() 
     272    def make_input(self, q_vectors): 
     273        """ 
     274        Make q input vectors available to the model. 
     275 
     276        Note that each model needs its own q vector even if the case of 
     277        mixture models because some models may be OpenCL, some may be 
     278        ctypes and some may be pure python. 
     279        """ 
     280        return GpuInput(q_vectors, dtype=self.dtype) 
    285281 
    286282# TODO: check that we don't need a destructor for buffers which go out of scope 
    … …  
    308304        self.nq = q_vectors[0].size 
    309305        self.dtype = np.dtype(dtype) 
    310         self.is_2d = (len(q_vectors) == 2) 
     306        self.is_2D = (len(q_vectors) == 2) 
    311307        # TODO: stretch input based on get_warp() 
    312308        # not doing it now since warp depends on kernel, which is not known 
    … …  
    321317 
    322318    def release(self): 
    323         """ 
    324         Free the memory. 
    325         """ 
    326319        for b in self.q_buffers: 
    327320            b.release() 
    328321        self.q_buffers = [] 
    329322 
    330     def __del__(self): 
    331         self.release() 
    332  
    333323class GpuKernel(object): 
    334324    """ 
    335325    Callable SAS kernel. 
    336326 
    337     *kernel* is the GpuKernel object to call 
     327    *kernel* is the GpuKernel object to call. 
    338328 
    339329    *info* is the module information 
    340330 
    341     *q_vectors* is the q vectors at which the kernel should be evaluated 
    342  
    343     *dtype* is the kernel precision 
     331    *q_input* is the DllInput q vectors at which the kernel should be 
     332    evaluated. 
    344333 
    345334    The resulting call method takes the *pars*, a list of values for 
    … …  
    351340    Call :meth:`release` when done with the kernel instance. 
    352341    """ 
    353     def __init__(self, kernel, info, q_vectors, dtype): 
    354         q_input = GpuInput(q_vectors, dtype) 
     342    def __init__(self, kernel, info, q_input): 
     343        self.q_input = q_input 
    355344        self.kernel = kernel 
    356345        self.info = info 
    357346        self.res = np.empty(q_input.nq, q_input.dtype) 
    358         dim = '2d' if q_input.is_2d else '1d' 
     347        dim = '2d' if q_input.is_2D else '1d' 
    359348        self.fixed_pars = info['partype']['fixed-' + dim] 
    360349        self.pd_pars = info['partype']['pd-' + dim] 
    … …  
    369358                                q_input.global_size[0] * q_input.dtype.itemsize) 
    370359                      for _ in env.queues] 
    371         self.q_input = q_input 
     360 
    372361 
    373362    def __call__(self, fixed_pars, pd_pars, cutoff=1e-5): 
    … …  
    410399 
    411400    def release(self): 
    412         """ 
    413         Release resources associated with the kernel. 
    414         """ 
    415401        for b in self.loops_b: 
    416402            b.release() 
    … …  
    419405            b.release() 
    420406        self.res_b = [] 
    421         self.q_input.release() 
    422407 
    423408    def __del__(self): 

  • sasmodels/kerneldll.py

    reafc9fa r74667d3  
    11r""" 
    2 DLL driver for C kernels 
     2C types wrapper for sasview models. 
    33 
    44The global attribute *ALLOW_SINGLE_PRECISION_DLLS* should be set to *True* if 
    … …  
    4444directory for this application, then OpenMP should be supported. 
    4545""" 
    46 from __future__ import print_function 
    4746 
    4847import sys 
    … …  
    123122    models are allowed as DLLs. 
    124123    """ 
    125     if callable(info.get('Iq', None)): 
     124    if callable(info.get('Iq',None)): 
    126125        return PyModel(info) 
    127126 
    … …  
    140139 
    141140    source = generate.convert_type(source, dtype) 
    142     source_files = generate.model_sources(info) + [info['filename']] 
    143     dll = dll_path(info, dtype) 
     141    source_files = generate.sources(info) + [info['filename']] 
     142    dll= dll_path(info, dtype) 
    144143    newest = max(os.path.getmtime(f) for f in source_files) 
    145     if not os.path.exists(dll) or os.path.getmtime(dll) < newest: 
     144    if not os.path.exists(dll) or os.path.getmtime(dll)<newest: 
    146145        # Replace with a proper temp file 
    147         fid, filename = tempfile.mkstemp(suffix=".c", prefix=tempfile_prefix) 
    148         os.fdopen(fid, "w").write(source) 
     146        fid, filename = tempfile.mkstemp(suffix=".c",prefix=tempfile_prefix) 
     147        os.fdopen(fid,"w").write(source) 
    149148        command = COMPILE%{"source":filename, "output":dll} 
    150149        print("Compile command: "+command) 
    … …  
    161160def load_dll(source, info, dtype="double"): 
    162161    """ 
    163     Create and load a dll corresponding to the source, info pair returned 
     162    Create and load a dll corresponding to the source,info pair returned 
    164163    from :func:`sasmodels.generate.make` compiled for the target precision. 
    165164 
    … …  
    200199        Npd2d = len(self.info['partype']['pd-2d']) 
    201200 
    202         #print("dll", self.dllpath) 
     201        #print("dll",self.dllpath) 
    203202        try: 
    204203            self.dll = ct.CDLL(self.dllpath) 
    … …  
    211210              else c_longdouble) 
    212211        pd_args_1d = [c_void_p, fp] + [c_int]*Npd1d if Npd1d else [] 
    213         pd_args_2d = [c_void_p, fp] + [c_int]*Npd2d if Npd2d else [] 
     212        pd_args_2d= [c_void_p, fp] + [c_int]*Npd2d if Npd2d else [] 
    214213        self.Iq = self.dll[generate.kernel_name(self.info, False)] 
    215214        self.Iq.argtypes = IQ_ARGS + pd_args_1d + [fp]*Nfixed1d 
    … …  
    219218 
    220219    def __getstate__(self): 
    221         return self.info, self.dllpath 
     220        return {'info': self.info, 'dllpath': self.dllpath, 'dll': None} 
    222221 
    223222    def __setstate__(self, state): 
    224         self.info, self.dllpath = state 
    225         self.dll = None 
    226  
    227     def __call__(self, q_vectors): 
    228         q_input = PyInput(q_vectors, self.dtype) 
     223        self.__dict__ = state 
     224 
     225    def __call__(self, q_input): 
     226        if self.dtype != q_input.dtype: 
     227            raise TypeError("data is %s kernel is %s" % (q_input.dtype, self.dtype)) 
    229228        if self.dll is None: self._load_dll() 
    230         kernel = self.Iqxy if q_input.is_2d else self.Iq 
     229        kernel = self.Iqxy if q_input.is_2D else self.Iq 
    231230        return DllKernel(kernel, self.info, q_input) 
    232231 
     232    # pylint: disable=no-self-use 
     233    def make_input(self, q_vectors): 
     234        """ 
     235        Make q input vectors available to the model. 
     236 
     237        Note that each model needs its own q vector even if the case of 
     238        mixture models because some models may be OpenCL, some may be 
     239        ctypes and some may be pure python. 
     240        """ 
     241        return PyInput(q_vectors, dtype=self.dtype) 
     242 
    233243    def release(self): 
    234         """ 
    235         Release any resources associated with the model. 
    236         """ 
    237244        pass # TODO: should release the dll 
    238245 
    … …  
    250257 
    251258    The resulting call method takes the *pars*, a list of values for 
    252     the fixed parameters to the kernel, and *pd_pars*, a list of (value, weight) 
     259    the fixed parameters to the kernel, and *pd_pars*, a list of (value,weight) 
    253260    vectors for the polydisperse parameters.  *cutoff* determines the 
    254261    integration limits: any points with combined weight less than *cutoff* 
    … …  
    262269        self.kernel = kernel 
    263270        self.res = np.empty(q_input.nq, q_input.dtype) 
    264         dim = '2d' if q_input.is_2d else '1d' 
     271        dim = '2d' if q_input.is_2D else '1d' 
    265272        self.fixed_pars = info['partype']['fixed-'+dim] 
    266273        self.pd_pars = info['partype']['pd-'+dim] 
    … …  
    292299 
    293300    def release(self): 
    294         """ 
    295         Release any resources associated with the kernel. 
    296         """ 
    297301        pass 

  • sasmodels/kernelpy.py

    reafc9fa r4c2c535  
    1 """ 
    2 Python driver for python kernels 
    3  
    4 Calls the kernel with a vector of $q$ values for a single parameter set. 
    5 Polydispersity is supported by looping over different parameter sets and 
    6 summing the results.  The interface to :class:`PyModel` matches those for 
    7 :class:`kernelcl.GpuModel` and :class:`kerneldll.DllModel`. 
    8 """ 
    91import numpy as np 
    102from numpy import pi, cos 
    … …  
    135 
    146class PyModel(object): 
    15     """ 
    16     Wrapper for pure python models. 
    17     """ 
    187    def __init__(self, info): 
    198        self.info = info 
    209 
    21     def __call__(self, q_vectors): 
    22         q_input = PyInput(q_vectors, dtype=F64) 
    23         kernel = self.info['Iqxy'] if q_input.is_2d else self.info['Iq'] 
     10    def __call__(self, q_input): 
     11        kernel = self.info['Iqxy'] if q_input.is_2D else self.info['Iq'] 
    2412        return PyKernel(kernel, self.info, q_input) 
    2513 
     14    # pylint: disable=no-self-use 
     15    def make_input(self, q_vectors): 
     16        return PyInput(q_vectors, dtype=F64) 
     17 
    2618    def release(self): 
    27         """ 
    28         Free resources associated with the model. 
    29         """ 
    3019        pass 
    3120 
    … …  
    5241        self.nq = q_vectors[0].size 
    5342        self.dtype = dtype 
    54         self.is_2d = (len(q_vectors) == 2) 
     43        self.is_2D = (len(q_vectors) == 2) 
    5544        self.q_vectors = [np.ascontiguousarray(q, self.dtype) for q in q_vectors] 
    5645        self.q_pointers = [q.ctypes.data for q in self.q_vectors] 
    5746 
    5847    def release(self): 
    59         """ 
    60         Free resources associated with the model inputs. 
    61         """ 
    6248        self.q_vectors = [] 
    6349 
    … …  
    8571        self.q_input = q_input 
    8672        self.res = np.empty(q_input.nq, q_input.dtype) 
    87         dim = '2d' if q_input.is_2d else '1d' 
     73        dim = '2d' if q_input.is_2D else '1d' 
    8874        # Loop over q unless user promises that the kernel is vectorized by 
    8975        # taggining it with vectorized=True 
    … …  
    9177            if dim == '2d': 
    9278                def vector_kernel(qx, qy, *args): 
    93                     """ 
    94                     Vectorized 2D kernel. 
    95                     """ 
    9679                    return np.array([kernel(qxi, qyi, *args) 
    9780                                     for qxi, qyi in zip(qx, qy)]) 
    9881            else: 
    9982                def vector_kernel(q, *args): 
    100                     """ 
    101                     Vectorized 1D kernel. 
    102                     """ 
    10383                    return np.array([kernel(qi, *args) 
    10484                                     for qi in q]) 
    … …  
    142122 
    143123    def release(self): 
    144         """ 
    145         Free resources associated with the kernel. 
    146         """ 
    147124        self.q_input = None 
    148125 

  • sasmodels/model_test.py

    r5c962df r9404dd3  
    99    if model1 is 'all', then all except the remaining models will be tested 
    1010 
    11 Each model is tested using the default parameters at q=0.1, (qx, qy)=(0.1, 0.1), 
     11Each model is tested using the default parameters at q=0.1, (qx,qy)=(0.1,0.1), 
    1212and the ER and VR are computed.  The return values at these points are not 
    1313considered.  The test is only to verify that the models run to completion, 
    … …  
    3030        [ {parameters}, (qx, qy), I(qx, Iqy)], 
    3131        [ {parameters}, [(qx1, qy1), (qx2, qy2), ...], 
    32                         [I(qx1, qy1), I(qx2, qy2), ...]], 
     32                        [I(qx1,qy1), I(qx2,qy2), ...]], 
    3333 
    3434        [ {parameters}, 'ER', ER(pars) ], 
    … …  
    4343Precision defaults to 5 digits (relative). 
    4444""" 
    45 from __future__ import print_function 
    4645 
    4746import sys 
    … …  
    5655 
    5756def make_suite(loaders, models): 
    58     """ 
    59     Construct the pyunit test suite. 
    60  
    61     *loaders* is the list of kernel drivers to use, which is one of 
    62     *["dll", "opencl"]*, *["dll"]* or *["opencl"]*.  For python models, 
    63     the python driver is always used. 
    64  
    65     *models* is the list of models to test, or *["all"]* to test all models. 
    66     """ 
    6757 
    6858    ModelTestCase = _hide_model_case_from_nosetests() 
    … …  
    8575        # don't try to call cl kernel since it will not be 
    8676        # available in some environmentes. 
    87         is_py = callable(getattr(model_definition, 'Iq', None)) 
     77        is_py = callable(getattr(model_definition,'Iq', None)) 
    8878 
    8979        if is_py:  # kernel implemented in python 
    … …  
    121111def _hide_model_case_from_nosetests(): 
    122112    class ModelTestCase(unittest.TestCase): 
    123         """ 
    124         Test suit for a particular model with a particular kernel driver. 
    125  
    126         The test suite runs a simple smoke test to make sure the model 
    127         functions, then runs the list of tests at the bottom of the model 
    128         description file. 
    129         """ 
    130113        def __init__(self, test_name, definition, test_method_name, 
    131114                     platform, dtype): 
    … …  
    140123        def _runTest(self): 
    141124            smoke_tests = [ 
    142                 [{}, 0.1, None], 
    143                 [{}, (0.1, 0.1), None], 
    144                 [{}, 'ER', None], 
    145                 [{}, 'VR', None], 
     125                [{},0.1,None], 
     126                [{},(0.1,0.1),None], 
     127                [{},'ER',None], 
     128                [{},'VR',None], 
    146129                ] 
    147130 
    … …  
    180163                actual = [call_VR(model.info, pars)] 
    181164            elif isinstance(x[0], tuple): 
    182                 Qx, Qy = zip(*x) 
     165                Qx,Qy = zip(*x) 
    183166                q_vectors = [np.array(Qx), np.array(Qy)] 
    184167                kernel = make_kernel(model, q_vectors) 
    … …  
    196179                    # smoke test --- make sure it runs and produces a value 
    197180                    self.assertTrue(np.isfinite(actual_yi), 
    198                                     'invalid f(%s): %s' % (xi, actual_yi)) 
     181                        'invalid f(%s): %s' % (xi, actual_yi)) 
    199182                else: 
    200183                    err = abs(yi - actual_yi) 
    201184                    nrm = abs(yi) 
    202185                    self.assertLess(err * 10**5, nrm, 
    203                                     'f(%s); expected:%s; actual:%s' 
    204                                     % (xi, yi, actual_yi)) 
     186                        'f(%s); expected:%s; actual:%s' % (xi, yi, actual_yi)) 
    205187 
    206188    return ModelTestCase 
    … …  
    255237    Run "nosetests sasmodels" on the command line to invoke it. 
    256238    """ 
    257     tests = make_suite(['opencl', 'dll'], ['all']) 
     239    tests = make_suite(['opencl','dll'],['all']) 
    258240    for test_i in tests: 
    259241        yield test_i._runTest 

  • sasmodels/resolution.py

    r5925e90 r190fc2b  
    1414           "pinhole_extend_q", "slit_extend_q", "bin_edges", 
    1515           "interpolate", "linear_extrapolation", "geometric_extrapolation", 
    16           ] 
     16           ] 
    1717 
    1818MINIMUM_RESOLUTION = 1e-8 
    … …  
    8080        # default to Perfect1D if the pinhole geometry is not defined. 
    8181        self.q, self.q_width = q, q_width 
    82         self.q_calc = (pinhole_extend_q(q, q_width, nsigma=nsigma) 
    83                        if q_calc is None else np.sort(q_calc)) 
    84         self.weight_matrix = pinhole_resolution( 
    85             self.q_calc, self.q, np.maximum(q_width, MINIMUM_RESOLUTION)) 
     82        self.q_calc = pinhole_extend_q(q, q_width, nsigma=nsigma) \ 
     83            if q_calc is None else np.sort(q_calc) 
     84        self.weight_matrix = pinhole_resolution(self.q_calc, 
     85                self.q, np.maximum(q_width, MINIMUM_RESOLUTION)) 
    8686 
    8787    def apply(self, theory): 
    … …  
    470470 
    471471def eval_form(q, form, pars): 
    472     """ 
    473     Return the SAS model evaluated at *q*. 
    474  
    475     *form* is the SAS model returned from :fun:`core.load_model`. 
    476  
    477     *pars* are the parameter values to use when evaluating. 
    478     """ 
    479472    from sasmodels import core 
    480473    kernel = core.make_kernel(form, [q]) 
    … …  
    485478 
    486479def gaussian(q, q0, dq): 
    487     """ 
    488     Return the Gaussian resolution function. 
    489  
    490     *q0* is the center, *dq* is the width and *q* are the points to evaluate. 
    491     """ 
    492480    from numpy import exp, pi 
    493481    return exp(-0.5*((q-q0)/dq)**2)/(sqrt(2*pi)*dq) 
    … …  
    571559 
    572560class ResolutionTest(unittest.TestCase): 
    573     """ 
    574     Test the resolution calculations. 
    575     """ 
    576561 
    577562    def setUp(self): 
    … …  
    683668 
    684669class IgorComparisonTest(unittest.TestCase): 
    685     """ 
    686     Test resolution calculations against those returned by Igor. 
    687     """ 
    688670 
    689671    def setUp(self): 
    … …  
    693675        self.model = core.load_model(sphere, dtype='double') 
    694676 
    695     def _eval_sphere(self, pars, resolution): 
     677    def Iq_sphere(self, pars, resolution): 
    696678        from sasmodels import core 
    697679        kernel = core.make_kernel(self.model, [resolution.q_calc]) 
    … …  
    701683        return result 
    702684 
    703     def _compare(self, q, output, answer, tolerance): 
     685    def compare(self, q, output, answer, tolerance): 
    704686        #err = (output - answer)/answer 
    705687        #idx = abs(err) >= tolerance 
    … …  
    718700        q, width, answer, _ = data 
    719701        resolution = Perfect1D(q) 
    720         output = self._eval_sphere(pars, resolution) 
    721         self._compare(q, output, answer, 1e-6) 
     702        output = self.Iq_sphere(pars, resolution) 
     703        self.compare(q, output, answer, 1e-6) 
    722704 
    723705    def test_pinhole(self): 
    … …  
    731713        q, q_width, answer = data 
    732714        resolution = Pinhole1D(q, q_width) 
    733         output = self._eval_sphere(pars, resolution) 
     715        output = self.Iq_sphere(pars, resolution) 
    734716        # TODO: relative error should be lower 
    735         self._compare(q, output, answer, 3e-4) 
     717        self.compare(q, output, answer, 3e-4) 
    736718 
    737719    def test_pinhole_romberg(self): 
    … …  
    754736        q_calc, tol = None, 0.01 
    755737        resolution = Pinhole1D(q, q_width, q_calc=q_calc) 
    756         output = self._eval_sphere(pars, resolution) 
     738        output = self.Iq_sphere(pars, resolution) 
    757739        if 0: # debug plot 
    758740            import matplotlib.pyplot as plt 
    759741            resolution = Perfect1D(q) 
    760             source = self._eval_sphere(pars, resolution) 
     742            source = self.Iq_sphere(pars, resolution) 
    761743            plt.loglog(q, source, '.') 
    762744            plt.loglog(q, answer, '-', hold=True) 
    763745            plt.loglog(q, output, '-', hold=True) 
    764746            plt.show() 
    765         self._compare(q, output, answer, tol) 
     747        self.compare(q, output, answer, tol) 
    766748 
    767749    def test_slit(self): 
    … …  
    775757        q, delta_qv, _, answer = data 
    776758        resolution = Slit1D(q, width=delta_qv, height=0) 
    777         output = self._eval_sphere(pars, resolution) 
     759        output = self.Iq_sphere(pars, resolution) 
    778760        # TODO: eliminate Igor test since it is too inaccurate to be useful. 
    779761        # This means we can eliminate the test data as well, and instead 
    780762        # use a generated q vector. 
    781         self._compare(q, output, answer, 0.5) 
     763        self.compare(q, output, answer, 0.5) 
    782764 
    783765    def test_slit_romberg(self): 
    … …  
    795777                               delta_qv[0], 0.) 
    796778        resolution = Slit1D(q, width=delta_qv, height=0, q_calc=q_calc) 
    797         output = self._eval_sphere(pars, resolution) 
     779        output = self.Iq_sphere(pars, resolution) 
    798780        # TODO: relative error should be lower 
    799         self._compare(q, output, answer, 0.025) 
     781        self.compare(q, output, answer, 0.025) 
    800782 
    801783    def test_ellipsoid(self): 
    … …  
    817799        # TODO: 10% is too much error; use better algorithm 
    818800        #print(np.max(abs(answer-output)/answer)) 
    819         self._compare(q, output, answer, 0.1) 
     801        self.compare(q, output, answer, 0.1) 
    820802 
    821803    #TODO: can sas q spacing be too sparse for the resolution calculation? 
    … …  
    831813        q, q_width, answer = data[:, ::20] # Take every nth point 
    832814        resolution = Pinhole1D(q, q_width) 
    833         output = self._eval_sphere(pars, resolution) 
    834         self._compare(q, output, answer, 1e-6) 
     815        output = self.Iq_sphere(pars, resolution) 
     816        self.compare(q, output, answer, 1e-6) 
    835817 
    836818 
    … …  
    10851067 
    10861068def demo_pinhole_1d(): 
    1087     """ 
    1088     Show example of pinhole smearing. 
    1089     """ 
    10901069    q = np.logspace(-4, np.log10(0.2), 400) 
    10911070    q_width = 0.1*q 
    … …  
    10941073 
    10951074def demo_slit_1d(): 
    1096     """ 
    1097     Show example of slit smearing. 
    1098     """ 
    10991075    q = np.logspace(-4, np.log10(0.2), 100) 
    11001076    w = h = 0. 
    … …  
    11071083 
    11081084def demo(): 
    1109     """ 
    1110     Run the resolution demos. 
    1111     """ 
    11121085    import matplotlib.pyplot as plt 
    11131086    plt.subplot(121) 

  • sasmodels/resolution2d.py

    r823e620 r841753c  
    143143        if self.coords == 'polar': 
    144144            q_r = sqrt(qx**2 + qy**2) 
    145             qx_res = ((dqx*cos(dphi) + q_r) * cos(-q_phi) 
    146                       + dqy*sin(dphi) * sin(-q_phi)) 
    147             qy_res = (-(dqx*cos(dphi) + q_r) * sin(-q_phi) 
    148                       + dqy*sin(dphi) * cos(-q_phi)) 
     145            qx_res = ((dqx*cos(dphi) + q_r) * cos(-q_phi) + 
     146                          dqy*sin(dphi) * sin(-q_phi)) 
     147            qy_res = (-(dqx*cos(dphi) + q_r) * sin(-q_phi) + 
     148                          dqy*sin(dphi) * cos(-q_phi)) 
    149149        else: 
    150             qx_res = qx + dqx*cos(dphi) 
    151             qy_res = qy + dqy*sin(dphi) 
     150            qx_res = qx +  dqx*cos(dphi) 
     151            qy_res = qy +  dqy*sin(dphi) 
    152152 
    153153 
    … …  
    162162            theory = np.reshape(theory, (nbins, nq)) 
    163163            ## Averaging with Gaussian weighting: normalization included. 
    164             value = np.average(theory, axis=0, weights=self.q_calc_weights) 
     164            value =np.average(theory, axis=0, weights=self.q_calc_weights) 
    165165            ## Return the smeared values in the range of self.index 
    166166            return value 

  • sasmodels/sasview_model.py

    reafc9fa r9404dd3  
    283283        """ 
    284284        q_vectors = [np.asarray(q) for q in args] 
    285         fn = self._model(q_vectors) 
     285        fn = self._model(self._model.make_input(q_vectors)) 
    286286        pars = [self.params[v] for v in fn.fixed_pars] 
    287287        pd_pars = [self._get_weights(p) for p in fn.pd_pars] 

  • sasmodels/weights.py

    r5c962df r3c56da87  
    2222 
    2323    def get_pars(self): 
    24         """ 
    25         Return the parameters to the disperser as a dictionary. 
    26         """ 
    2724        pars = {'type': self.type} 
    2825        pars.update(self.__dict__) 
    … …  
    3128    # pylint: disable=no-self-use 
    3229    def set_weights(self, values, weights): 
    33         """ 
    34         Set the weights on the disperser if it is :class:`ArrayDispersion`. 
    35         """ 
    3630        raise RuntimeError("set_weights is only available for ArrayDispersion") 
    3731 
    … …  
    4236        *center* is the center of the distribution 
    4337 
    44         *lb*, *ub* are the min and max allowed values 
     38        *lb*,*ub* are the min and max allowed values 
    4539 
    4640        *relative* is True if the distribution width is proportional to the 
    … …  
    6963 
    7064class GaussianDispersion(Dispersion): 
    71     r""" 
    72     Gaussian dispersion, with 1-\ $\sigma$ width. 
    73  
    74     .. math:: 
    75  
    76         w = \exp\left(-\tfrac12 (x - c)^2/\sigma^2\right) 
    77     """ 
    7865    type = "gaussian" 
    7966    default = dict(npts=35, width=0, nsigmas=3) 
    … …  
    8572 
    8673class RectangleDispersion(Dispersion): 
    87     r""" 
    88     Uniform dispersion, with width $\sqrt{3}\sigma$. 
    89  
    90     .. math:: 
    91  
    92         w = 1 
    93     """ 
    9474    type = "rectangle" 
    9575    default = dict(npts=35, width=0, nsigmas=1.70325) 
    … …  
    10181 
    10282class LogNormalDispersion(Dispersion): 
    103     r""" 
    104     log Gaussian dispersion, with 1-\ $\sigma$ width. 
    105  
    106     .. math:: 
    107  
    108         w = \frac{\exp\left(-\tfrac12 (\ln x - c)^2/\sigma^2\right)}{x\sigma} 
    109     """ 
    11083    type = "lognormal" 
    11184    default = dict(npts=80, width=0, nsigmas=8) 
    11285    def _weights(self, center, sigma, lb, ub): 
    113         x = self._linspace(center, sigma, max(lb, 1e-8), max(ub, 1e-8)) 
    114         px = np.exp(-0.5*(np.log(x)-center)**2/sigma**2)/(x*sigma) 
     86        x = self._linspace(center, sigma, max(lb,1e-8), max(ub,1e-8)) 
     87        px = np.exp(-0.5*(np.log(x)-center)**2)/sigma**2/(x*sigma) 
    11588        return x, px 
    11689 
    11790 
    11891class SchulzDispersion(Dispersion): 
    119     r""" 
    120     Schultz dispersion, with 1-\ $\sigma$ width. 
    121  
    122     .. math:: 
    123  
    124         w = \frac{z^z\,R^{z-1}}{e^{Rz}\,c \Gamma(z)} 
    125  
    126     where $c$ is the center of the distribution, $R = x/c$ and $z=(c/\sigma)^2$. 
    127  
    128     This is evaluated using logarithms as 
    129  
    130     .. math:: 
    131  
    132         w = \exp\left(z \ln z + (z-1)\ln R - Rz - \ln c - \ln \Gamma(z) \right) 
    133     """ 
    13492    type = "schulz" 
    13593    default = dict(npts=80, width=0, nsigmas=8) 
    13694    def _weights(self, center, sigma, lb, ub): 
    137         x = self._linspace(center, sigma, max(lb, 1e-8), max(ub, 1e-8)) 
    138         R = x/center 
     95        x = self._linspace(center, sigma, max(lb,1e-8), max(ub,1e-8)) 
     96        R= x/center 
    13997        z = (center/sigma)**2 
    14098        arg = z*np.log(z) + (z-1)*np.log(R) - R*z - np.log(center) - gammaln(z) 
    … …  
    144102 
    145103class ArrayDispersion(Dispersion): 
    146     r""" 
    147     Empirical dispersion curve. 
    148  
    149     Use :meth:`set_weights` to set $w = f(x)$. 
    150     """ 
    151104    type = "array" 
    152105    default = dict(npts=35, width=0, nsigmas=1) 
    … …  
    157110 
    158111    def set_weights(self, values, weights): 
    159         """ 
    160         Set the weights for the given x values. 
    161         """ 
    162112        self.values = np.ascontiguousarray(values, 'd') 
    163113        self.weights = np.ascontiguousarray(weights, 'd') 
    … …  
    167117        # TODO: interpolate into the array dispersion using npts 
    168118        x = center + self.values*sigma 
    169         idx = (x >= lb) & (x <= ub) 
     119        idx = (x>=lb)&(x<=ub) 
    170120        x = x[idx] 
    171121        px = self.weights[idx] 
    … …  
    174124 
    175125# dispersion name -> disperser lookup table. 
    176 models = dict((d.type, d) for d in ( 
     126models = dict((d.type,d) for d in ( 
    177127    GaussianDispersion, RectangleDispersion, 
    178128    ArrayDispersion, SchulzDispersion, LogNormalDispersion 
    … …  
    199149    of the parameter, and false if it is an absolute width. 
    200150 
    201     Returns *(value, weight)*, where *value* and *weight* are vectors. 
     151    Returns *(value,weight)*, where *value* and *weight* are vectors. 
    202152    """ 
    203153    cls = models[disperser] 
    204154    obj = cls(n, width, nsigmas) 
    205     v, w = obj.get_weights(value, limits[0], limits[1], relative) 
    206     return v, w 
     155    v,w = obj.get_weights(value, limits[0], limits[1], relative) 
     156    return v,w 
    207157 
    208158# Hack to allow sasview dispersion objects to interoperate with sasmodels 
    209 dispersers = dict((v.__name__, k) for k, v in models.items()) 
     159dispersers = dict((v.__name__,k) for k,v in models.items()) 
    210160dispersers['DispersionModel'] = RectangleDispersion.type 
    211161 
Note: See TracChangeset for help on using the changeset viewer.