diff --git a/src/utils.py b/src/utils.py
index f375ea6..f5b8067 100644
--- a/src/utils.py
+++ b/src/utils.py
@@ -32,11 +32,20 @@ def crossings(profile_x, profile_z, constant_z):
     indicies = np.where(z[:-1] * z[1:] < 0)[0]
 
     # Use linear interpolation to find intersample crossings.
-    return [
+    x_crossings = [
         profile_x[i] - (profile_x[i] - profile_x[i + 1]) / (z[i] - z[i + 1]) * (z[i])
         for i in indicies
     ]
 
+    # Also need to check the end points as the above will not include them if they are close
+    x_crossings += [
+        profile_x[i]
+        for i in [0, len(profile_z) - 1]
+        if np.isclose(constant_z, profile_z[i])
+    ]
+
+    return sorted(x_crossings)
+
 
 # TODO Think of a better way to do this than having to manually specify the coordinate systems
 def convert_coord_systems(