Fix estimate_parameters for binned and non uniform axis

hyperspy/_components/doniach.py

@@ -19,6 +19,7 @@
 import math
 import numpy as np
 
+from hyperspy.component import _get_scaling_factor
 from hyperspy._components.expression import Expression
 from hyperspy._components.gaussian import _estimate_gaussian_parameters
 from hyperspy.misc.utils import is_binned # remove in v2.0
@@ -145,8 +146,8 @@ def estimate_parameters(self, signal, x1, x2, only_current=False):
         axis = signal.axes_manager.signal_axes[0]
         centre, height, sigma = _estimate_gaussian_parameters(signal, x1, x2,
                                                               only_current)
-        scaling_factor = axis.scale if axis.is_uniform \
-                         else np.gradient(axis.axis)[axis.value2index(centre)]
+        scaling_factor = _get_scaling_factor(signal, axis, centre)
+
         if only_current is True:
             self.centre.value = centre
             self.sigma.value = sigma

hyperspy/_components/exponential.py

@@ -52,7 +52,7 @@ class Exponential(Expression):
     """
 
     def __init__(self, A=1., tau=1., module="numexpr", **kwargs):
-        super(Exponential, self).__init__(
+        super().__init__(
             expression="A * exp(-x / tau)",
             name="Exponential",
             A=A,
@@ -86,7 +86,7 @@ def estimate_parameters(self, signal, x1, x2, only_current=False):
         bool
 
         """
-        super(Exponential, self)._estimate_parameters(signal)
+        super()._estimate_parameters(signal)
         axis = signal.axes_manager.signal_axes[0]
         i1, i2 = axis.value_range_to_indices(x1, x2)
         if i1 + 1 == i2:

hyperspy/_components/gaussian.py

@@ -21,6 +21,7 @@
 import numpy as np
 import dask.array as da
 
+from hyperspy.component import _get_scaling_factor
 from hyperspy._components.expression import Expression
 from hyperspy.misc.utils import is_binned # remove in v2.0
 
@@ -160,12 +161,12 @@ def estimate_parameters(self, signal, x1, x2, only_current=False):
         >>> g.estimate_parameters(s, -10, 10, False)
         """
 
-        super(Gaussian, self)._estimate_parameters(signal)
+        super()._estimate_parameters(signal)
         axis = signal.axes_manager.signal_axes[0]
         centre, height, sigma = _estimate_gaussian_parameters(signal, x1, x2,
                                                               only_current)
-        scaling_factor = axis.scale if axis.is_uniform \
-                         else np.gradient(axis.axis)[axis.value2index(centre)]
+        scaling_factor = _get_scaling_factor(signal, axis, centre)
+
         if only_current is True:
             self.centre.value = centre
             self.sigma.value = sigma

hyperspy/_components/gaussianhf.py

@@ -17,10 +17,10 @@
 # along with  HyperSpy.  If not, see <http://www.gnu.org/licenses/>.
 
 import math
-import numpy as np
 
 from hyperspy._components.expression import Expression
 from hyperspy._components.gaussian import _estimate_gaussian_parameters
+from hyperspy.component import _get_scaling_factor
 from hyperspy.misc.utils import is_binned # remove in v2.0
 
 sqrt2pi = math.sqrt(2 * math.pi)
@@ -136,12 +136,12 @@ def estimate_parameters(self, signal, x1, x2, only_current=False):
         >>> g.estimate_parameters(s, -10, 10, False)
         """
 
-        super(GaussianHF, self)._estimate_parameters(signal)
+        super()._estimate_parameters(signal)
         axis = signal.axes_manager.signal_axes[0]
         centre, height, sigma = _estimate_gaussian_parameters(signal, x1, x2,
                                                               only_current)
-        scaling_factor = axis.scale if axis.is_uniform \
-                         else np.gradient(axis.axis)[axis.value2index(centre)]
+        scaling_factor = _get_scaling_factor(signal, axis, centre)
+
         if only_current is True:
             self.centre.value = centre
             self.fwhm.value = sigma * sigma2fwhm

hyperspy/_components/lorentzian.py

@@ -19,6 +19,7 @@
 import numpy as np
 import dask.array as da
 
+from hyperspy.component import _get_scaling_factor
 from hyperspy._components.expression import Expression
 from hyperspy.misc.utils import is_binned # remove in v2.0
 
@@ -160,12 +161,12 @@ def estimate_parameters(self, signal, x1, x2, only_current=False):
         >>> g.estimate_parameters(s, -10, 10, False)
         """
 
-        super(Lorentzian, self)._estimate_parameters(signal)
+        super()._estimate_parameters(signal)
         axis = signal.axes_manager.signal_axes[0]
         centre, height, gamma = _estimate_lorentzian_parameters(signal, x1, x2,
                                                               only_current)
-        scaling_factor = axis.scale if axis.is_uniform \
-                         else np.gradient(axis.axis)[axis.value2index(centre)]
+        scaling_factor = _get_scaling_factor(signal, axis, centre)
+
         if only_current is True:
             self.centre.value = centre
             self.gamma.value = gamma

hyperspy/_components/offset.py

@@ -29,19 +29,19 @@ class Offset(Component):
     r"""Component to add a constant value in the y-axis.
 
     .. math::
-    
+
         f(x) = k
 
     ============ =============
-    Variable      Parameter 
+    Variable      Parameter
     ============ =============
-    :math:`k`     offset   
+    :math:`k`     offset
     ============ =============
 
     Parameters
     -----------
-    offset : float 
-        
+    offset : float
+
     """
 
     def __init__(self, offset=0.):
@@ -86,13 +86,17 @@ def estimate_parameters(self, signal, x1, x2, only_current=False):
         bool
 
         """
-        super(Offset, self)._estimate_parameters(signal)
+        super()._estimate_parameters(signal)
         axis = signal.axes_manager.signal_axes[0]
         i1, i2 = axis.value_range_to_indices(x1, x2)
-        # using the mean of the gradient for non-uniform axes is a best guess
-        # to the scaling of binned signals for the estimation
-        scaling_factor = axis.scale if axis.is_uniform \
-                         else np.mean(np.gradient(axis.axis))
+        if is_binned(signal):
+        # in v2 replace by
+        #if axis.is_binned:
+            # using the mean of the gradient for non-uniform axes is a best
+            # guess to the scaling of binned signals for the estimation
+            scaling_factor = axis.scale if axis.is_uniform \
+                             else np.mean(np.gradient(axis.axis), axis=-1)
+
         if only_current is True:
             self.offset.value = signal()[i1:i2].mean()
             if is_binned(signal):

hyperspy/_components/pes_voigt.py

@@ -19,7 +19,7 @@
 import numpy as np
 import math
 
-from hyperspy.component import Component
+from hyperspy.component import Component, _get_scaling_factor
 from hyperspy._components.gaussian import _estimate_gaussian_parameters
 from hyperspy.misc.utils import is_binned # remove in v2.0
 
@@ -292,12 +292,12 @@ def estimate_parameters(self, signal, E1, E2, only_current=False):
         >>> g.estimate_parameters(s, -10, 10, False)
 
         """
-        super(PESVoigt, self)._estimate_parameters(signal)
+        super()._estimate_parameters(signal)
         axis = signal.axes_manager.signal_axes[0]
         centre, height, sigma = _estimate_gaussian_parameters(signal, E1, E2,
                                                               only_current)
-        scaling_factor = axis.scale if axis.is_uniform \
-                         else np.gradient(axis.axis)[axis.value2index(centre)]
+        scaling_factor = _get_scaling_factor(signal, axis, centre)
+
         if only_current is True:
             self.centre.value = centre
             self.FWHM.value = sigma * sigma2fwhm

hyperspy/_components/polynomial.py

@@ -31,7 +31,7 @@ class Polynomial(Expression):
     """n-order polynomial component.
 
     Polynomial component consisting of order + 1 parameters.
-    The parameters are named "a" followed by the corresponding order, 
+    The parameters are named "a" followed by the corresponding order,
     i.e.
 
     .. math::
@@ -58,10 +58,10 @@ def __init__(self, order=2, module="numexpr", **kwargs):
             raise ValueError("Polynomial of order 0 is not supported.")
         coeff_list = ['{}'.format(o).zfill(len(list(str(order)))) for o in
                       range(order, -1, -1)]
-        expr = "+".join(["a{}*x**{}".format(c, o) for c, o in 
+        expr = "+".join(["a{}*x**{}".format(c, o) for c, o in
                          zip(coeff_list, range(order, -1, -1))])
         name = "{} order Polynomial".format(ordinal(order))
-        super().__init__(expression=expr, name=name, module=module, 
+        super().__init__(expression=expr, name=name, module=module,
              autodoc=False, **kwargs)
         self._id_name = "eab91275-88db-4855-917a-cdcbe7209592"
 
@@ -92,10 +92,15 @@ def estimate_parameters(self, signal, x1, x2, only_current=False):
         super()._estimate_parameters(signal)
         axis = signal.axes_manager.signal_axes[0]
         i1, i2 = axis.value_range_to_indices(x1, x2)
-        # using the mean of the gradient for non-uniform axes is a best guess
-        # to the scaling of binned signals for the estimation
-        scaling_factor = axis.scale if axis.is_uniform \
-                         else np.mean(np.gradient(axis.axis))
+
+        if is_binned(signal):
+        # in v2 replace by
+        #if axis.is_binned:
+            # using the mean of the gradient for non-uniform axes is a best
+            # guess to the scaling of binned signals for the estimation
+            scaling_factor = axis.scale if axis.is_uniform \
+                             else np.mean(np.gradient(axis.axis), axis=-1)
+
         if only_current is True:
             estimation = np.polyfit(axis.axis[i1:i2],
                                     signal()[i1:i2],
@@ -147,7 +152,7 @@ def convert_to_polynomial(poly_dict):
     """
     _logger.info("Converting the polynomial to the new definition")
     poly_order = poly_dict['order']
-    coeff_list = ['{}'.format(o).zfill(len(list(str(poly_dict['order'])))) 
+    coeff_list = ['{}'.format(o).zfill(len(list(str(poly_dict['order']))))
                   for o in range(poly_dict['order'], -1, -1)]
     poly2_dict = dict(poly_dict)
     coefficient_dict = poly_dict['parameters'][0]

hyperspy/_components/polynomial_deprecated.py

@@ -126,13 +126,18 @@ def estimate_parameters(self, signal, x1, x2, only_current=False):
         -------
         bool
         """
-        super(Polynomial, self)._estimate_parameters(signal)
+        super()._estimate_parameters(signal)
         axis = signal.axes_manager.signal_axes[0]
         i1, i2 = axis.value_range_to_indices(x1, x2)
-        # using the mean of the gradient for non-uniform axes is a best guess
-        # to the scaling of binned signals for the estimation
-        scaling_factor = axis.scale if axis.is_uniform \
-                         else np.mean(np.gradient(axis.axis))
+
+        if is_binned(signal):
+        # in v2 replace by
+        #if axis.is_binned:
+            # using the mean of the gradient for non-uniform axes is a best
+            # guess to the scaling of binned signals for the estimation
+            scaling_factor = axis.scale if axis.is_uniform \
+                             else np.mean(np.gradient(axis.axis), axis=-1)
+
         if only_current is True:
             estimation = np.polyfit(axis.axis[i1:i2],
                                     signal()[i1:i2],

hyperspy/_components/power_law.py

@@ -111,7 +111,7 @@ def estimate_parameters(self, signal, x1, x2, only_current=False,
         {bool, tuple of values}
 
         """
-        super(PowerLaw, self)._estimate_parameters(signal)
+        super()._estimate_parameters(signal)
         axis = signal.axes_manager.signal_axes[0]
         i1, i2 = axis.value_range_to_indices(x1, x2)
         if not (i2 + i1) % 2 == 0:

hyperspy/_components/skew_normal.py

@@ -20,6 +20,7 @@
 import dask.array as da
 import sympy
 
+from hyperspy.component import _get_scaling_factor
 from hyperspy._components.expression import Expression
 from hyperspy.misc.utils import is_binned # remove in v2.0
 
@@ -205,12 +206,13 @@ def estimate_parameters(self, signal, x1, x2, only_current=False):
         >>> g.estimate_parameters(s, -10, 10, False)
         """
 
-        super(SkewNormal, self)._estimate_parameters(signal)
+        super()._estimate_parameters(signal)
         axis = signal.axes_manager.signal_axes[0]
-        x0, height, scale, shape = _estimate_skewnormal_parameters(signal, x1,
-                                                                   x2, only_current)
-        scaling_factor = axis.scale if axis.is_uniform \
-                         else np.gradient(axis.axis)[axis.value2index(x0)]
+        x0, height, scale, shape = _estimate_skewnormal_parameters(
+            signal, x1, x2, only_current
+            )
+        scaling_factor = _get_scaling_factor(signal, axis, x0)
+
         if only_current is True:
             self.x0.value = x0
             self.A.value = height * sqrt2pi

hyperspy/_components/split_voigt.py

@@ -18,7 +18,7 @@
 
 import numpy as np
 
-from hyperspy.component import Component
+from hyperspy.component import Component, _get_scaling_factor
 from hyperspy._components.gaussian import _estimate_gaussian_parameters
 from hyperspy.docstrings.parameters import FUNCTION_ND_DOCSTRING
 from hyperspy.misc.utils import is_binned # remove in v2.0
@@ -190,12 +190,12 @@ def estimate_parameters(self, signal, x1, x2, only_current=False):
         >>> g.estimate_parameters(s, -10,10, False)
 
         """
-        super(SplitVoigt, self)._estimate_parameters(signal)
+        super()._estimate_parameters(signal)
         axis = signal.axes_manager.signal_axes[0]
         centre, height, sigma = _estimate_gaussian_parameters(signal, x1, x2,
                                                               only_current)
-        scaling_factor = axis.scale if axis.is_uniform \
-                         else np.gradient(axis.axis)[axis.value2index(centre)]
+        scaling_factor = _get_scaling_factor(signal, axis, centre)
+
         if only_current is True:
             self.centre.value = centre
             self.sigma1.value = sigma

hyperspy/_components/voigt.py

@@ -20,6 +20,7 @@
 import sympy
 import numpy as np
 
+from hyperspy.component import _get_scaling_factor
 from hyperspy._components.expression import Expression
 from hyperspy._components.gaussian import _estimate_gaussian_parameters
 from hyperspy.misc.utils import is_binned # remove in v2.0
@@ -157,12 +158,12 @@ def estimate_parameters(self, signal, x1, x2, only_current=False):
         >>> g.estimate_parameters(s, -10, 10, False)
 
         """
-        super(Voigt, self)._estimate_parameters(signal)
+        super()._estimate_parameters(signal)
         axis = signal.axes_manager.signal_axes[0]
         centre, height, sigma = _estimate_gaussian_parameters(signal, x1, x2,
                                                               only_current)
-        scaling_factor = axis.scale if axis.is_uniform \
-                         else np.gradient(axis.axis)[axis.value2index(centre)]
+        scaling_factor = _get_scaling_factor(signal, axis, centre)
+
         if only_current is True:
             self.centre.value = centre
             self.sigma.value = sigma

hyperspy/axes.py

@@ -559,22 +559,33 @@ def value2index(self, value, rounding=round):
         ------
         ValueError
             If any value is out of the axis limits.
-
+            If value is NaN
+            If value is an array containing at least 1 NaN value
+            If value is an array containing at least 1 out-of-bound value
         """
         if value is None:
             return None
-        if self.low_value <= value <= self.high_value:
-            return (np.abs(self.axis - value)).argmin()
+        #Should evaluate on both arrays and scalars. Raises error if there are
+        #nan values in array
+        if np.all((value >= self.low_value)*(value <= self.high_value)):
+            #initialise the index same dimension as input, force type to int
+            index = np.empty_like(value,dtype=int)
+            #assign on flat, iterate on flat.
+            for i,v in enumerate(np.asarray(value).flat):
+                index.flat[i] = (np.abs(self.axis - v)).argmin()
+            #Squeezing to get a scalar out if scalar in. See squeeze doc
+            return np.squeeze(index)[()]
         else:
-            index = int(value)
-            if self.size > index >= 0:
-                return index
-            else:
-                raise ValueError(
-                    f'The value {value} is out of the limits '
-                    f'[{self.low_value:.3g}-{self.high_value:.3g}] of the '
-                    f'"{self._get_name()}" axis.'
-                    )
+            ### I commented out this because it seemed dangerous
+            # index = int(value)
+            # if self.size > index >= 0:
+            #     return index
+            # else:
+            raise ValueError(
+                f'The value {value} is out of the limits '
+                f'[{self.low_value:.3g}-{self.high_value:.3g}] of the '
+                f'"{self._get_name()}" axis.'
+                )
 
     def index2value(self, index):
         if isinstance(index, da.Array):
@@ -623,7 +634,6 @@ def value_range_to_indices(self, v1, v2):
                 i1 = self.value2index(v1)
             if v2 is not None and self.low_value < v2 <= self.high_value:
                 i2 = self.value2index(v2)
-
         return i1, i2
 
     def update_from(self, axis, attributes):