Add unit test for numerical derivative

src/mvesuvio/analysis_reduction.py

@@ -12,7 +12,7 @@
                             VesuvioThickness, Integration, Divide, Multiply, DeleteWorkspaces, \
                             CreateWorkspace, CreateSampleWorkspace
 
-from mvesuvio.util.analysis_helpers import loadConstants, numericalThirdDerivative, extend_range_of_array
+from mvesuvio.util.analysis_helpers import loadConstants, numerical_third_derivativ, extend_range_of_array
 
 
 
@@ -651,7 +651,7 @@ def _neutron_compton_profiles(self, pars):
         JOfY = scipy.special.voigt_profile(y_space_arrays_extended - centers, totalGaussWidth, lorzRes)
 
         FSE = (
-            -numericalThirdDerivative(y_space_arrays_extended, JOfY)
+            -numerical_third_derivative(y_space_arrays_extended, JOfY)
             * widths**4
             / deltaQ
             * 0.72

src/mvesuvio/util/analysis_helpers.py

@@ -255,7 +255,7 @@ def extend_range_of_array(arr, n_columns):
     return np.concatenate([left_extend, arr, right_extend], axis=-1)
 
 
-def numericalThirdDerivative(x, y):
+def numerical_third_derivative(x, y):
     k6 = (- y[:, 12:] + y[:, :-12]) * 1
     k5 = (+ y[:, 11:-1] - y[:, 1:-11]) * 24
     k4 = (- y[:, 10:-2] + y[:, 2:-10]) * 192

tests/unit/analysis/test_analysis_helpers.py

@@ -1,11 +1,13 @@
 import unittest
 import numpy as np
+import scipy
 import dill
 import numpy.testing as nptest
 from mock import MagicMock
 from mvesuvio.util.analysis_helpers import extractWS, _convert_dict_to_table,  \
-    fix_profile_parameters, calculate_h_ratio, numericalThirdDerivative, extend_range_of_array
+    fix_profile_parameters, calculate_h_ratio, numerical_third_derivative, extend_range_of_array
 from mantid.simpleapi import CreateWorkspace, DeleteWorkspace
+import matplotlib.pyplot as plt
 
 
 class TestAnalysisHelpers(unittest.TestCase):
@@ -139,5 +141,15 @@ def test_extend_range_of_array_for_decreasing_range(self):
         x_extended = extend_range_of_array(x, 5)
         np.testing.assert_array_equal(x_extended, np.vstack([np.linspace(-7.5, 7.5, 31), np.linspace(-15, 15, 31)]))
 
+
+    def test_numerical_third_derivative(self):
+        x= np.linspace(-20, 20, 300)    # Workspaces are about 300 points of range
+        x = np.vstack([x, 2*x])
+        y = scipy.special.voigt_profile(x, 5, 5)
+        numerical_derivative = numerical_third_derivative(x, y)
+        expected_derivative = np.array([np.gradient(np.gradient(np.gradient(y_i, x_i), x_i), x_i)[6: -6] for y_i, x_i in zip(y, x) ])
+        np.testing.assert_allclose(numerical_derivative, expected_derivative, atol=1e-6)
+
+
 if __name__ == "__main__":
     unittest.main()