remove double valued energies at edges and adding test for interpolation

roentgen/absorption/material.py

@@ -309,6 +309,8 @@ def __init__(self, material):
         self.energy = u.Quantity(data[:, 0] * 1000, "keV")
         self.data = u.Quantity(data[:, 1], "cm^2/g")
 
+        self._remove_double_vals_from_data()
+
         data_energy_kev = np.log10(self.energy.value)
         data_attenuation_coeff = np.log10(self.data.value)
         self._f = interpolate.interp1d(
@@ -318,3 +320,14 @@ def __init__(self, material):
         self.func = lambda x: u.Quantity(
             10 ** self._f(np.log10(x.to("keV").value)), "cm^2/g"
         )
+
+    def _remove_double_vals_from_data(self):
+        """Remove double-values energy values. Edges are represented with
+        the same energy index and at the bottom and top value of the edge. This
+        must be removed to enable correct interpolation."""
+        uniq, count = np.unique(self.energy, return_counts=True)
+        duplicates = uniq[count > 1]
+        for this_dup in duplicates:
+            ind = (self.energy == this_dup).nonzero()
+            # shift the first instance of the energy, the bottom of the edge
+            self.energy[ind[0][0]] -= 1e-3 * u.eV

roentgen/tests/test_massatten.py

@@ -0,0 +1,18 @@
+import pytest
+
+from numpy.testing import assert_allclose
+
+import astropy.units as u
+
+from roentgen.absorption.material import MassAttenuationCoefficient
+
+
+def test_interpolate_matches_at_data():
+    for this_element in ['Te', 'Si', 'Ge', 'cdte']:
+        te = MassAttenuationCoefficient(this_element)
+        assert_allclose(te.data, te.func(te.energy))
+
+        # now change one point and make sure it no longer works
+        te.energy[0] = 10 * u.keV
+        with pytest.raises(AssertionError):
+            assert_allclose(te.data, te.func(te.energy))