source: sasview/src/sas/sascalc/data_util/qsmearing.py @ d3911e3

ESS_GUIESS_GUI_DocsESS_GUI_batch_fittingESS_GUI_bumps_abstractionESS_GUI_iss1116ESS_GUI_iss879ESS_GUI_iss959ESS_GUI_openclESS_GUI_orderingESS_GUI_sync_sascalccostrafo411magnetic_scattrelease-4.1.1release-4.1.2release-4.2.2ticket-1009ticket-1094-headlessticket-1242-2d-resolutionticket-1243ticket-1249ticket885unittest-saveload
Last change on this file since d3911e3 was d3911e3, checked in by Paul Kienzle <pkienzle@…>, 3 years ago

Reproduce sasview smear 2D wrapper so 2D fits work. Fixes #811, #825.

  • Property mode set to 100644
File size: 7.0 KB
Line 
1"""
2    Handle Q smearing
3"""
4#####################################################################
5#This software was developed by the University of Tennessee as part of the
6#Distributed Data Analysis of Neutron Scattering Experiments (DANSE)
7#project funded by the US National Science Foundation.
8#See the license text in license.txt
9#copyright 2008, University of Tennessee
10######################################################################
11import numpy
12import math
13import logging
14import sys
15
16from sasmodels.resolution import Slit1D, Pinhole1D
17from sasmodels.resolution2d import Pinhole2D
18
19def smear_selection(data, model = None):
20    """
21    Creates the right type of smearer according
22    to the data.
23    The canSAS format has a rule that either
24    slit smearing data OR resolution smearing data
25    is available.
26
27    For the present purpose, we choose the one that
28    has none-zero data. If both slit and resolution
29    smearing arrays are filled with good data
30    (which should not happen), then we choose the
31    resolution smearing data.
32
33    :param data: Data1D object
34    :param model: sas.model instance
35    """
36    # Sanity check. If we are not dealing with a SAS Data1D
37    # object, just return None
38    if  data.__class__.__name__ not in ['Data1D', 'Theory1D']:
39        if data == None:
40            return None
41        elif data.dqx_data == None or data.dqy_data == None:
42            return None
43        return PySmear2D(data, model)
44
45    if  not hasattr(data, "dx") and not hasattr(data, "dxl")\
46         and not hasattr(data, "dxw"):
47        return None
48
49    # Look for resolution smearing data
50    _found_resolution = False
51    if data.dx is not None and len(data.dx) == len(data.x):
52
53        # Check that we have non-zero data
54        if data.dx[0] > 0.0:
55            _found_resolution = True
56            #print "_found_resolution",_found_resolution
57            #print "data1D.dx[0]",data1D.dx[0],data1D.dxl[0]
58    # If we found resolution smearing data, return a QSmearer
59    if _found_resolution == True:
60         return pinhole_smear(data, model)
61
62    # Look for slit smearing data
63    _found_slit = False
64    if data.dxl is not None and len(data.dxl) == len(data.x) \
65        and data.dxw is not None and len(data.dxw) == len(data.x):
66
67        # Check that we have non-zero data
68        if data.dxl[0] > 0.0 or data.dxw[0] > 0.0:
69            _found_slit = True
70
71        # Sanity check: all data should be the same as a function of Q
72        for item in data.dxl:
73            if data.dxl[0] != item:
74                _found_resolution = False
75                break
76
77        for item in data.dxw:
78            if data.dxw[0] != item:
79                _found_resolution = False
80                break
81    # If we found slit smearing data, return a slit smearer
82    if _found_slit == True:
83        return slit_smear(data, model)
84    return None
85
86
87class PySmear(object):
88    """
89    Wrapper for pure python sasmodels resolution functions.
90    """
91    def __init__(self, resolution, model):
92        self.model = model
93        self.resolution = resolution
94        self.offset = numpy.searchsorted(self.resolution.q_calc, self.resolution.q[0])
95
96    def apply(self, iq_in, first_bin=0, last_bin=None):
97        """
98        Apply the resolution function to the data.
99        Note that this is called with iq_in matching data.x, but with
100        iq_in[first_bin:last_bin] set to theory values for these bins,
101        and the remainder left undefined.  The first_bin, last_bin values
102        should be those returned from get_bin_range.
103        The returned value is of the same length as iq_in, with the range
104        first_bin:last_bin set to the resolution smeared values.
105        """
106        if last_bin is None: last_bin = len(iq_in)
107        start, end = first_bin + self.offset, last_bin + self.offset
108        q_calc = self.resolution.q_calc
109        iq_calc = numpy.empty_like(q_calc)
110        if start > 0:
111            iq_calc[:start] = self.model.evalDistribution(q_calc[:start])
112        if end+1 < len(q_calc):
113            iq_calc[end+1:] = self.model.evalDistribution(q_calc[end+1:])
114        iq_calc[start:end+1] = iq_in[first_bin:last_bin+1]
115        smeared = self.resolution.apply(iq_calc)
116        return smeared
117    __call__ = apply
118
119    def get_bin_range(self, q_min=None, q_max=None):
120        """
121        For a given q_min, q_max, find the corresponding indices in the data.
122        Returns first, last.
123        Note that these are indexes into q from the data, not the q_calc
124        needed by the resolution function.  Note also that these are the
125        indices, not the range limits.  That is, the complete range will be
126        q[first:last+1].
127        """
128        q = self.resolution.q
129        first = numpy.searchsorted(q, q_min)
130        last = numpy.searchsorted(q, q_max)
131        return first, min(last,len(q)-1)
132
133def slit_smear(data, model=None):
134    q = data.x
135    width = data.dxw if data.dxw is not None else 0
136    height = data.dxl if data.dxl is not None else 0
137    # TODO: width and height seem to be reversed
138    return PySmear(Slit1D(q, height, width), model)
139
140def pinhole_smear(data, model=None):
141    q = data.x
142    width = data.dx if data.dx is not None else 0
143    return PySmear(Pinhole1D(q, width), model)
144
145
146class PySmear2D(object):
147    """
148    Q smearing class for SAS 2d pinhole data
149    """
150
151    def __init__(self, data=None, model=None):
152        self.data = data
153        self.model = model
154        self.accuracy = 'Low'
155        self.limit = 3.0
156        self.index = None
157        self.coords = 'polar'
158        self.smearer = True
159
160    def set_accuracy(self, accuracy='Low'):
161        """
162        Set accuracy.
163
164        :param accuracy:  string
165        """
166        self.accuracy = accuracy
167
168    def set_smearer(self, smearer=True):
169        """
170        Set whether or not smearer will be used
171
172        :param smearer: smear object
173
174        """
175        self.smearer = smearer
176
177    def set_data(self, data=None):
178        """
179        Set data.
180
181        :param data: DataLoader.Data_info type
182        """
183        self.data = data
184
185    def set_model(self, model=None):
186        """
187        Set model.
188
189        :param model: sas.models instance
190        """
191        self.model = model
192
193    def set_index(self, index=None):
194        """
195        Set index.
196
197        :param index: 1d arrays
198        """
199        self.index = index
200
201    def get_value(self):
202        """
203        Over sampling of r_nbins times phi_nbins, calculate Gaussian weights,
204        then find smeared intensity
205        """
206        if self.smearer:
207            res = Pinhole2D(data=self.data, index=self.index,
208                            nsigma=3.0, accuracy=self.accuracy,
209                            coords=self.coords)
210            val = self.model.evalDistribution(res.q_calc)
211            return res.apply(val)
212        else:
213            index = self.index if self.index is not None else slice(None)
214            qx_data = self.data.qx_data[index]
215            qy_data = self.data.qy_data[index]
216            q_calc = [qx_data, qy_data]
217            val = self.model.evalDistribution(q_calc)
218            return val
219
Note: See TracBrowser for help on using the repository browser.