diff --git a/src/analysis/forecast_twl.py b/src/analysis/forecast_twl.py
index 11ae0ad..e15df6e 100644
--- a/src/analysis/forecast_twl.py
+++ b/src/analysis/forecast_twl.py
@@ -152,9 +152,11 @@ def foreshore_slope_from_profile(profile_x, profile_z, tide, runup_function, **k
         return None
 
     # Initalize estimates
-    max_number_iterations = 20
+    max_number_iterations = 30
     iteration_count = 0
-    min_accuracy = 0.001
+    averaged_accuracy = 0.03  # if slopes within this amount, average after max number of iterations
+    acceptable_accuracy = 0.01  # if slopes within this amount, accept after max number of iterations
+    preferred_accuracy = 0.001  # if slopes within this amount, accept
     beta = 0.05
 
     while True:
@@ -173,12 +175,17 @@ def foreshore_slope_from_profile(profile_x, profile_z, tide, runup_function, **k
             return None
 
         # If slopes do not change much between interactions, return the slope
-        if abs(beta_new - beta) < min_accuracy:
+        if abs(beta_new - beta) < preferred_accuracy:
             return beta
 
         # If we can't converge a solution, return None
         if iteration_count > max_number_iterations:
-            return None
+            if abs(beta_new - beta) < acceptable_accuracy:
+                return beta
+            elif abs(beta_new - beta) < averaged_accuracy:
+                return (beta_new + beta) / 2
+            else:
+                return None
 
         beta = beta_new
         iteration_count += 1
@@ -230,7 +237,7 @@ def slope_from_profile(profile_x, profile_z, top_elevation, btm_elevation, metho
         return -(z_top - z_btm) / (x_top - x_btm)
 
     elif method == "least_squares":
-        profile_mask = [True if end_points["top"]["x"] < pts < end_points["btm"]["x"] else False for pts in x]
+        profile_mask = [True if end_points["top"]["x"] < pts < end_points["btm"]["x"] else False for pts in profile_x]
         slope_x = np.array(profile_x)[profile_mask].tolist()
         slope_z = np.array(profile_z)[profile_mask].tolist()
         slope, _, _, _, _ = stats.linregress(slope_x, slope_z)