Index: DamEngine/trunk/src/Deltares.DamEngine.Calculators.Tests/KernelWrappers/MacroStabilityCommon/PlLinesToWaternetConverterTests.cs
===================================================================
diff -u -r4838 -r4846
--- DamEngine/trunk/src/Deltares.DamEngine.Calculators.Tests/KernelWrappers/MacroStabilityCommon/PlLinesToWaternetConverterTests.cs	(.../PlLinesToWaternetConverterTests.cs)	(revision 4838)
+++ DamEngine/trunk/src/Deltares.DamEngine.Calculators.Tests/KernelWrappers/MacroStabilityCommon/PlLinesToWaternetConverterTests.cs	(.../PlLinesToWaternetConverterTests.cs)	(revision 4846)
@@ -35,6 +35,8 @@
 public class PlLinesToWaternetConverterTests
 {
     private const double precision3Decimals = 0.00051;
+    private const double precision6Decimals = 0.00000051;
+    
     private static readonly Random random = new(21);
     private static readonly double leftCoordinate = random.NextDouble();
     private static readonly double rightCoordinate = leftCoordinate + 1;
@@ -537,6 +539,7 @@
         AssertGeometry([
             new GeometryPoint(leftCoordinate, -10),
             new GeometryPoint(middleXCoordinate, -10),
+            new GeometryPoint(middleXCoordinate, -15),
             new GeometryPoint(rightCoordinate, -15)
         ], pl4WaternetLine.Points);
     }
@@ -643,6 +646,7 @@
         AssertGeometry(new[]
         {
             bottomLeftUpperLayer,
+            bottomIntermediateLowerUpperLayer,
             bottomIntermediateUpperUpperLayer,
             bottomRightUpperLayer
         }, pl4WaternetLine.Points);
@@ -775,6 +779,7 @@
         {
             bottomLeftUpperLayer,
             bottomIntermediateUpperUpperLayer,
+            bottomIntermediateLowerUpperLayer,
             bottomRightUpperLayer
         }, pl4WaternetLine.Points);
     }
@@ -905,6 +910,7 @@
         AssertGeometry(new[]
         {
             bottomLeftUpperLayer,
+            bottomIntermediateUpperUpperLayer,
             bottomIntermediateLowerUpperLayer,
             bottomRightUpperLayer
         }, pl4WaternetLine.Points);
@@ -1343,6 +1349,7 @@
         yield return new TestCaseData(soilProfileRightSoilLayerBottomInRange, new[]
         {
             new GeometryPoint(leftCoordinate, -10),
+            new GeometryPoint(middleXCoordinate, -10),
             new GeometryPoint(middleXCoordinate, -5),
             new GeometryPoint(rightCoordinate, -5)
         }).SetName("Right aquifer only bottom in range");
@@ -1375,6 +1382,7 @@
         {
             new GeometryPoint(leftCoordinate, -10),
             new GeometryPoint(middleXCoordinate, -10),
+            new GeometryPoint(middleXCoordinate, -15),
             new GeometryPoint(rightCoordinate, -15)
         }).SetName("Right aquifer only top in range");
         
@@ -1404,6 +1412,7 @@
         yield return new TestCaseData(soilProfileRightAquiferLayerFullyEnvelopsLeft, new[]
         {
             new GeometryPoint(leftCoordinate, -10),
+            new GeometryPoint(middleXCoordinate, -10),
             new GeometryPoint(middleXCoordinate, -5),
             new GeometryPoint(rightCoordinate, -5)
         }).SetName("Right aquifer fully envelopes left aquifer");
@@ -1436,6 +1445,7 @@
         {
             new GeometryPoint(leftCoordinate, -10),
             new GeometryPoint(middleXCoordinate, -10),
+            new GeometryPoint(middleXCoordinate, -15),
             new GeometryPoint(rightCoordinate, -15)
         }).SetName("Right aquifer fully enveloped by left aquifer");
     }
@@ -1498,8 +1508,8 @@
             GeometryPoint actualPoint = actualPoints.ElementAt(i);
             Assert.Multiple(() =>
             {
-                Assert.That(actualPoint.X, Is.EqualTo(expectedPoint.X));
-                Assert.That(actualPoint.Z, Is.EqualTo(expectedPoint.Z));
+                Assert.That(actualPoint.X, Is.EqualTo(expectedPoint.X).Within(precision6Decimals));
+                Assert.That(actualPoint.Z, Is.EqualTo(expectedPoint.Z).Within(precision6Decimals));
             });
         }
     }
Index: DamEngine/trunk/src/Deltares.DamEngine.Calculators/KernelWrappers/MacroStabilityCommon/PlLinesToWaternetConverter.cs
===================================================================
diff -u -r4540 -r4846
--- DamEngine/trunk/src/Deltares.DamEngine.Calculators/KernelWrappers/MacroStabilityCommon/PlLinesToWaternetConverter.cs	(.../PlLinesToWaternetConverter.cs)	(revision 4540)
+++ DamEngine/trunk/src/Deltares.DamEngine.Calculators/KernelWrappers/MacroStabilityCommon/PlLinesToWaternetConverter.cs	(.../PlLinesToWaternetConverter.cs)	(revision 4846)
@@ -34,6 +34,8 @@
 
 public static class PlLinesToWaternetConverter
 {
+    private const double epsilon = 1e-9;
+    
     /// <summary>
     /// Converts the pl lines to a waternet.
     /// </summary>
@@ -121,7 +123,7 @@
         ThrowWhenPlLinesIsNull(plLines);
         ThrowWhenSoilProfileIsNull(soilProfile);
 
-        // Get all the xCoordinates to make cross sections
+        // Get all the xCoordinates to make cross-sections
         IEnumerable<Point2D> points = soilProfile.Surfaces.SelectMany(surf => surf.GeometrySurface.OuterLoop.CalcPoints);
         double[] xCoordinates = points.Select(point => point.X).OrderBy(x => x).Distinct().ToArray();
 
@@ -219,32 +221,44 @@
         return line;
     }
 
-    private static IEnumerable<Point2D> GetAquiferCoordinates(Func<SoilProfile1D, SoilLayer1D> getSoilLayerFunc, IEnumerable<double> xCoordinates,
+    private static IEnumerable<Point2D> GetAquiferCoordinates(Func<SoilProfile1D, SoilLayer1D> getSoilLayerFunc, double[] xCoordinates,
                                                               SoilProfile2D soilProfile)
     {
         SoilLayer1D previousAquiferLayer = null;
         var coordinates = new List<Point2D>();
+        var xCoordinatesAll = new List<double>();
         foreach (double xCoordinate in xCoordinates)
         {
-            SoilProfile1D crossSection = soilProfile.GetSoilProfile1D(xCoordinate);
-
-            // Determine if the bottom aquifer layer is in range of the previous aquifer layer
-            // If not, return empty coordinates, because the aquifer layer is interrupted
-            SoilLayer1D currentAquifer = getSoilLayerFunc(crossSection);
-            if (previousAquiferLayer != null &&
-                currentAquifer != null
-                && !AreHorizontallyConnected(previousAquiferLayer, currentAquifer))
+            if (HasAquiferAVerticalPartAtGivenX(getSoilLayerFunc, xCoordinate, soilProfile))
             {
-                return Enumerable.Empty<Point2D>();
+                xCoordinatesAll.Add(xCoordinate - epsilon);
+                xCoordinatesAll.Add(xCoordinate + epsilon);
             }
-
-            if (currentAquifer != null)
+            else
             {
-                previousAquiferLayer = currentAquifer;
-                coordinates.Add(new Point2D(xCoordinate, currentAquifer.TopLevel));
+                xCoordinatesAll.Add(xCoordinate);
             }
         }
 
+        foreach (double xCoordinate in xCoordinatesAll)
+        {
+            SoilProfile1D crossSection = soilProfile.GetSoilProfile1D(xCoordinate);
+                // Determine if the aquifer layer is in range of the previous aquifer layer
+                // If not, return empty coordinates, because the aquifer layer is interrupted
+                SoilLayer1D currentAquifer = getSoilLayerFunc(crossSection);
+                if (previousAquiferLayer != null && currentAquifer != null
+                                                 && !AreHorizontallyConnected(previousAquiferLayer, currentAquifer))
+                {
+                    return Enumerable.Empty<Point2D>();
+                }
+
+                if (currentAquifer != null)
+                {
+                    previousAquiferLayer = currentAquifer;
+                    coordinates.Add(new Point2D(xCoordinate, currentAquifer.TopLevel));
+                }
+        }
+
         // Perform a short validation that the coordinates are fully defined from the beginning to the end 
         // of the profile. If not, the aquifer is not continuous.
         if (!coordinates.Any() || coordinates.First().X != xCoordinates.First() || coordinates.Last().X != xCoordinates.Last())
@@ -255,6 +269,22 @@
         return coordinates;
     }
 
+    private static bool HasAquiferAVerticalPartAtGivenX(Func<SoilProfile1D, SoilLayer1D> getSoilLayerFunc, double xCoordinate, SoilProfile2D soilProfile)
+    {
+        SoilProfile1D crossSectionLeft = soilProfile.GetSoilProfile1D(xCoordinate - epsilon);
+        SoilProfile1D crossSectionRight = soilProfile.GetSoilProfile1D(xCoordinate + epsilon);
+        if (crossSectionLeft != null && crossSectionRight != null)
+        {
+            SoilLayer1D currentAquiferLeft = getSoilLayerFunc(crossSectionLeft);
+            SoilLayer1D currentAquiferRight = getSoilLayerFunc(crossSectionRight);
+            if (currentAquiferLeft != null && currentAquiferRight != null)
+            {
+                return Math.Abs(currentAquiferLeft.TopLevel - currentAquiferRight.TopLevel) > GeometryConstants.Accuracy;
+            }
+        }
+        return false;
+    }
+
     private static bool AreHorizontallyConnected(SoilLayer1D leftSoilLayer, SoilLayer1D rightSoilLayer)
     {
         // The layers are identical