source: sasview/src/sas/sasgui/perspectives/fitting/model_thread.py @ 5251ec6

magnetic_scattrelease-4.2.2ticket-1009ticket-1249
Last change on this file since 5251ec6 was 5251ec6, checked in by Paul Kienzle <pkienzle@…>, 13 months ago

improved support for py37 in sasgui

  • Property mode set to 100644
File size: 10.0 KB
Line 
1"""
2Calculation thread for modeling
3"""
4
5import time
6import math
7
8import numpy as np
9
10from sas.sascalc.data_util.calcthread import CalcThread
11from sas.sascalc.fit.MultiplicationModel import MultiplicationModel
12
13class Calc2D(CalcThread):
14    """
15    Compute 2D model
16    This calculation assumes a 2-fold symmetry of the model
17    where points are computed for one half of the detector
18    and I(qx, qy) = I(-qx, -qy) is assumed.
19    """
20    def __init__(self, data, model, smearer, qmin, qmax, page_id,
21                 state=None,
22                 weight=None,
23                 fid=None,
24                 toggle_mode_on=False,
25                 completefn=None,
26                 updatefn=None,
27                 update_chisqr=True,
28                 source='model',
29                 yieldtime=0.04,
30                 worktime=0.04,
31                 exception_handler=None,
32                ):
33        CalcThread.__init__(self, completefn, updatefn, yieldtime, worktime,
34                            exception_handler=exception_handler)
35        self.qmin = qmin
36        self.qmax = qmax
37        self.weight = weight
38        self.fid = fid
39        #self.qstep = qstep
40        self.toggle_mode_on = toggle_mode_on
41        self.data = data
42        self.page_id = page_id
43        self.state = None
44        # the model on to calculate
45        self.model = model
46        self.smearer = smearer
47        self.starttime = 0
48        self.update_chisqr = update_chisqr
49        self.source = source
50
51    def compute(self):
52        """
53        Compute the data given a model function
54        """
55        self.starttime = time.time()
56        # Determine appropriate q range
57        if self.qmin is None:
58            self.qmin = 0
59        if self.qmax is None:
60            if self.data is not None:
61                newx = max(math.fabs(self.data.xmax), math.fabs(self.data.xmin))
62                newy = max(math.fabs(self.data.ymax), math.fabs(self.data.ymin))
63                self.qmax = math.sqrt(newx**2 + newy**2)
64
65        if self.data is None:
66            msg = "Compute Calc2D receive data = %s.\n" % str(self.data)
67            raise ValueError(msg)
68
69        # Define matrix where data will be plotted
70        radius = np.sqrt(self.data.qx_data**2 + self.data.qy_data**2)
71
72        # For theory, qmax is based on 1d qmax
73        # so that must be mulitified by sqrt(2) to get actual max for 2d
74        index_model = (self.qmin <= radius) & (radius <= self.qmax)
75        index_model &= self.data.mask
76        index_model &= np.isfinite(self.data.data)
77
78        if self.smearer is not None:
79            # Set smearer w/ data, model and index.
80            fn = self.smearer
81            fn.set_model(self.model)
82            fn.set_index(index_model)
83            # Calculate smeared Intensity
84            #(by Gaussian averaging): DataLoader/smearing2d/Smearer2D()
85            value = fn.get_value()
86        else:
87            # calculation w/o smearing
88            value = self.model.evalDistribution([
89                self.data.qx_data[index_model],
90                self.data.qy_data[index_model]
91            ])
92        # Initialize output to NaN so masked elements do not get plotted.
93        output = np.empty_like(self.data.qx_data)
94        # output default is None
95        # This method is to distinguish between masked
96        #point(nan) and data point = 0.
97        output[:] = np.NaN
98        # set value for self.mask==True, else still None to Plottools
99        output[index_model] = value
100        elapsed = time.time() - self.starttime
101        self.complete(image=output,
102                      data=self.data,
103                      page_id=self.page_id,
104                      model=self.model,
105                      state=self.state,
106                      toggle_mode_on=self.toggle_mode_on,
107                      elapsed=elapsed,
108                      index=index_model,
109                      fid=self.fid,
110                      qmin=self.qmin,
111                      qmax=self.qmax,
112                      weight=self.weight,
113                      #qstep=self.qstep,
114                      update_chisqr=self.update_chisqr,
115                      source=self.source)
116
117
118class Calc1D(CalcThread):
119    """
120    Compute 1D data
121    """
122    def __init__(self, model,
123                 page_id,
124                 data,
125                 fid=None,
126                 qmin=None,
127                 qmax=None,
128                 weight=None,
129                 smearer=None,
130                 toggle_mode_on=False,
131                 state=None,
132                 completefn=None,
133                 update_chisqr=True,
134                 source='model',
135                 updatefn=None,
136                 yieldtime=0.01,
137                 worktime=0.01,
138                 exception_handler=None,
139                ):
140        """
141        """
142        CalcThread.__init__(self, completefn, updatefn, yieldtime, worktime,
143                            exception_handler=exception_handler)
144        self.fid = fid
145        self.data = data
146        self.qmin = qmin
147        self.qmax = qmax
148        self.model = model
149        self.weight = weight
150        self.toggle_mode_on = toggle_mode_on
151        self.state = state
152        self.page_id = page_id
153        self.smearer = smearer
154        self.starttime = 0
155        self.update_chisqr = update_chisqr
156        self.source = source
157        self.out = None
158        self.index = None
159
160    def compute(self):
161        """
162        Compute model 1d value given qmin , qmax , x value
163        """
164        self.starttime = time.time()
165        output = np.zeros((len(self.data.x)))
166        index = (self.qmin <= self.data.x) & (self.data.x <= self.qmax)
167
168        # If we use a smearer, also return the unsmeared model
169        unsmeared_output = None
170        unsmeared_data = None
171        unsmeared_error = None
172        ##smearer the ouput of the plot
173        if self.smearer is not None:
174            first_bin, last_bin = self.smearer.get_bin_range(self.qmin,
175                                                             self.qmax)
176            mask = self.data.x[first_bin:last_bin+1]
177            unsmeared_output = np.zeros((len(self.data.x)))
178            unsmeared_output[first_bin:last_bin+1] = self.model.evalDistribution(mask)
179            self.smearer.model = self.model
180            output = self.smearer(unsmeared_output, first_bin, last_bin)
181
182            # Rescale data to unsmeared model
183            # Check that the arrays are compatible. If we only have a model but no data,
184            # the length of data.y will be zero.
185            if isinstance(self.data.y, np.ndarray) and output.shape == self.data.y.shape:
186                unsmeared_data = np.zeros((len(self.data.x)))
187                unsmeared_error = np.zeros((len(self.data.x)))
188                unsmeared_data[first_bin:last_bin+1] = self.data.y[first_bin:last_bin+1]\
189                                                        * unsmeared_output[first_bin:last_bin+1]\
190                                                        / output[first_bin:last_bin+1]
191                unsmeared_error[first_bin:last_bin+1] = self.data.dy[first_bin:last_bin+1]\
192                                                        * unsmeared_output[first_bin:last_bin+1]\
193                                                        / output[first_bin:last_bin+1]
194                unsmeared_output = unsmeared_output[index]
195                unsmeared_data = unsmeared_data[index]
196                unsmeared_error = unsmeared_error
197        else:
198            output[index] = self.model.evalDistribution(self.data.x[index])
199
200        x=self.data.x[index]
201        y=output[index]
202        sq_values = None
203        pq_values = None
204        if isinstance(self.model, MultiplicationModel):
205            s_model = self.model.s_model
206            p_model = self.model.p_model
207            sq_values = s_model.evalDistribution(x)
208            pq_values = p_model.evalDistribution(x)
209        elif hasattr(self.model, "calc_composition_models"):
210            results = self.model.calc_composition_models(x)
211            if results is not None:
212                pq_values, sq_values = results
213
214
215        elapsed = time.time() - self.starttime
216
217        self.complete(x=x, y=y,
218                      page_id=self.page_id,
219                      state=self.state,
220                      weight=self.weight,
221                      fid=self.fid,
222                      toggle_mode_on=self.toggle_mode_on,
223                      elapsed=elapsed, index=index, model=self.model,
224                      data=self.data,
225                      update_chisqr=self.update_chisqr,
226                      source=self.source,
227                      unsmeared_model=unsmeared_output,
228                      unsmeared_data=unsmeared_data,
229                      unsmeared_error=unsmeared_error,
230                      pq_model=pq_values,
231                      sq_model=sq_values)
232
233    def results(self):
234        """
235        Send resuts of the computation
236        """
237        return [self.out, self.index]
238
239"""
240Example: ::
241
242    class CalcCommandline:
243        def __init__(self, n=20000):
244            #print(thread.get_ident())
245
246            from sasmodels.sasview_model import _make_standard_model
247            cylinder = _make_standard_model('cylinder')
248            model = cylinder()
249
250            print(model.runXY([0.01, 0.02]))
251
252            qmax = 0.01
253            qstep = 0.0001
254            self.done = False
255
256            x = numpy.arange(-qmax, qmax+qstep*0.01, qstep)
257            y = numpy.arange(-qmax, qmax+qstep*0.01, qstep)
258
259            calc_thread_2D = Calc2D(x, y, None, model.clone(),None,
260                                    -qmax, qmax,qstep,
261                                    completefn=self.complete,
262                                    updatefn=self.update ,
263                                    yieldtime=0.0)
264
265            calc_thread_2D.queue()
266            calc_thread_2D.ready(2.5)
267
268            while not self.done:
269                time.sleep(1)
270
271        def update(self,output):
272            print("update")
273
274        def complete(self, image, data, model, elapsed, qmin, qmax,index, qstep ):
275            print("complete")
276            self.done = True
277
278    if __name__ == "__main__":
279        CalcCommandline()
280"""
Note: See TracBrowser for help on using the repository browser.