Changes in / [ba51e00:fbaef04] in sasmodels


Ignore:
Files:
11 edited

Legend:

Unmodified
Added
Removed
  • doc/guide/scripting.rst

    rbd7630d r4aa5dce  
    1010The key functions are :func:`sasmodels.core.load_model` for loading the 
    1111model definition and compiling the kernel and 
    12 :func:`sasmodels.data.load_data` for calling sasview to load the data. 
     12:func:`sasmodels.data.load_data` for calling sasview to load the data. Need 
     13the data because that defines the resolution function and the q values to 
     14evaluate. If there is no data, then use :func:`sasmodels.data.empty_data1D` 
     15or :func:`sasmodels.data.empty_data2D` to create some data with a given $q$. 
    1316 
    14 Preparing data 
    15 ============== 
     17Using sasmodels through bumps 
     18============================= 
    1619 
    17 Usually you will load data via the sasview loader, with the 
    18 :func:`sasmodels.data.load_data` function.  For example:: 
    19  
    20     from sasmodels.data import load_data 
    21     data = load_data("sasmodels/example/093191_201.dat") 
    22  
    23 You may want to apply a data mask, such a beam stop, and trim high $q$:: 
    24  
    25     from sasmodels.data import set_beam_stop 
    26     set_beam_stop(data, qmin, qmax) 
    27  
    28 The :func:`sasmodels.data.set_beam_stop` method simply sets the *mask* 
    29 attribute for the data. 
    30  
    31 The data defines the resolution function and the q values to evaluate, so 
    32 even if you simulating experiments prior to making measurements, you still 
    33 need a data object for reference. Use :func:`sasmodels.data.empty_data1D` 
    34 or :func:`sasmodels.data.empty_data2D` to create a container with a 
    35 given $q$ and $\Delta q/q$.  For example:: 
    36  
    37     import numpy as np 
    38     from sasmodels.data import empty_data1D 
    39  
    40     # 120 points logarithmically spaced from 0.005 to 0.2, with dq/q = 5% 
    41     q = np.logspace(np.log10(5e-3), np.log10(2e-1), 120) 
    42     data = empty_data1D(q, resolution=0.05) 
    43  
    44 To use a more realistic model of resolution, or to load data from a file 
    45 format not understood by SasView, you can use :class:`sasmodels.data.Data1D` 
    46 or :class:`sasmodels.data.Data2D` directly.  The 1D data uses 
    47 *x*, *y*, *dx* and *dy* for $x = q$ and $y = I(q)$, and 2D data uses 
    48 *x*, *y*, *z*, *dx*, *dy*, *dz* for $x, y = qx, qy$ and $z = I(qx, qy)$. 
    49 [Note: internally, the Data2D object uses SasView conventions, 
    50 *qx_data*, *qy_data*, *data*, *dqx_data*, *dqy_data*, and *err_data*.] 
    51  
    52 For USANS data, use 1D data, but set *dxl* and *dxw* attributes to 
    53 indicate slit resolution:: 
    54  
    55     data.dxl = 0.117 
    56  
    57 See :func:`sasmodels.resolution.slit_resolution` for details. 
    58  
    59 SESANS data is more complicated; if your SESANS format is not supported by 
    60 SasView you need to define a number of attributes beyond *x*, *y*.  For 
    61 example:: 
    62  
    63     SElength = np.linspace(0, 2400, 61) # [A] 
    64     data = np.ones_like(SElength) 
    65     err_data = np.ones_like(SElength)*0.03 
    66  
    67     class Source: 
    68         wavelength = 6 # [A] 
    69         wavelength_unit = "A" 
    70     class Sample: 
    71         zacceptance = 0.1 # [A^-1] 
    72         thickness = 0.2 # [cm] 
    73  
    74     class SESANSData1D: 
    75         #q_zmax = 0.23 # [A^-1] 
    76         lam = 0.2 # [nm] 
    77         x = SElength 
    78         y = data 
    79         dy = err_data 
    80         sample = Sample() 
    81     data = SESANSData1D() 
    82  
    83     x, y = ... # create or load sesans 
    84     data = smd.Data 
    85  
    86 The *data* module defines various data plotters as well. 
    87  
    88 Using sasmodels directly 
    89 ======================== 
    90  
    91 Once you have a computational kernel and a data object, you can evaluate 
    92 the model for various parameters using 
    93 :class:`sasmodels.direct_model.DirectModel`.  The resulting object *f* 
    94 will be callable as *f(par=value, ...)*, returning the $I(q)$ for the $q$ 
    95 values in the data.  For example:: 
    96  
    97     import numpy as np 
    98     from sasmodels.data import empty_data1D 
    99     from sasmodels.core import load_model 
    100     from sasmodels.direct_model import DirectModel 
    101  
    102     # 120 points logarithmically spaced from 0.005 to 0.2, with dq/q = 5% 
    103     q = np.logspace(np.log10(5e-3), np.log10(2e-1), 120) 
    104     data = empty_data1D(q, resolution=0.05) 
    105     kernel = load_model("ellipsoid) 
    106     f = DirectModel(data, kernel) 
    107     Iq = f(radius_polar=100) 
    108  
    109 Polydispersity information is set with special parameter names: 
    110  
    111     * *par_pd* for polydispersity width, $\Delta p/p$, 
    112     * *par_pd_n* for the number of points in the distribution, 
    113     * *par_pd_type* for the distribution type (as a string), and 
    114     * *par_pd_nsigmas* for the limits of the distribution. 
    115  
    116 Using sasmodels through the bumps optimizer 
    117 =========================================== 
    118  
    119 Like DirectModel, you can wrap data and a kernel in a *bumps* model with 
     20With the data and the model, you can wrap it in a *bumps* model with 
    12021class:`sasmodels.bumps_model.Model` and create an 
    121 class:`sasmodels.bumps_model.Experiment` that you can fit with the *bumps* 
     22class:`sasmodels.bump_model.Experiment` that you can fit with the *bumps* 
    12223interface. Here is an example from the *example* directory such as 
    12324*example/model.py*:: 
     
    17475    SasViewCom bumps.cli example/model.py --preview 
    17576 
    176 Calling the computation kernel 
    177 ============================== 
     77Using sasmodels directly 
     78======================== 
     79 
     80Bumps has a notion of parameter boxes in which you can set and retrieve 
     81values.  Instead of using bumps, you can create a directly callable function 
     82with :class:`sasmodels.direct_model.DirectModel`.  The resulting object *f* 
     83will be callable as *f(par=value, ...)*, returning the $I(q)$ for the $q$ 
     84values in the data.  Polydisperse parameters use the same naming conventions 
     85as in the bumps model, with e.g., radius_pd being the polydispersity associated 
     86with radius. 
    17887 
    17988Getting a simple function that you can call on a set of q values and return 
  • sasmodels/compare.py

    rbd7630d raf7a97c  
    12881288 
    12891289    if opts['datafile'] is not None: 
    1290         data0 = load_data(os.path.expanduser(opts['datafile'])) 
    1291         data = data0, data0 
     1290        data = load_data(os.path.expanduser(opts['datafile'])) 
    12921291    else: 
    12931292        # Hack around the fact that make_data doesn't take a pair of resolutions 
  • sasmodels/data.py

    rbd7630d rc036ddb  
    183183    *x* is spin echo length and *y* is polarization (P/P0). 
    184184    """ 
    185     isSesans = True 
    186185    def __init__(self, **kw): 
    187186        Data1D.__init__(self, **kw) 
     
    302301        self.wavelength_unit = "A" 
    303302 
    304 class Sample(object): 
    305     """ 
    306     Sample attributes. 
    307     """ 
    308     def __init__(self): 
    309         # type: () -> None 
    310         pass 
    311303 
    312304def empty_data1D(q, resolution=0.0, L=0., dL=0.): 
  • sasmodels/direct_model.py

    r7b9e4dd r1a8c11c  
    3131from . import resolution2d 
    3232from .details import make_kernel_args, dispersion_mesh 
    33 from .modelinfo import DEFAULT_BACKGROUND 
    3433 
    3534# pylint: disable=unused-import 
     
    350349 
    351350        # Need to pull background out of resolution for multiple scattering 
    352         background = pars.get('background', DEFAULT_BACKGROUND) 
     351        background = pars.get('background', 0.) 
    353352        pars = pars.copy() 
    354353        pars['background'] = 0. 
  • sasmodels/generate.py

    r5399809 r5399809  
    991991    docs = model_info.docs if model_info.docs is not None else "" 
    992992    docs = convert_section_titles_to_boldface(docs) 
    993     if model_info.structure_factor: 
    994         pars = model_info.parameters.kernel_parameters 
    995     else: 
    996         pars = model_info.parameters.COMMON + model_info.parameters.kernel_parameters 
    997     partable = make_partable(pars) 
     993    pars = make_partable(model_info.parameters.COMMON 
     994                         + model_info.parameters.kernel_parameters) 
    998995    subst = dict(id=model_info.id.replace('_', '-'), 
    999996                 name=model_info.name, 
    1000997                 title=model_info.title, 
    1001                  parameters=partable, 
     998                 parameters=pars, 
    1002999                 returns=Sq_units if model_info.structure_factor else Iq_units, 
    10031000                 docs=docs) 
  • sasmodels/kernel_iq.c

    r70530778 rc57ee9e  
    192192    QACRotation *rotation, 
    193193    double qx, double qy, 
    194     double *qab_out, double *qc_out) 
     194    double *qa_out, double *qc_out) 
    195195{ 
     196    const double dqc = rotation->R31*qx + rotation->R32*qy; 
    196197    // Indirect calculation of qab, from qab^2 = |q|^2 - qc^2 
    197     const double dqc = rotation->R31*qx + rotation->R32*qy; 
    198     const double dqab_sq = -dqc*dqc + qx*qx + qy*qy; 
    199     //*qab_out = sqrt(fabs(dqab_sq)); 
    200     *qab_out = dqab_sq > 0.0 ? sqrt(dqab_sq) : 0.0; 
     198    const double dqa = sqrt(-dqc*dqc + qx*qx + qy*qy); 
     199 
     200    *qa_out = dqa; 
    201201    *qc_out = dqc; 
    202202} 
  • sasmodels/modelinfo.py

    r5399809 r5399809  
    4545# Note that scale and background cannot be coordinated parameters whose value 
    4646# depends on the some polydisperse parameter with the current implementation 
    47 DEFAULT_BACKGROUND = 1e-3 
    4847COMMON_PARAMETERS = [ 
    4948    ("scale", "", 1, (0.0, np.inf), "", "Source intensity"), 
    50     ("background", "1/cm", DEFAULT_BACKGROUND, (-np.inf, np.inf), "", "Source background"), 
     49    ("background", "1/cm", 1e-3, (-np.inf, np.inf), "", "Source background"), 
    5150] 
    5251assert (len(COMMON_PARAMETERS) == 2 
  • sasmodels/models/ellipsoid.py

    rd277229 rd277229  
    126126from numpy import inf, sin, cos, pi 
    127127 
    128 try: 
    129     from numpy import cbrt 
    130 except ImportError: 
    131     def cbrt(x): return x ** (1.0/3.0) 
    132  
    133128name = "ellipsoid" 
    134129title = "Ellipsoid of revolution with uniform scattering length density." 
  • sasmodels/models/spinodal.py

    r475ff58 r71b751d  
    33---------- 
    44 
    5 This model calculates the SAS signal of a phase separating system  
    6 undergoing spinodal decomposition. The scattering intensity $I(q)$ is calculated  
    7 as  
     5This model calculates the SAS signal of a phase separating solution 
     6under spinodal decomposition. The scattering intensity $I(q)$ is calculated as 
    87 
    98.. math:: 
    109    I(q) = I_{max}\frac{(1+\gamma/2)x^2}{\gamma/2+x^{2+\gamma}}+B 
    1110 
    12 where $x=q/q_0$, $q_0$ is the peak position, $I_{max}$ is the intensity  
    13 at $q_0$ (parameterised as the $scale$ parameter), and $B$ is a flat  
    14 background. The spinodal wavelength is given by $2\pi/q_0$.  
    15  
    16 The exponent $\gamma$ is equal to $d+1$ for off-critical concentration  
    17 mixtures (smooth interfaces) and $2d$ for critical concentration mixtures  
    18 (entangled interfaces), where $d$ is the dimensionality (ie, 1, 2, 3) of the  
    19 system. Thus 2 <= $\gamma$ <= 6. A transition from $\gamma=d+1$ to $\gamma=2d$  
    20 is expected near the percolation threshold.  
    21  
    22 As this function tends to zero as $q$ tends to zero, in practice it may be  
    23 necessary to combine it with another function describing the low-angle  
    24 scattering, or to simply omit the low-angle scattering from the fit. 
     11where $x=q/q_0$ and $B$ is a flat background. The characteristic structure 
     12length scales with the correlation peak at $q_0$. The exponent $\gamma$ is 
     13equal to $d+1$ with d the dimensionality of the off-critical concentration 
     14mixtures. A transition to $\gamma=2d$ is seen near the percolation threshold 
     15into the critical concentration regime. 
    2516 
    2617References 
     
    3122Physica A 123,497 (1984). 
    3223 
    33 Revision History 
    34 ---------------- 
     24Authorship and Verification 
     25---------------------------- 
    3526 
    36 * **Author:**  Dirk Honecker **Date:** Oct 7, 2016 
    37 * **Revised:** Steve King    **Date:** Sep 7, 2018 
     27* **Author:** Dirk Honecker **Date:** Oct 7, 2016 
    3828""" 
    3929 
     
    4434title = "Spinodal decomposition model" 
    4535description = """\ 
    46       I(q) = Imax ((1+gamma/2)x^2)/(gamma/2+x^(2+gamma)) + background 
     36      I(q) = scale ((1+gamma/2)x^2)/(gamma/2+x^(2+gamma))+background 
    4737 
    4838      List of default parameters: 
    49        
    50       Imax = correlation peak intensity at q_0 
    51       background = incoherent background 
    52       gamma = exponent (see model documentation) 
     39      scale = scaling 
     40      gamma = exponent 
     41      x = q/q_0 
    5342      q_0 = correlation peak position [1/A] 
    54       x = q/q_0""" 
    55        
     43      background = Incoherent background""" 
    5644category = "shape-independent" 
    5745 
  • sasmodels/product.py

    r33d7be3 r33d7be3  
    4242    pars = [] 
    4343    if p_info.have_Fq: 
    44         par = parse_parameter( 
    45                 "structure_factor_mode", 
    46                 "", 
    47                 0, 
    48                 [["P*S","P*(1+beta*(S-1))"]], 
    49                 "", 
    50                 "Structure factor calculation") 
     44        par = parse_parameter("structure_factor_mode", "", 0, [["P*S","P*(1+beta*(S-1))"]], "", 
     45                        "Structure factor calculation") 
    5146        pars.append(par) 
    5247    if p_info.effective_radius_type is not None: 
    53         par = parse_parameter( 
    54                 "radius_effective_mode", 
    55                 "", 
    56                 0, 
    57                 [["unconstrained"] + p_info.effective_radius_type], 
    58                 "", 
    59                 "Effective radius calculation") 
     48        par = parse_parameter("radius_effective_mode", "", 0, [["unconstrained"] + p_info.effective_radius_type], "", 
     49                        "Effective radius calculation") 
    6050        pars.append(par) 
    6151    return pars 
     
    179169    )) 
    180170 
    181 def _intermediates_beta(F1,              # type: np.ndarray 
    182                         F2,              # type: np.ndarray 
    183                         S,               # type: np.ndarray 
    184                         scale,           # type: np.ndarray 
    185                         bg,              # type: np.ndarray 
    186                         effective_radius # type: float 
    187                         ): 
    188     # type: (...) -> OrderedDict[str, Union[np.ndarray, float]] 
     171def _intermediates_beta(F1, F2, S, scale, bg, effective_radius): 
     172    # type: (np.ndarray, np.ndarray, np.ndarray, np.ndarray, np.ndarray, float) -> OrderedDict[str, Union[np.ndarray, float]] 
    189173    """ 
    190174    Returns intermediate results for beta approximation-enabled product. The result may be an array or a float. 
  • sasmodels/sasview_model.py

    r84f2962 rd32de68  
    543543            # pylint: disable=unpacking-non-sequence 
    544544            q, phi = x 
    545             result, _ = self.calculate_Iq([q*math.cos(phi)], [q*math.sin(phi)]) 
    546             return result[0] 
    547         else: 
    548             result, _ = self.calculate_Iq([x]) 
    549             return result[0] 
     545            return self.calculate_Iq([q*math.cos(phi)], [q*math.sin(phi)])[0] 
     546        else: 
     547            return self.calculate_Iq([x])[0] 
    550548 
    551549 
     
    562560        """ 
    563561        if isinstance(x, (list, tuple)): 
    564             result, _ = self.calculate_Iq([x[0]], [x[1]]) 
    565             return result[0] 
    566         else: 
    567             result, _ = self.calculate_Iq([x]) 
    568             return result[0] 
     562            return self.calculate_Iq([x[0]], [x[1]])[0] 
     563        else: 
     564            return self.calculate_Iq([x])[0] 
    569565 
    570566    def evalDistribution(self, qdist): 
     
    600596            # Check whether we have a list of ndarrays [qx,qy] 
    601597            qx, qy = qdist 
    602             result, _ = self.calculate_Iq(qx, qy) 
    603             return result 
     598            return self.calculate_Iq(qx, qy) 
    604599 
    605600        elif isinstance(qdist, np.ndarray): 
    606601            # We have a simple 1D distribution of q-values 
    607             result, _ = self.calculate_Iq(qdist) 
    608             return result 
     602            return self.calculate_Iq(qdist) 
    609603 
    610604        else: 
     
    646640        if composition and composition[0] == 'product': # only P*S for now 
    647641            with calculation_lock: 
    648                 _, lazy_results = self._calculate_Iq(qx) 
    649                 # for compatibility with sasview 4.x 
    650                 results = lazy_results() 
    651                 pq_data = results.get("P(Q)") 
    652                 sq_data = results.get("S(Q)") 
    653                 return pq_data, sq_data 
     642                self._calculate_Iq(qx) 
     643                # for compatibility with sasview 4.3 
     644                results = self._intermediate_results() 
     645                return results["P(Q)"], results["S(Q)"] 
    654646        else: 
    655647            return None 
    656648 
    657     def calculate_Iq(self, 
    658                      qx,     # type: Sequence[float] 
    659                      qy=None # type: Optional[Sequence[float]] 
    660                      ): 
    661         # type: (...) -> Tuple[np.ndarray, Callable[[], collections.OrderedDict[str, np.ndarray]]] 
     649    def calculate_Iq(self, qx, qy=None): 
     650        # type: (Sequence[float], Optional[Sequence[float]]) -> np.ndarray 
    662651        """ 
    663652        Calculate Iq for one set of q with the current parameters. 
     
    667656        This should NOT be used for fitting since it copies the *q* vectors 
    668657        to the card for each evaluation. 
    669  
    670         The returned tuple contains the scattering intensity followed by a 
    671         callable which returns a dictionary of intermediate data from 
    672         ProductKernel. 
    673658        """ 
    674659        ## uncomment the following when trying to debug the uncoordinated calls 
     
    703688        result = calculator(call_details, values, cutoff=self.cutoff, 
    704689                            magnetic=is_magnetic) 
    705         lazy_results = getattr(calculator, 'results', 
    706                                lambda: collections.OrderedDict()) 
    707690        #print("result", result) 
    708  
     691        self._intermediate_results = getattr(calculator, 'results', None) 
    709692        calculator.release() 
    710693        #self._model.release() 
    711  
    712         return result, lazy_results 
     694        return result 
    713695 
    714696    def calculate_ER(self): 
     
    899881    q = np.linspace(-0.35, 0.35, 500) 
    900882    qx, qy = np.meshgrid(q, q) 
    901     result, _ = model.calculate_Iq(qx.flatten(), qy.flatten()) 
     883    result = model.calculate_Iq(qx.flatten(), qy.flatten()) 
    902884    result = result.reshape(qx.shape) 
    903885 
Note: See TracChangeset for help on using the changeset viewer.