SasView / sasmodels

Package for calculation of small angle scattering models using OpenCL.
BSD 3-Clause "New" or "Revised" License
15 stars 27 forks source link

New plotting option iq #491

Closed Caddy-Jones closed 2 years ago

Caddy-Jones commented 2 years ago

Required for the new_plotting_option_IQ branch in sasView

butlerpd commented 2 years ago

After discussions with Lionel and at the fortnightly meeting, it was agreed that this is not ready to merge. In fact we should separate this into two separate issues. The first addressing the ability to fit data in rescaled plots (e.g. I*Q^4 vs Q^2) so that the fit line is also plotted in that scale. The second would be to add more weighting functions beyond I, root I etc to include things like Q and Q^2 etc. Alternatively change the weighting choice from a finite selection to a free form?

In the meantime we should close this request but keep the branch which has started the work.

gonzalezma commented 2 years ago

Concerning the plotting of the data, I think that the change scale menu does the job, e.g: IQ2 The only issues that I see are:

  1. The functionality is hidden in the contextual menu that only appears when clicking the right-mouse button in the figure, so perhaps many users are not aware? Perhaps we will need a graph options menu in the principal menu, but this is a matter of taste and there are good reasons to keep the principal menu as simple as possible.
  2. The scaling of the axis for the rescaled figure is not good (at least for y), so one needs to set the y scale manually (either calling Set Graph Range or the graph navigation menu and then using the zoom, available always in the contextual menu).
butlerpd commented 2 years ago

Interesting @gonzalezma. I'm guessing you did the fit normally and then did change scale? Maybe @Caddy-Jones was saying that if you use change scale on the data and then do the fit that the fit line did not respect the new axis?

That may in fact be the issue now that I think about it because SasView probably does not know about the "scale" of the existing plot when it sends a new date set to an existing graph? Also I think that is what Lionel is trying to achieve: do fitting, including just changing parameters before doing a full fit, and see how closely the curve matches the data (in the new scale)? I suppose we may need to understand the use case better?

But it sounds like at least a good chunk of what is needed is already there?

pkienzle commented 1 year ago

I'm deleting stale branches, including this one. The PR will still have "restore branch" if you need it. In case I'm wrong, here's the patch:

diff --git a/sasmodels/sasview_model.py b/sasmodels/sasview_model.py
index ff122fa9..71c78f27 100644
--- a/sasmodels/sasview_model.py
+++ b/sasmodels/sasview_model.py
@@ -9,6 +9,8 @@ create a sasview model class to run that kernel as follows::
 """
 from __future__ import print_function

+import ast
+import operator as op
 import math
 from copy import deepcopy
 import collections
@@ -414,6 +416,7 @@ class SasviewModel(object):
         self.params = collections.OrderedDict()
         self.dispersion = collections.OrderedDict()
         self.details = {}
+        self.expression = "I"
         for p in self._model_info.parameters.user_parameters({}, is2d=True):
             if p.name in hidden:
                 continue
@@ -430,6 +433,14 @@ class SasviewModel(object):
                     'type': 'gaussian',
                 }

+    def set_expression(self, expression):
+        """Define the expression attribute"""
+        self.expression = expression
+
+    def get_expression(self):
+        """Return the expression attribute"""
+        return self.expression
+
     def __get_state__(self):
         # type: () -> Dict[str, Any]
         state = self.__dict__.copy()
@@ -740,14 +751,62 @@ class SasviewModel(object):
                             magnetic=is_magnetic)
         lazy_results = getattr(calculator, 'results',
                                lambda: collections.OrderedDict())
-        #print("result", result)

         calculator.release()
         #self._model.release()

+        # Calculate the expression string like a line of code
+        try:
+            result = self.safeEval(ast.parse(self.expression, mode='eval').body, qx, result)
+        except TypeError as e:
+            logging.error(f"Unknown symbol used in Y-axis Data, please only use I, Q and python mathematical operators"
+                          f"\n{e}")
+
         return result, lazy_results

+    def safeEval(self, node, q_vals, I_vals):
+        """ Calculate a string as a line of code
+
+        A safe form of the eval command, only calculates mathematical operators, can accept I, Q, q as variables.
+
+        :param node: expression to calculate
+        :type node: string
+        :param q_vals: q values
+        :type q_vals: np.array
+        :param I_vals: Calculated I values
+        :type I_vals: np.array
+        :return: The calculation outlined in the node string
+        :rtype: np.array
+        :raises TypeError: Exception isn't a variable named I, Q, q; or a mathematical operators.
+        """
+
+        # Supported operators
+        operators = {ast.Add: op.add, ast.Sub: op.sub, ast.Mult: op.mul, ast.Div: op.truediv, ast.Pow: op.pow,
+                     ast.BitXor: op.xor, ast.USub: op.neg}
+
+        # If node is just a number
+        if isinstance(node, ast.Num):
+            return node.n
+        # If mathematical operator eg, <num> <operator> <num>
+        elif isinstance(node, ast.BinOp):
+            return operators[type(node.op)](self.safeEval(node.left, q_vals, I_vals),
+                                            self.safeEval(node.right, q_vals, I_vals))
+        # Special case eg -1
+        elif isinstance(node, ast.UnaryOp):
+            return operators[type(node.op)](self.safeEval(node.operand, q_vals, I_vals))
+        # Accepted variable name
+        elif isinstance(node, ast.Name):
+            if node.id == "I":
+                return I_vals
+            elif node.id == "Q":
+                return q_vals
+            elif node.id == "q":
+                return q_vals
+        else:
+            raise TypeError(node)
+
+
     def calculate_ER(self, mode=1):
         # type: (int) -> float
         """