Changes in / [fb5914f:7a1143b] in sasmodels


Ignore:
Files:
1 added
14 edited

Legend:

Unmodified
Added
Removed

  • example/fit.py

    r7cf2cfd r1182da5  
    3131    model = Model(kernel, 
    3232        scale=0.08, 
    33         rpolar=15, requatorial=800, 
    34         sld=.291, solvent_sld=7.105, 
     33        r_polar=15, r_equatorial=800, 
     34        sld=.291, sld_solvent=7.105, 
    3535        background=0, 
    3636        theta=90, phi=phi, 
    3737        theta_pd=15, theta_pd_n=40, theta_pd_nsigma=3, 
    38         rpolar_pd=0.222296, rpolar_pd_n=1, rpolar_pd_nsigma=0, 
    39         requatorial_pd=.000128, requatorial_pd_n=1, requatorial_pd_nsigma=0, 
     38        r_polar_pd=0.222296, r_polar_pd_n=1, r_polar_pd_nsigma=0, 
     39        r_equatorial_pd=.000128, r_equatorial_pd_n=1, r_equatorial_pd_nsigma=0, 
    4040        phi_pd=0, phi_pd_n=20, phi_pd_nsigma=3, 
    4141        ) 
    … …  
    4343 
    4444    # SET THE FITTING PARAMETERS 
    45     model.rpolar.range(15, 1000) 
    46     model.requatorial.range(15, 1000) 
     45    model.r_polar.range(15, 1000) 
     46    model.r_equatorial.range(15, 1000) 
    4747    model.theta_pd.range(0, 360) 
    4848    model.background.range(0,1000) 
    … …  
    8181    pars = dict( 
    8282        scale=.01, background=35, 
    83         sld=.291, solvent_sld=5.77,  
     83        sld=.291, sld_solvent=5.77, 
    8484        radius=250, length=178,  
    8585        theta=90, phi=phi, 
    … …  
    105105    model = Model(kernel, 
    106106        scale= .031, radius=19.5, thickness=30, length=22, 
    107         core_sld=7.105, shell_sld=.291, solvent_sld=7.105, 
     107        sld_core=7.105, sld_shell=.291, sdl_solvent=7.105, 
    108108        background=0, theta=0, phi=phi, 
    109109 
    … …  
    132132    model = Model(kernel, 
    133133        scale=.08, radius=20, cap_radius=40, length=400, 
    134         sld_capcyl=1, sld_solv=6.3, 
     134        sld=1, sld_solvent=6.3, 
    135135        background=0, theta=0, phi=phi, 
    136136        radius_pd=.1, radius_pd_n=5, radius_pd_nsigma=3, 
    … …  
    147147    model = Model(kernel, 
    148148        scale=0.08, req_minor=15, req_major=20, rpolar=500, 
    149         sldEll=7.105, solvent_sld=.291, 
     149        sld=7.105, solvent_sld=.291, 
    150150        background=5, theta=0, phi=phi, psi=0, 
    151151        theta_pd=20, theta_pd_n=40, theta_pd_nsigma=3, 

  • example/model.py

    r15d2285 r1182da5  
    1717model = Model(kernel, 
    1818    scale=0.08, 
    19     rpolar=15, requatorial=800, 
    20     sld=.291, solvent_sld=7.105, 
     19    r_polar=15, r_equatorial=800, 
     20    sld=.291, sld_solvent=7.105, 
    2121    background=0, 
    2222    theta=90, phi=0, 
    2323    theta_pd=15, theta_pd_n=40, theta_pd_nsigma=3, 
    24     rpolar_pd=0.222296, rpolar_pd_n=1, rpolar_pd_nsigma=0, 
    25     requatorial_pd=.000128, requatorial_pd_n=1, requatorial_pd_nsigma=0, 
     24    r_polar_pd=0.222296, r_polar_pd_n=1, r_polar_pd_nsigma=0, 
     25    r_equatorial_pd=.000128, r_equatorial_pd_n=1, r_equatorial_pd_nsigma=0, 
    2626    phi_pd=0, phi_pd_n=20, phi_pd_nsigma=3, 
    2727    ) 
    2828 
    2929# SET THE FITTING PARAMETERS 
    30 model.rpolar.range(15, 1000) 
    31 model.requatorial.range(15, 1000) 
     30model.r_polar.range(15, 1000) 
     31model.r_equatorial.range(15, 1000) 
    3232model.theta_pd.range(0, 360) 
    3333model.background.range(0,1000) 

  • sasmodels/bumps_model.py

    rd19962c rea75043  
    219219        """ 
    220220        data, theory, resid = self._data, self.theory(), self.residuals() 
    221         plot_theory(data, theory, resid, view) 
     221        plot_theory(data, theory, resid, view, Iq_calc = self.Iq_calc) 
    222222 
    223223    def simulate_data(self, noise=None): 

  • sasmodels/compare.py

    rd19962c rea75043  
    468468        data = empty_data2D(np.linspace(-qmax, qmax, nq), resolution=res) 
    469469        data.accuracy = opts['accuracy'] 
    470         set_beam_stop(data, 0.004) 
     470        set_beam_stop(data, 0.0004) 
    471471        index = ~data.mask 
    472472    else: 
    … …  
    784784    n1 = int(args[1]) if len(args) > 1 else 1 
    785785    n2 = int(args[2]) if len(args) > 2 else 1 
     786    use_sasview = any(engine=='sasview' and count>0 
     787                      for engine, count in zip(engines, [n1, n2])) 
    786788 
    787789    # Get demo parameters from model definition, or use default parameters 
    … …  
    811813    #import pprint; pprint.pprint(model_info) 
    812814    constrain_pars(model_info, pars) 
    813     constrain_new_to_old(model_info, pars) 
     815    if use_sasview: 
     816        constrain_new_to_old(model_info, pars) 
    814817    if opts['show_pars']: 
    815818        print(str(parlist(model_info, pars, opts['is2d']))) 

  • sasmodels/convert.py

    rce896fd r2a3e1f5  
    138138    namelist = name.split('*') if '*' in name else [name] 
    139139    for name in namelist: 
    140         if name == 'pearl_necklace': 
     140        if name == 'stacked_disks': 
     141            _remove_pd(oldpars, 'n_stacking', name) 
     142        elif name == 'pearl_necklace': 
    141143            _remove_pd(oldpars, 'num_pearls', name) 
    142144            _remove_pd(oldpars, 'thick_string', name) 

  • sasmodels/data.py

    r092cb3c rd6f5da6  
    4141    Load data using a sasview loader. 
    4242    """ 
    43     from sas.dataloader.loader import Loader 
     43    from sas.sascalc.dataloader.loader import Loader 
    4444    loader = Loader() 
    4545    data = loader.load(filename) 
    … …  
    322322 
    323323def plot_theory(data, theory, resid=None, view='log', 
    324                 use_data=True, limits=None): 
     324                use_data=True, limits=None, Iq_calc=None): 
    325325    """ 
    326326    Plot theory calculation. 
    … …  
    343343        _plot_result2D(data, theory, resid, view, use_data, limits=limits) 
    344344    else: 
    345         _plot_result1D(data, theory, resid, view, use_data, limits=limits) 
     345        _plot_result1D(data, theory, resid, view, use_data, 
     346                       limits=limits, Iq_calc=Iq_calc) 
    346347 
    347348 
    … …  
    366367 
    367368@protect 
    368 def _plot_result1D(data, theory, resid, view, use_data, limits=None): 
     369def _plot_result1D(data, theory, resid, view, use_data, 
     370                   limits=None, Iq_calc=None): 
    369371    """ 
    370372    Plot the data and residuals for 1D data. 
    … …  
    376378    use_theory = theory is not None 
    377379    use_resid = resid is not None 
    378     num_plots = (use_data or use_theory) + use_resid 
     380    use_calc = use_theory and Iq_calc is not None 
     381    num_plots = (use_data or use_theory) + use_calc + use_resid 
     382 
    379383 
    380384    scale = data.x**4 if view == 'q4' else 1.0 
    381385 
    382386    if use_data or use_theory: 
     387        if num_plots > 1: 
     388            plt.subplot(1, num_plots, 1) 
     389 
    383390        #print(vmin, vmax) 
    384391        all_positive = True 
    … …  
    399406            if view is 'log': 
    400407                mtheory[mtheory <= 0] = masked 
    401             plt.plot(data.x/10, scale*mtheory, '-', hold=True) 
     408            plt.plot(data.x, scale*mtheory, '-', hold=True) 
    402409            all_positive = all_positive and (mtheory > 0).all() 
    403410            some_present = some_present or (mtheory.count() > 0) 
    … …  
    406413            plt.ylim(*limits) 
    407414 
    408         if num_plots > 1: 
    409             plt.subplot(1, num_plots, 1) 
    410415        plt.xscale('linear' if not some_present else view) 
    411416        plt.yscale('linear' 
    … …  
    415420        plt.ylabel('$I(q)$') 
    416421 
     422    if use_calc: 
     423        # Only have use_calc if have use_theory 
     424        plt.subplot(1, num_plots, 2) 
     425        qx, qy, Iqxy = Iq_calc 
     426        plt.pcolormesh(qx, qy[qy>0], np.log10(Iqxy[qy>0,:])) 
     427        plt.xlabel("$q_x$/A$^{-1}$") 
     428        plt.xlabel("$q_y$/A$^{-1}$") 
     429        plt.xscale('log') 
     430        plt.yscale('log') 
     431        #plt.axis('equal') 
     432 
    417433    if use_resid: 
    418434        mresid = masked_array(resid, data.mask.copy()) 
    … …  
    421437 
    422438        if num_plots > 1: 
    423             plt.subplot(1, num_plots, (use_data or use_theory) + 1) 
     439            plt.subplot(1, num_plots, use_calc + 2) 
    424440        plt.plot(data.x, mresid, '-') 
    425441        plt.xlabel("$q$/A$^{-1}$") 
    … …  
    561577    plottable = masked_array(image, ~valid | data.mask) 
    562578    # Divide range by 10 to convert from angstroms to nanometers 
    563     xmin, xmax = min(data.qx_data)/10, max(data.qx_data)/10 
    564     ymin, ymax = min(data.qy_data)/10, max(data.qy_data)/10 
     579    xmin, xmax = min(data.qx_data), max(data.qx_data) 
     580    ymin, ymax = min(data.qy_data), max(data.qy_data) 
    565581    if view == 'log': 
    566582        vmin, vmax = np.log10(vmin), np.log10(vmax) 
    567583    plt.imshow(plottable.reshape(len(data.x_bins), len(data.y_bins)), 
    568                interpolation='nearest', aspect=1, origin='upper', 
     584               interpolation='nearest', aspect=1, origin='lower', 
    569585               extent=[xmin, xmax, ymin, ymax], vmin=vmin, vmax=vmax) 
    570     plt.xlabel("$q_x$/nm$^{-1}$") 
    571     plt.ylabel("$q_y$/nm$^{-1}$") 
     586    plt.xlabel("$q_x$/A$^{-1}$") 
     587    plt.ylabel("$q_y$/A$^{-1}$") 
    572588    return vmin, vmax 
    573589 

  • sasmodels/direct_model.py

    r60eab2a rea75043  
    6262        elif hasattr(data, 'qx_data'): 
    6363            self.data_type = 'Iqxy' 
     64        elif getattr(data, 'oriented', False): 
     65            self.data_type = 'Iq-oriented' 
    6466        else: 
    6567            self.data_type = 'Iq' 
    … …  
    113115                  and getattr(data, 'dxw', None) is not None): 
    114116                res = resolution.Slit1D(data.x[index], 
    115                                         width=data.dxh[index], 
    116                                         height=data.dxw[index]) 
     117                                        qx_width=data.dxw[index], 
     118                                        qy_width=data.dxl[index]) 
    117119            else: 
    118120                res = resolution.Perfect1D(data.x[index]) 
    … …  
    120122            #self._theory = np.zeros_like(self.Iq) 
    121123            q_vectors = [res.q_calc] 
     124            q_mono = [] 
     125        elif self.data_type == 'Iq-oriented': 
     126            index = (data.x >= data.qmin) & (data.x <= data.qmax) 
     127            if data.y is not None: 
     128                index &= ~np.isnan(data.y) 
     129                Iq = data.y[index] 
     130                dIq = data.dy[index] 
     131            else: 
     132                Iq, dIq = None, None 
     133            if (getattr(data, 'dxl', None) is None 
     134                or getattr(data, 'dxw', None) is None): 
     135                raise ValueError("oriented sample with 1D data needs slit resolution") 
     136 
     137            res = resolution2d.Slit2D(data.x[index], 
     138                                      qx_width=data.dxw[index], 
     139                                      qy_width=data.dxl[index]) 
     140            q_vectors = res.q_calc 
    122141            q_mono = [] 
    123142        else: 
    … …  
    139158        y = Iq + np.random.randn(*dy.shape) * dy 
    140159        self.Iq = y 
    141         if self.data_type == 'Iq': 
     160        if self.data_type in ('Iq', 'Iq-oriented'): 
    142161            self._data.dy[self.index] = dy 
    143162            self._data.y[self.index] = y 
    … …  
    152171        if self._kernel is None: 
    153172            self._kernel = self._model.make_kernel(self._kernel_inputs) 
    154             self._kernel_mono = ( 
    155                 self._model.make_kernel(self._kernel_mono_inputs) 
    156                 if self._kernel_mono_inputs else None 
    157             ) 
     173            self._kernel_mono = (self._model.make_kernel(self._kernel_mono_inputs) 
     174                                 if self._kernel_mono_inputs else None) 
    158175 
    159176        Iq_calc = call_kernel(self._kernel, pars, cutoff=cutoff) 
    160         Iq_mono = (call_kernel(self._kernel_mono, pars, mono=True) 
    161                    if self._kernel_mono_inputs else None) 
     177        # TODO: may want to plot the raw Iq for other than oriented usans 
     178        self.Iq_calc = None 
    162179        if self.data_type == 'sesans': 
     180            Iq_mono = (call_kernel(self._kernel_mono, pars, mono=True) 
     181                       if self._kernel_mono_inputs else None) 
    163182            result = sesans.transform(self._data, 
    164183                                   self._kernel_inputs[0], Iq_calc,  
    … …  
    166185        else: 
    167186            result = self.resolution.apply(Iq_calc) 
     187            if hasattr(self.resolution, 'nx'): 
     188                self.Iq_calc = ( 
     189                    self.resolution.qx_calc, self.resolution.qy_calc, 
     190                    np.reshape(Iq_calc, (self.resolution.ny, self.resolution.nx)) 
     191                ) 
    168192        return result         
    169193 

  • sasmodels/models/_spherepy.py

    r49da079 rd7028dc  
    11r""" 
    2 For information about polarised and magnetic scattering, click here_. 
    3  
    4 .. _here: polar_mag_help.html 
     2For information about polarised and magnetic scattering, see  
     3the :doc:`magnetic help <../sasgui/perspectives/fitting/mag_help>` documentation. 
    54 
    65Definition 

  • sasmodels/models/fuzzy_sphere.py

    r0cc31e1 rd7028dc  
    11r""" 
    2 For information about polarised and magnetic scattering, click here_. 
    3  
    4 .. _here: polar_mag_help.html 
     2For information about polarised and magnetic scattering, see  
     3the :doc:`magnetic help <../sasgui/perspectives/fitting/mag_help>` documentation. 
    54 
    65Definition 

  • sasmodels/models/multilayer_vesicle.py

    rc6ca41e rd7028dc  
    2929    is used as the effective radius for *S(Q)* when $P(Q) * S(Q)$ is applied. 
    3030 
    31  
    32  
    33 For information about polarised and magnetic scattering, click here_. 
    34  
    35 .. _here: polar_mag_help.html 
     31For information about polarised and magnetic scattering, see  
     32the :doc:`magnetic help <../sasgui/perspectives/fitting/mag_help>` documentation. 
    3633 
    3734This code is based on the form factor calculations implemented in the NIST 

  • sasmodels/models/sphere.py

    r49da079 rd7028dc  
    11r""" 
    2 For information about polarised and magnetic scattering, click here_. 
    3  
    4 .. _here: polar_mag_help.html 
     2For information about polarised and magnetic scattering, see  
     3the :doc:`magnetic help <../sasgui/perspectives/fitting/mag_help>` documentation. 
    54 
    65Definition 

  • sasmodels/models/stacked_disks.py

    re664a11 r53215cf  
    100100J S Higgins and H C Benoit, *Polymers and Neutron Scattering*, Clarendon, Oxford, 1994 
    101101 
     102**Author:** NIST IGOR/DANSE **on:** pre 2010 
     103 
     104**Last Modified by:** Piotr Rozyczko **on:** February 18, 2016 
     105 
     106**Last Reviewed by:** Richard Heenan **on:** March 22, 2016 
    102107""" 
    103108 
    … …  
    157162 
    158163tests = [ 
    159     # Accuracy tests based on content in test/utest_extra_models.py 
     164    # Accuracy tests based on content in test/utest_extra_models.py. 
     165    # Added 2 tests with n_stacked = 5 using SasView 3.1.2 - PDB 
    160166    [{'core_thick': 10.0, 
    161167      'layer_thick': 15.0, 
    … …  
    171177      'background': 0.001, 
    172178     }, 0.001, 5075.12], 
     179 
     180    [{'core_thick': 10.0, 
     181      'layer_thick': 15.0, 
     182      'radius': 3000.0, 
     183      'n_stacking': 5.0, 
     184      'sigma_d': 0.0, 
     185      'sld_core': 4.0, 
     186      'sld_layer': -0.4, 
     187      'solvent_sd': 5.0, 
     188      'theta': 0.0, 
     189      'phi': 0.0, 
     190      'scale': 0.01, 
     191      'background': 0.001, 
     192     }, 0.001, 5065.12793824], 
     193 
     194    [{'core_thick': 10.0, 
     195      'layer_thick': 15.0, 
     196      'radius': 3000.0, 
     197      'n_stacking': 5.0, 
     198      'sigma_d': 0.0, 
     199      'sld_core': 4.0, 
     200      'sld_layer': -0.4, 
     201      'solvent_sd': 5.0, 
     202      'theta': 0.0, 
     203      'phi': 0.0, 
     204      'scale': 0.01, 
     205      'background': 0.001, 
     206     }, 0.164, 0.0127673597265], 
    173207 
    174208    [{'core_thick': 10.0, 

  • sasmodels/resolution.py

    r48fbd50 rea75043  
    8282        self.q_calc = (pinhole_extend_q(q, q_width, nsigma=nsigma) 
    8383                       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)) 
     84        self.weight_matrix = pinhole_resolution(self.q_calc, self.q, 
     85                                np.maximum(q_width, MINIMUM_RESOLUTION)) 
    8686 
    8787    def apply(self, theory): 
    … …  
    9191class Slit1D(Resolution): 
    9292    """ 
    93     Slit aperture with a complicated resolution function. 
     93    Slit aperture with resolution function. 
    9494 
    9595    *q* points at which the data is measured. 
    9696 
    97     *qx_width* slit width 
    98  
    99     *qy_height* slit height 
     97    *dqx* slit width in qx 
     98 
     99    *dqy* slit height in qy 
    100100 
    101101    *q_calc* is the list of points to calculate, or None if this should 
    … …  
    104104    The *weight_matrix* is computed by :func:`slit1d_resolution` 
    105105    """ 
    106     def __init__(self, q, width, height=0., q_calc=None): 
    107         # Remember what width/height was used even though we won't need them 
     106    def __init__(self, q, qx_width, qy_width=0., q_calc=None): 
     107        # Remember what width/dqy was used even though we won't need them 
    108108        # after the weight matrix is constructed 
    109         self.width, self.height = width, height 
     109        self.qx_width, self.qy_width = qx_width, qy_width 
    110110 
    111111        # Allow independent resolution on each point even though it is not 
    112112        # needed in practice. 
    113         if np.isscalar(width): 
    114             width = np.ones(len(q))*width 
     113        if np.isscalar(qx_width): 
     114            qx_width = np.ones(len(q))*qx_width 
    115115        else: 
    116             width = np.asarray(width) 
    117         if np.isscalar(height): 
    118             height = np.ones(len(q))*height 
     116            qx_width = np.asarray(qx_width) 
     117        if np.isscalar(qy_width): 
     118            qy_width = np.ones(len(q))*qy_width 
    119119        else: 
    120             height = np.asarray(height) 
     120            qy_width = np.asarray(qy_width) 
    121121 
    122122        self.q = q.flatten() 
    123         self.q_calc = slit_extend_q(q, width, height) \ 
     123        self.q_calc = slit_extend_q(q, qx_width, qy_width) \ 
    124124            if q_calc is None else np.sort(q_calc) 
    125125        self.weight_matrix = \ 
    126             slit_resolution(self.q_calc, self.q, width, height) 
     126            slit_resolution(self.q_calc, self.q, qx_width, qy_width) 
    127127 
    128128    def apply(self, theory): 
    … …  
    396396    """ 
    397397    q = np.sort(q) 
    398     if q_min < q[0]: 
     398    if q_min + 2*MINIMUM_RESOLUTION < q[0]: 
    399399        if q_min <= 0: q_min = q_min*MIN_Q_SCALE_FOR_NEGATIVE_Q_EXTRAPOLATION 
    400400        n_low = np.ceil((q[0]-q_min) / (q[1]-q[0])) if q[1] > q[0] else 15 
    … …  
    402402    else: 
    403403        q_low = [] 
    404     if q_max > q[-1]: 
     404    if q_max - 2*MINIMUM_RESOLUTION > q[-1]: 
    405405        n_high = np.ceil((q_max-q[-1]) / (q[-1]-q[-2])) if q[-1] > q[-2] else 15 
    406406        q_high = np.linspace(q[-1], q_max, n_high+1)[1:] 
    … …  
    596596        Slit smearing with perfect resolution. 
    597597        """ 
    598         resolution = Slit1D(self.x, width=0, height=0, q_calc=self.x) 
     598        resolution = Slit1D(self.x, qx_width=0, qy_width=0, q_calc=self.x) 
    599599        theory = self.Iq(resolution.q_calc) 
    600600        output = resolution.apply(theory) 
    … …  
    606606        Slit smearing with height 0.005 
    607607        """ 
    608         resolution = Slit1D(self.x, width=0, height=0.005, q_calc=self.x) 
     608        resolution = Slit1D(self.x, qx_width=0, qy_width=0.005, q_calc=self.x) 
    609609        theory = self.Iq(resolution.q_calc) 
    610610        output = resolution.apply(theory) 
    … …  
    621621        """ 
    622622        q = np.logspace(-4, -1, 10) 
    623         resolution = Slit1D(q, width=0.2, height=np.inf) 
     623        resolution = Slit1D(q, qx_width=0.2, qy_width=np.inf) 
    624624        theory = 1000*self.Iq(resolution.q_calc**4) 
    625625        output = resolution.apply(theory) 
    … …  
    635635        Slit smearing with width 0.0002 
    636636        """ 
    637         resolution = Slit1D(self.x, width=0.0002, height=0, q_calc=self.x) 
     637        resolution = Slit1D(self.x, qx_width=0.0002, qy_width=0, q_calc=self.x) 
    638638        theory = self.Iq(resolution.q_calc) 
    639639        output = resolution.apply(theory) 
    … …  
    648648        Slit smearing with width > 100*height. 
    649649        """ 
    650         resolution = Slit1D(self.x, width=0.0002, height=0.000001, 
     650        resolution = Slit1D(self.x, qx_width=0.0002, qy_width=0.000001, 
    651651                            q_calc=self.x) 
    652652        theory = self.Iq(resolution.q_calc) 
    … …  
    773773        data = np.loadtxt(data_string.split('\n')).T 
    774774        q, delta_qv, _, answer = data 
    775         resolution = Slit1D(q, width=delta_qv, height=0) 
     775        resolution = Slit1D(q, qx_width=delta_qv, qy_width=0) 
    776776        output = self._eval_sphere(pars, resolution) 
    777777        # TODO: eliminate Igor test since it is too inaccurate to be useful. 
    … …  
    793793        q_calc = slit_extend_q(interpolate(q, 2*np.pi/radius/20), 
    794794                               delta_qv[0], 0.) 
    795         resolution = Slit1D(q, width=delta_qv, height=0, q_calc=q_calc) 
     795        resolution = Slit1D(q, qx_width=delta_qv, qy_width=0, q_calc=q_calc) 
    796796        output = self._eval_sphere(pars, resolution) 
    797797        # TODO: relative error should be lower 
    … …  
    811811        q = np.logspace(log10(4e-5), log10(2.5e-2), 68) 
    812812        width, height = 0.117, 0. 
    813         resolution = Slit1D(q, width=width, height=height) 
     813        resolution = Slit1D(q, qx_width=width, qy_width=height) 
    814814        answer = romberg_slit_1d(q, width, height, form, pars) 
    815815        output = resolution.apply(eval_form(resolution.q_calc, form, pars)) 
    … …  
    10681068 
    10691069    if isinstance(resolution, Slit1D): 
    1070         width, height = resolution.width, resolution.height 
     1070        width, height = resolution.dqx, resolution.dqy 
    10711071        Iq_romb = romberg_slit_1d(resolution.q, width, height, model, pars) 
    10721072    else: 

  • sasmodels/resolution2d.py

    r823e620 rd6f5da6  
    1010from numpy import pi, cos, sin, sqrt 
    1111 
     12from . import resolution 
    1213from .resolution import Resolution 
    1314 
    … …  
    2021NR = {'xhigh':10, 'high':5, 'med':5, 'low':3} 
    2122NPHI = {'xhigh':20, 'high':12, 'med':6, 'low':4} 
     23 
     24## Defaults 
     25N_SLIT_PERP = {'xhigh':1000, 'high':500, 'med':200, 'low':50} 
     26N_SLIT_PERP_DOC = ", ".join("%s=%d"%(name,value) for value,name in 
     27                            sorted((2*v+1,k) for k,v in N_SLIT_PERP.items())) 
    2228 
    2329class Pinhole2D(Resolution): 
    … …  
    167173        else: 
    168174            return theory 
     175 
     176 
     177class Slit2D(Resolution): 
     178    """ 
     179    Slit aperture with resolution function on an oriented sample. 
     180 
     181    *q* points at which the data is measured. 
     182 
     183    *qx_width* slit width in qx 
     184 
     185    *qy_width* slit height in qy; current implementation requires a fixed 
     186    qy_width for all q points. 
     187 
     188    *q_calc* is the list of q points to calculate, or None if this 
     189    should be estimated from the *q* and *qx_width*. 
     190 
     191    *accuracy* determines the number of *qy* points to compute for each *q*. 
     192    The values are stored in sasmodels.resolution2d.N_SLIT_PERP.  The default 
     193    values are: %s 
     194    """ 
     195    __doc__ = __doc__%N_SLIT_PERP_DOC 
     196    def __init__(self, q, qx_width, qy_width=0., q_calc=None, accuracy='low'): 
     197        # Remember what q and width was used even though we won't need them 
     198        # after the weight matrix is constructed 
     199        self.q, self.qx_width, self.qy_width = q, qx_width, qy_width 
     200 
     201        # Allow independent resolution on each qx point even though it is not 
     202        # needed in practice.  Set qy_width to the maximum qy width. 
     203        if np.isscalar(qx_width): 
     204            qx_width = np.ones(len(q))*qx_width 
     205        else: 
     206            qx_width = np.asarray(qx_width) 
     207        if not np.isscalar(qy_width): 
     208            qy_width = np.max(qy_width) 
     209 
     210        # Build grid of qx, qy points 
     211        if q_calc is not None: 
     212            qx_calc = np.sort(q_calc) 
     213        else: 
     214            qx_calc = resolution.pinhole_extend_q(q, qx_width, nsigma=3) 
     215        qy_min, qy_max = np.log10(np.min(q)), np.log10(qy_width) 
     216        qy_calc = np.logspace(qy_min, qy_max, N_SLIT_PERP[accuracy]) 
     217        qy_calc = np.hstack((-qy_calc[::-1], 0, qy_calc)) 
     218        self.q_calc = [v.flatten() for v in np.meshgrid(qx_calc, qy_calc)] 
     219        self.qx_calc, self.qy_calc = qx_calc, qy_calc 
     220        self.nx, self.ny = len(qx_calc), len(qy_calc) 
     221        self.dy = 2*qy_width/self.ny 
     222 
     223        # Build weight matrix for resolution integration 
     224        if np.any(qx_width > 0): 
     225            self.weights = resolution.pinhole_resolution(qx_calc, q, 
     226                    np.maximum(qx_width, resolution.MINIMUM_RESOLUTION)) 
     227        elif len(qx_calc)==len(q) and np.all(qx_calc == q): 
     228            self.weights = None 
     229        else: 
     230            raise ValueError("Slit2D fails with q_calc != q") 
     231 
     232    def apply(self, theory): 
     233        Iq = np.trapz(theory.reshape(self.ny, self.nx), axis=0, x=self.qy_calc) 
     234        if self.weights is not None: 
     235            Iq = resolution.apply_resolution_matrix(self.weights, Iq) 
     236        return Iq 
     237 
Note: See TracChangeset for help on using the changeset viewer.