source: sasview/src/sas/qtgui/Perspectives/Inversion/InversionLogic.py @ 318b353e

ESS_GUIESS_GUI_DocsESS_GUI_batch_fittingESS_GUI_bumps_abstractionESS_GUI_iss1116ESS_GUI_iss879ESS_GUI_iss959ESS_GUI_openclESS_GUI_orderingESS_GUI_sync_sascalc
Last change on this file since 318b353e was 318b353e, checked in by krzywon, 4 years ago

Automatically plot P(r) fits when calculated and small tweaks.

  • Property mode set to 100644
File size: 8.0 KB
Line 
1import math
2import pylab
3import numpy as np
4
5from sas.qtgui.Plotting.PlotterData import Data1D
6
7PR_FIT_LABEL = r"$P_{fit}(r)$"
8PR_LOADED_LABEL = r"$P_{loaded}(r)$"
9IQ_DATA_LABEL = r"$I_{obs}(q)$"
10IQ_FIT_LABEL = r"$I_{fit}(q)$"
11IQ_SMEARED_LABEL = r"$I_{smeared}(q)$"
12GROUP_ID_IQ_DATA = r"$I_{obs}(q)$"
13GROUP_ID_PR_FIT = r"$P_{fit}(r)$"
14
15
16class InversionLogic(object):
17    """
18    All the data-related logic. This class deals exclusively with Data1D/2D
19    No QStandardModelIndex here.
20    """
21
22    # TODO: Add way to change this value
23    _pr_n_pts = 51
24
25    def __init__(self, data=None):
26        self._data = data
27        self.data_is_loaded = False
28        if data is not None:
29            self.data_is_loaded = True
30
31    @property
32    def data(self):
33        return self._data
34
35    @data.setter
36    def data(self, value):
37        """ data setter """
38        self._data = value
39        self.data_is_loaded = (self._data is not None)
40
41    def isLoadedData(self):
42        """ accessor """
43        return self.data_is_loaded
44
45    def new1DPlot(self, out, pr, q=None):
46        """
47        Create a new 1D data instance based on fitting results
48        """
49
50        qtemp = pr.x
51        if q is not None:
52            qtemp = q
53
54        # Make a plot
55        maxq = max(qtemp)
56
57        minq = min(qtemp)
58
59        # Check for user min/max
60        if pr.q_min is not None and maxq >= pr.q_min >= minq:
61            minq = pr.q_min
62        if pr.q_max is not None and maxq >= pr.q_max >= minq:
63            maxq = pr.q_max
64
65        x = pylab.arange(minq, maxq, maxq / 301.0)
66        y = np.zeros(len(x))
67        err = np.zeros(len(x))
68        for i in range(len(x)):
69            value = pr.iq(out, x[i])
70            y[i] = value
71            try:
72                err[i] = math.sqrt(math.fabs(value))
73            except:
74                err[i] = 1.0
75                print(("Error getting error", value, x[i]))
76
77        new_plot = Data1D(x, y)
78        new_plot.name = IQ_FIT_LABEL
79        new_plot.xaxis("\\rm{Q}", 'A^{-1}')
80        new_plot.yaxis("\\rm{Intensity} ", "cm^{-1}")
81        title = "I(q)"
82        new_plot.title = title
83
84        # If we have a group ID, use it
85        if 'plot_group_id' in pr.info:
86            new_plot.group_id = pr.info["plot_group_id"]
87        new_plot.id = IQ_FIT_LABEL
88
89        # If we have used slit smearing, plot the smeared I(q) too
90        if pr.slit_width > 0 or pr.slit_height > 0:
91            x = pylab.arange(minq, maxq, maxq / 301.0)
92            y = np.zeros(len(x))
93            err = np.zeros(len(x))
94            for i in range(len(x)):
95                value = pr.iq_smeared(pr.out, x[i])
96                y[i] = value
97                try:
98                    err[i] = math.sqrt(math.fabs(value))
99                except:
100                    err[i] = 1.0
101                    print(("Error getting error", value, x[i]))
102
103            new_plot = Data1D(x, y)
104            new_plot.name = IQ_SMEARED_LABEL
105            new_plot.xaxis("\\rm{Q}", 'A^{-1}')
106            new_plot.yaxis("\\rm{Intensity} ", "cm^{-1}")
107            # If we have a group ID, use it
108            if 'plot_group_id' in pr.info:
109                new_plot.group_id = pr.info["plot_group_id"]
110            new_plot.id = IQ_SMEARED_LABEL
111            new_plot.title = title
112
113        return new_plot
114
115    def update1DPlot(self, plot, out, pr, q=None):
116        """
117        Create a new 1D data instance based on fitting results
118        """
119
120        qtemp = pr.x
121        if q is not None:
122            qtemp = q
123
124        # Make a plot
125        maxq = max(qtemp)
126
127        minq = min(qtemp)
128
129        # Check for user min/max
130        if pr.q_min is not None and maxq >= pr.q_min >= minq:
131            minq = pr.q_min
132        if pr.q_max is not None and maxq >= pr.q_max >= minq:
133            maxq = pr.q_max
134
135        x = pylab.arange(minq, maxq, maxq / 301.0)
136        y = np.zeros(len(x))
137        err = np.zeros(len(x))
138        for i in range(len(x)):
139            value = pr.iq(out, x[i])
140            y[i] = value
141            try:
142                err[i] = math.sqrt(math.fabs(value))
143            except:
144                err[i] = 1.0
145                print(("Error getting error", value, x[i]))
146
147        plot.x = x
148        plot.y = y
149
150        # If we have used slit smearing, plot the smeared I(q) too
151        if pr.slit_width > 0 or pr.slit_height > 0:
152            x = pylab.arange(minq, maxq, maxq / 301.0)
153            y = np.zeros(len(x))
154            err = np.zeros(len(x))
155            for i in range(len(x)):
156                value = pr.iq_smeared(pr.out, x[i])
157                y[i] = value
158                try:
159                    err[i] = math.sqrt(math.fabs(value))
160                except:
161                    err[i] = 1.0
162                    print(("Error getting error", value, x[i]))
163
164            plot.x = x
165            plot.y = y
166
167        return plot
168
169    def newPRPlot(self, out, pr, cov=None):
170        """
171        """
172        # Show P(r)
173        x = pylab.arange(0.0, pr.d_max, pr.d_max / self._pr_n_pts)
174
175        y = np.zeros(len(x))
176        dy = np.zeros(len(x))
177
178        total = 0.0
179        pmax = 0.0
180        cov2 = np.ascontiguousarray(cov)
181
182        for i in range(len(x)):
183            if cov2 is None:
184                value = pr.pr(out, x[i])
185            else:
186                (value, dy[i]) = pr.pr_err(out, cov2, x[i])
187            total += value * pr.d_max / len(x)
188
189            # keep track of the maximum P(r) value
190            if value > pmax:
191                pmax = value
192
193            y[i] = value
194
195        # if self._normalize_output == True:
196        #     y = y / total
197        #     dy = dy / total
198        # elif self._scale_output_unity == True:
199        #     y = y / pmax
200        #     dy = dy / pmax
201
202        if cov2 is None:
203            new_plot = Data1D(x, y)
204        else:
205            new_plot = Data1D(x, y, dy=dy)
206        new_plot.name = PR_FIT_LABEL
207        new_plot.xaxis("\\rm{r}", 'A')
208        new_plot.yaxis("\\rm{P(r)} ", "cm^{-3}")
209        new_plot.title = "P(r) fit"
210        new_plot.id = PR_FIT_LABEL
211        new_plot.scale = "linear"
212        new_plot.group_id = GROUP_ID_PR_FIT
213
214        return new_plot
215
216    def updatePRPlot(self, plot, out, pr, cov=None):
217        x = pylab.arange(0.0, pr.d_max, pr.d_max / self._pr_n_pts)
218
219        y = np.zeros(len(x))
220        dy = np.zeros(len(x))
221
222        total = 0.0
223        pmax = 0.0
224        cov2 = np.ascontiguousarray(cov)
225
226        for i in range(len(x)):
227            if cov2 is None:
228                value = pr.pr(out, x[i])
229            else:
230                (value, dy[i]) = pr.pr_err(out, cov2, x[i])
231            total += value * pr.d_max / len(x)
232
233            # keep track of the maximum P(r) value
234            if value > pmax:
235                pmax = value
236
237            y[i] = value
238
239        # if self._normalize_output == True:
240        #     y = y / total
241        #     dy = dy / total
242        # elif self._scale_output_unity == True:
243        #     y = y / pmax
244        #     dy = dy / pmax
245        plot.x = x
246        plot.y = y
247
248        if cov2 is not None:
249            plot.dy = dy
250
251        return plot
252
253    def computeDataRange(self):
254        """
255        Wrapper for calculating the data range based on local dataset
256        """
257        return self.computeRangeFromData(self.data)
258
259    def computeRangeFromData(self, data):
260        """
261        Compute the minimum and the maximum range of the data
262        return the npts contains in data
263        """
264        qmin, qmax = None, None
265        if isinstance(data, Data1D):
266            try:
267                qmin = min(data.x)
268                qmax = max(data.x)
269            except (ValueError, TypeError):
270                msg = "Unable to find min/max/length of \n data named %s" % \
271                            self.data.filename
272                raise ValueError(msg)
273
274        else:
275            qmin = 0
276            try:
277                x = max(np.fabs(data.xmin), np.fabs(data.xmax))
278                y = max(np.fabs(data.ymin), np.fabs(data.ymax))
279            except (ValueError, TypeError):
280                msg = "Unable to find min/max of \n data named %s" % \
281                            self.data.filename
282                raise ValueError(msg)
283            qmax = np.sqrt(x * x + y * y)
284        return qmin, qmax
Note: See TracBrowser for help on using the repository browser.