source: sasview/src/sas/qtgui/Perspectives/Fitting/FittingUtilities.py @ b5cc06e

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 b5cc06e was b5cc06e, checked in by Piotr Rozyczko <rozyczko@…>, 4 years ago

iqxy_parameters are no longer available

  • Property mode set to 100644
File size: 15.6 KB
Line 
1import copy
2
3from PyQt5 import QtCore
4from PyQt5 import QtGui
5from PyQt5 import QtWidgets
6
7import numpy
8
9from sas.qtgui.Plotting.PlotterData import Data1D
10from sas.qtgui.Plotting.PlotterData import Data2D
11
12model_header_captions = ['Parameter', 'Value', 'Min', 'Max', 'Units']
13
14model_header_tooltips = ['Select parameter for fitting',
15                         'Enter parameter value',
16                         'Enter minimum value for parameter',
17                         'Enter maximum value for parameter',
18                         'Unit of the parameter']
19
20poly_header_captions = ['Parameter', 'PD[ratio]', 'Min', 'Max', 'Npts', 'Nsigs',
21                        'Function', 'Filename']
22
23poly_header_tooltips = ['Select parameter for fitting',
24                        'Enter polydispersity ratio (STD/mean). '
25                        'STD: standard deviation from the mean value',
26                        'Enter minimum value for parameter',
27                        'Enter maximum value for parameter',
28                        'Enter number of points for parameter',
29                        'Enter number of sigmas parameter',
30                        'Select distribution function',
31                        'Select filename with user-definable distribution']
32
33error_tooltip = 'Error value for fitted parameter'
34header_error_caption = 'Error'
35
36def replaceShellName(param_name, value):
37    """
38    Updates parameter name from <param_name>[n_shell] to <param_name>value
39    """
40    assert '[' in param_name
41    return param_name[:param_name.index('[')]+str(value)
42
43def getIterParams(model):
44    """
45    Returns a list of all multi-shell parameters in 'model'
46    """
47    return list([par for par in model.iq_parameters if "[" in par.name])
48
49def getMultiplicity(model):
50    """
51    Finds out if 'model' has multishell parameters.
52    If so, returns the name of the counter parameter and the number of shells
53    """
54    iter_params = getIterParams(model)
55    param_name = ""
56    param_length = 0
57    if iter_params:
58        param_length = iter_params[0].length
59        param_name = iter_params[0].length_control
60        if param_name is None and '[' in iter_params[0].name:
61            param_name = iter_params[0].name[:iter_params[0].name.index('[')]
62    return (param_name, param_length)
63
64def addParametersToModel(parameters, kernel_module, is2D):
65    """
66    Update local ModelModel with sasmodel parameters
67    """
68    multishell_parameters = getIterParams(parameters)
69    multishell_param_name, _ = getMultiplicity(parameters)
70    if is2D:
71        params = [p for p in parameters.kernel_parameters if p.type != 'magnetic']
72    else:
73        params = parameters.iq_parameters
74    item = []
75    for param in params:
76        # don't include shell parameters
77        if param.name == multishell_param_name:
78            continue
79        # Modify parameter name from <param>[n] to <param>1
80        item_name = param.name
81        if param in multishell_parameters:
82            continue
83        #    item_name = replaceShellName(param.name, 1)
84
85        item1 = QtGui.QStandardItem(item_name)
86        item1.setCheckable(True)
87        item1.setEditable(False)
88        # item_err = QtGui.QStandardItem()
89        # check for polydisp params
90        if param.polydisperse:
91            poly_item = QtGui.QStandardItem("Polydispersity")
92            poly_item.setEditable(False)
93            item1_1 = QtGui.QStandardItem("Distribution")
94            item1_1.setEditable(False)
95            # Find param in volume_params
96            for p in parameters.form_volume_parameters:
97                if p.name != param.name:
98                    continue
99                width = kernel_module.getParam(p.name+'.width')
100                type = kernel_module.getParam(p.name+'.type')
101
102                item1_2 = QtGui.QStandardItem(str(width))
103                item1_2.setEditable(False)
104                item1_3 = QtGui.QStandardItem()
105                item1_3.setEditable(False)
106                item1_4 = QtGui.QStandardItem()
107                item1_4.setEditable(False)
108                item1_5 = QtGui.QStandardItem(type)
109                item1_5.setEditable(False)
110                poly_item.appendRow([item1_1, item1_2, item1_3, item1_4, item1_5])
111                break
112            # Add the polydisp item as a child
113            item1.appendRow([poly_item])
114        # Param values
115        item2 = QtGui.QStandardItem(str(param.default))
116        # TODO: the error column.
117        # Either add a proxy model or a custom view delegate
118        #item_err = QtGui.QStandardItem()
119        item3 = QtGui.QStandardItem(str(param.limits[0]))
120        item4 = QtGui.QStandardItem(str(param.limits[1]))
121        item5 = QtGui.QStandardItem(param.units)
122        item5.setEditable(False)
123        item.append([item1, item2, item3, item4, item5])
124    return item
125
126def addSimpleParametersToModel(parameters, is2D):
127    """
128    Update local ModelModel with sasmodel parameters
129    """
130    if is2D:
131        params = [p for p in parameters.kernel_parameters if p.type != 'magnetic']
132    else:
133        params = parameters.iq_parameters
134    item = []
135    for param in params:
136        # Create the top level, checkable item
137        item_name = param.name
138        item1 = QtGui.QStandardItem(item_name)
139        item1.setCheckable(True)
140        item1.setEditable(False)
141        # Param values
142        # TODO: add delegate for validation of cells
143        item2 = QtGui.QStandardItem(str(param.default))
144        item4 = QtGui.QStandardItem(str(param.limits[0]))
145        item5 = QtGui.QStandardItem(str(param.limits[1]))
146        item6 = QtGui.QStandardItem(param.units)
147        item6.setEditable(False)
148        item.append([item1, item2, item4, item5, item6])
149    return item
150
151def addCheckedListToModel(model, param_list):
152    """
153    Add a QItem to model. Makes the QItem checkable
154    """
155    assert isinstance(model, QtGui.QStandardItemModel)
156    item_list = [QtGui.QStandardItem(item) for item in param_list]
157    item_list[0].setCheckable(True)
158    model.appendRow(item_list)
159
160def addHeadersToModel(model):
161    """
162    Adds predefined headers to the model
163    """
164    for i, item in enumerate(model_header_captions):
165        model.setHeaderData(i, QtCore.Qt.Horizontal, item)
166
167    model.header_tooltips = copy.copy(model_header_tooltips)
168
169def addErrorHeadersToModel(model):
170    """
171    Adds predefined headers to the model
172    """
173    model_header_error_captions = copy.copy(model_header_captions)
174    model_header_error_captions.insert(2, header_error_caption)
175    for i, item in enumerate(model_header_error_captions):
176        model.setHeaderData(i, QtCore.Qt.Horizontal, item)
177
178    model_header_error_tooltips = copy.copy(model_header_tooltips)
179    model_header_error_tooltips.insert(2, error_tooltip)
180    model.header_tooltips = copy.copy(model_header_error_tooltips)
181
182def addPolyHeadersToModel(model):
183    """
184    Adds predefined headers to the model
185    """
186    for i, item in enumerate(poly_header_captions):
187        model.setHeaderData(i, QtCore.Qt.Horizontal, item)
188
189    model.header_tooltips = copy.copy(poly_header_tooltips)
190
191
192def addErrorPolyHeadersToModel(model):
193    """
194    Adds predefined headers to the model
195    """
196    poly_header_error_captions = copy.copy(poly_header_captions)
197    poly_header_error_captions.insert(2, header_error_caption)
198    for i, item in enumerate(poly_header_error_captions):
199        model.setHeaderData(i, QtCore.Qt.Horizontal, item)
200
201    poly_header_error_tooltips = copy.copy(poly_header_tooltips)
202    poly_header_error_tooltips.insert(2, error_tooltip)
203    model.header_tooltips = copy.copy(poly_header_error_tooltips)
204
205def addShellsToModel(parameters, model, index):
206    """
207    Find out multishell parameters and update the model with the requested number of them
208    """
209    multishell_parameters = getIterParams(parameters)
210
211    for i in range(index):
212        for par in multishell_parameters:
213            # Create the name: <param>[<i>], e.g. "sld1" for parameter "sld[n]"
214            param_name = replaceShellName(par.name, i+1)
215            item1 = QtGui.QStandardItem(param_name)
216            item1.setCheckable(True)
217            # check for polydisp params
218            if par.polydisperse:
219                poly_item = QtGui.QStandardItem("Polydispersity")
220                item1_1 = QtGui.QStandardItem("Distribution")
221                # Find param in volume_params
222                for p in parameters.form_volume_parameters:
223                    if p.name != par.name:
224                        continue
225                    item1_2 = QtGui.QStandardItem(str(p.default))
226                    item1_3 = QtGui.QStandardItem(str(p.limits[0]))
227                    item1_4 = QtGui.QStandardItem(str(p.limits[1]))
228                    item1_5 = QtGui.QStandardItem(p.units)
229                    poly_item.appendRow([item1_1, item1_2, item1_3, item1_4, item1_5])
230                    break
231                item1.appendRow([poly_item])
232
233            item2 = QtGui.QStandardItem(str(par.default))
234            item3 = QtGui.QStandardItem(str(par.limits[0]))
235            item4 = QtGui.QStandardItem(str(par.limits[1]))
236            item5 = QtGui.QStandardItem(par.units)
237            model.appendRow([item1, item2, item3, item4, item5])
238
239def calculateChi2(reference_data, current_data):
240    """
241    Calculate Chi2 value between two sets of data
242    """
243
244    # WEIGHING INPUT
245    #from sas.sasgui.perspectives.fitting.utils import get_weight
246    #flag = self.get_weight_flag()
247    #weight = get_weight(data=self.data, is2d=self._is_2D(), flag=flag)
248    chisqr = None
249    if reference_data is None:
250        return chisqr
251
252    # temporary default values for index and weight
253    index = None
254    weight = None
255
256    # Get data: data I, theory I, and data dI in order
257    if isinstance(reference_data, Data2D):
258        if index is None:
259            index = numpy.ones(len(current_data.data), dtype=bool)
260        if weight is not None:
261            current_data.err_data = weight
262        # get rid of zero error points
263        index = index & (current_data.err_data != 0)
264        index = index & (numpy.isfinite(current_data.data))
265        fn = current_data.data[index]
266        gn = reference_data.data[index]
267        en = current_data.err_data[index]
268    else:
269        # 1 d theory from model_thread is only in the range of index
270        if index is None:
271            index = numpy.ones(len(current_data.y), dtype=bool)
272        if weight is not None:
273            current_data.dy = weight
274        if current_data.dy is None or current_data.dy == []:
275            dy = numpy.ones(len(current_data.y))
276        else:
277            ## Set consistently w/AbstractFitengine:
278            # But this should be corrected later.
279            dy = copy.deepcopy(current_data.dy)
280            dy[dy == 0] = 1
281        fn = current_data.y[index]
282        gn = reference_data.y
283        en = dy[index]
284    # Calculate the residual
285    try:
286        res = (fn - gn) / en
287    except ValueError:
288        #print "Chi2 calculations: Unmatched lengths %s, %s, %s" % (len(fn), len(gn), len(en))
289        return None
290
291    residuals = res[numpy.isfinite(res)]
292    chisqr = numpy.average(residuals * residuals)
293
294    return chisqr
295
296def residualsData1D(reference_data, current_data):
297    """
298    Calculate the residuals for difference of two Data1D sets
299    """
300    # temporary default values for index and weight
301    index = None
302    weight = None
303
304    # 1d theory from model_thread is only in the range of index
305    if current_data.dy is None or current_data.dy == []:
306        dy = numpy.ones(len(current_data.y))
307    else:
308        dy = weight if weight is not None else numpy.ones(len(current_data.y))
309        dy[dy == 0] = 1
310    fn = current_data.y[index][0]
311    gn = reference_data.y
312    en = dy[index][0]
313    # build residuals
314    residuals = Data1D()
315    if len(fn) == len(gn):
316        y = (fn - gn)/en
317        residuals.y = -y
318    else:
319        # TODO: fix case where applying new data from file on top of existing model data
320        try:
321            y = (fn - gn[index][0]) / en
322            residuals.y = y
323        except ValueError:
324            # value errors may show up every once in a while for malformed columns,
325            # just reuse what's there already
326            pass
327
328    residuals.x = current_data.x[index][0]
329    residuals.dy = numpy.ones(len(residuals.y))
330    residuals.dx = None
331    residuals.dxl = None
332    residuals.dxw = None
333    residuals.ytransform = 'y'
334    # For latter scale changes
335    residuals.xaxis('\\rm{Q} ', 'A^{-1}')
336    residuals.yaxis('\\rm{Residuals} ', 'normalized')
337
338    return residuals
339
340def residualsData2D(reference_data, current_data):
341    """
342    Calculate the residuals for difference of two Data2D sets
343    """
344    # temporary default values for index and weight
345    # index = None
346    weight = None
347
348    # build residuals
349    residuals = Data2D()
350    # Not for trunk the line below, instead use the line above
351    current_data.clone_without_data(len(current_data.data), residuals)
352    residuals.data = None
353    fn = current_data.data
354    gn = reference_data.data
355    en = current_data.err_data if weight is None else weight
356    residuals.data = (fn - gn) / en
357    residuals.qx_data = current_data.qx_data
358    residuals.qy_data = current_data.qy_data
359    residuals.q_data = current_data.q_data
360    residuals.err_data = numpy.ones(len(residuals.data))
361    residuals.xmin = min(residuals.qx_data)
362    residuals.xmax = max(residuals.qx_data)
363    residuals.ymin = min(residuals.qy_data)
364    residuals.ymax = max(residuals.qy_data)
365    residuals.q_data = current_data.q_data
366    residuals.mask = current_data.mask
367    residuals.scale = 'linear'
368    # check the lengths
369    if len(residuals.data) != len(residuals.q_data):
370        return None
371    return residuals
372
373def plotResiduals(reference_data, current_data):
374    """
375    Create Data1D/Data2D with residuals, ready for plotting
376    """
377    data_copy = copy.deepcopy(current_data)
378    # Get data: data I, theory I, and data dI in order
379    method_name = current_data.__class__.__name__
380    residuals_dict = {"Data1D": residualsData1D,
381                      "Data2D": residualsData2D}
382
383    residuals = residuals_dict[method_name](reference_data, data_copy)
384
385    theory_name = str(current_data.name.split()[0])
386    residuals.name = "Residuals for " + str(theory_name) + "[" + \
387                    str(reference_data.filename) + "]"
388    residuals.title = residuals.name
389    residuals.ytransform = 'y'
390
391    # when 2 data have the same id override the 1 st plotted
392    # include the last part if keeping charts for separate models is required
393    residuals.id = "res" + str(reference_data.id) # + str(theory_name)
394    # group_id specify on which panel to plot this data
395    group_id = reference_data.group_id
396    residuals.group_id = "res" + str(group_id)
397
398    # Symbol
399    residuals.symbol = 0
400    residuals.hide_error = False
401
402    return residuals
403
404def binary_encode(i, digits):
405    return [i >> d & 1 for d in range(digits)]
406
407def getWeight(data, is2d, flag=None):
408    """
409    Received flag and compute error on data.
410    :param flag: flag to transform error of data.
411    """
412    weight = None
413    if is2d:
414        dy_data = data.err_data
415        data = data.data
416    else:
417        dy_data = data.dy
418        data = data.y
419
420    if flag == 0:
421        weight = numpy.ones_like(data)
422    elif flag == 1:
423        weight = dy_data
424    elif flag == 2:
425        weight = numpy.sqrt(numpy.abs(data))
426    elif flag == 3:
427        weight = numpy.abs(data)
428    return weight
429
430def updateKernelWithResults(kernel, results):
431    """
432    Takes model kernel and applies results dict to its parameters,
433    returning the modified (deep) copy of the kernel.
434    """
435    assert(isinstance(results, dict))
436    local_kernel = copy.deepcopy(kernel)
437
438    for parameter in results.keys():
439        # Update the parameter value - note: this supports +/-inf as well
440        local_kernel.setParam(parameter, results[parameter][0])
441
442    return local_kernel
443
444
Note: See TracBrowser for help on using the repository browser.