Index: DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometrySurface.cs =================================================================== diff -u -r5061 -r5069 --- DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometrySurface.cs (.../GeometrySurface.cs) (revision 5061) +++ DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometrySurface.cs (.../GeometrySurface.cs) (revision 5069) @@ -35,6 +35,9 @@ { private GeometryLoop outerLoop = new GeometryLoop(); + private readonly List previousInnerLoops = new List(); + private GeometryLoop previousOuterLoop = new GeometryLoop(); + /// /// Empty constructor /// @@ -63,13 +66,33 @@ outerLoop = value; } } - + /// + /// Sets the outer loop to loop whilst storing the "old" value as previousOuterLoop. + /// + /// + public void SetOuterLoop(GeometryLoop loop) + { + if (loop == OuterLoop) + return; + previousOuterLoop = OuterLoop; + OuterLoop = loop; + } + + /// /// All internal loops encompassed by . /// public List InnerLoops { get; } = new List(); /// + /// Placeholder for the previous outer loop. + /// + public GeometryLoop PreviousOuterLoop => this.previousOuterLoop; + + + public List PreviousInnerLoops => this.previousInnerLoops; + + /// /// Add to if it isn't already /// in that collection. /// @@ -80,6 +103,18 @@ InnerLoops.Add(aLoop); } } + + /// + /// Remove all inner loops from the surface whilst updating the previous inner loops list. + /// + public void RemoveAllInnerLoops() + { + previousInnerLoops.Clear(); + if (InnerLoops.Count <= 0) + return; + previousInnerLoops.AddRange((IEnumerable) InnerLoops); + InnerLoops.Clear(); + } /// /// determine top of geometry surface outerloop as GeometryPointString Index: DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometryLoop.cs =================================================================== diff -u -r5061 -r5069 --- DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometryLoop.cs (.../GeometryLoop.cs) (revision 5061) +++ DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometryLoop.cs (.../GeometryLoop.cs) (revision 5069) @@ -252,7 +252,8 @@ GeometryPointString clonedGeometryPointString = base.Clone(); GeometryLoop clonedGeometryLoop = new GeometryLoop(); clonedGeometryLoop.Points.AddRange(clonedGeometryPointString.Points); - clonedGeometryLoop.CurveList.AddRange(CurveList); + clonedGeometryLoop.SyncCalcPoints(); //Always sync the calcPoints after adding points or preferably add the calcPoints to the clonedGeometryPointString and synchronize the points (using SyncPoints) + clonedGeometryLoop.CurveList.AddRange(CurveList); // Wrong: The points in the new curves must be retrieved by location in the new pointlist in order to properly clone the curves return clonedGeometryLoop; } } \ No newline at end of file Index: DamEngine/trunk/src/Deltares.DamEngine.Data.Tests/Geotechnics/SoilProfile2DSurfaceLineHelperTests.cs =================================================================== diff -u -r5053 -r5069 --- DamEngine/trunk/src/Deltares.DamEngine.Data.Tests/Geotechnics/SoilProfile2DSurfaceLineHelperTests.cs (.../SoilProfile2DSurfaceLineHelperTests.cs) (revision 5053) +++ DamEngine/trunk/src/Deltares.DamEngine.Data.Tests/Geotechnics/SoilProfile2DSurfaceLineHelperTests.cs (.../SoilProfile2DSurfaceLineHelperTests.cs) (revision 5069) @@ -314,7 +314,7 @@ Name = "dikemat" }; SoilProfile2D newSoilProfile2D = SoilProfile2DSurfaceLineHelper.CombineSurfaceLineWithSoilProfile2D( - testCaseSurfaceLine.SurfaceLine, soilProfile2D, defaultSoil, 0); + testCaseSurfaceLine.SurfaceLine.Geometry, soilProfile2D, defaultSoil, 0); Assert.That(newSoilProfile2D, Is.Not.Null); Assert.That(newSoilProfile2D.Surfaces, Has.Count.EqualTo(testCaseSurfaceLine.SurfaceCount)); } @@ -331,7 +331,7 @@ }; // When SoilProfile2D newSoilProfile2D = SoilProfile2DSurfaceLineHelper.CombineSurfaceLineWithSoilProfile2D( - testCaseSurfaceLine.GivenZigZagSurfaceLine, soilProfile2D, defaultSoil, 0); + testCaseSurfaceLine.GivenZigZagSurfaceLine.Geometry, soilProfile2D, defaultSoil, 0); // Then if (testCaseSurfaceLine.ExpectedSurfaceCount == 0) { @@ -445,7 +445,7 @@ return surfaceLine; } - private void CheckSoilProfileContainsSoilLayer(SoilProfile2D soilProfile2D, SoilLayer2D expectedSoilLayer) + private static void CheckSoilProfileContainsSoilLayer(SoilProfile2D soilProfile2D, SoilLayer2D expectedSoilLayer) { if (expectedSoilLayer != null) { Index: DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometryGenerator.cs =================================================================== diff -u -r5034 -r5069 --- DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometryGenerator.cs (.../GeometryGenerator.cs) (revision 5034) +++ DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometryGenerator.cs (.../GeometryGenerator.cs) (revision 5069) @@ -166,8 +166,7 @@ if (aLoop.IsLoop() && aLoop.HasArea()) { var newSurface = new GeometrySurface(); - //newSurface.SetOuterLoop(aLoop); #Bka: replaced with - newSurface.OuterLoop = aLoop; + newSurface.SetOuterLoop(aLoop); geometryData.Surfaces.Add(newSurface); return newSurface; } @@ -844,8 +843,7 @@ int newPointCount = newLoop.CurveList.Count; List newPolygon = newLoop.CalcPoints; - newSurface.InnerLoops.Clear(); - + newSurface.RemoveAllInnerLoops(); for (var index = 0; index < surfaceCount; index++) { if (newSurface == geometryData.Surfaces[index]) @@ -1010,37 +1008,37 @@ private void MergePoints() { //try + //{ + int count = geometryData.Points.Count; + List aConnectedAtHeadCurveList = new List(); + List aConnectedAtEndCurveList = new List(); + List list = geometryData.NewlyEffectedPoints.Select((Func) + (newPoint => geometryData.Points.IndexOf(newPoint))).ToList(); + list.Reverse(); + foreach (int index1 in list) { - int count = geometryData.Points.Count; - List aConnectedAtHeadCurveList = new List(); - List aConnectedAtEndCurveList = new List(); - List list = geometryData.NewlyEffectedPoints.Select((Func) (newPoint => geometryData.Points.IndexOf(newPoint))).ToList(); - list.Reverse(); - foreach (int index1 in list) + if (index1 >= 0 && index1 < geometryData.Points.Count) { - if (index1 >= 0 && index1 < geometryData.Points.Count) + for (int index2 = count - 1; index2 >= 0; --index2) { - for (int index2 = count - 1; index2 >= 0; --index2) + if ((index2 != index1) && (Routines2D.DetermineIfPointsCoincide( + geometryData.Points[index1].X, geometryData.Points[index1].Z, + geometryData.Points[index2].X, geometryData.Points[index2].Z, 0.001))) { - if (index2 != index1) - { - if (Routines2D.DetermineIfPointsCoincide(geometryData.Points[index1].X, geometryData.Points[index1].Z, - geometryData.Points[index2].X, geometryData.Points[index2].Z, 0.001)) - { - aConnectedAtHeadCurveList.Clear(); - aConnectedAtEndCurveList.Clear(); - GetConnectedCurves(this.geometryData.Points[index1], ref aConnectedAtHeadCurveList, ref aConnectedAtEndCurveList); - SetPointInCurves(aConnectedAtHeadCurveList, geometryData.Points[index2], true); - SetPointInCurves(aConnectedAtEndCurveList, geometryData.Points[index2], false); - --count; - geometryData.DeletePoint(geometryData.Points[index1]); - break; - } - } + aConnectedAtHeadCurveList.Clear(); + aConnectedAtEndCurveList.Clear(); + GetConnectedCurves(this.geometryData.Points[index1], ref aConnectedAtHeadCurveList, ref aConnectedAtEndCurveList); + SetPointInCurves(aConnectedAtHeadCurveList, geometryData.Points[index2], true); + SetPointInCurves(aConnectedAtEndCurveList, geometryData.Points[index2], false); + --count; + geometryData.DeletePoint(geometryData.Points[index1]); + break; } + } } } + //} // catch (Exception ex) // { // LogManager.Add(new LogMessage(LogMessageType.FatalError, (object) this, "MergePoints: " + ex.Message)); @@ -1062,23 +1060,7 @@ double num2 = !pointXz3.Z.IsNearEqual(pointXz4.Z) ? (pointXz3.X - pointXz4.X) / (pointXz3.Z - pointXz4.Z) : double.MaxValue; if (Math.Abs(Math.Abs(num1) - Math.Abs(num2)) > 5E-12) { - if (flag3 & flag4) - { - if (Routines2D.CalculateDistanceToLine(pointXz1.X, pointXz1.Z, pointXz3.X, pointXz3.Z, pointXz4.X, pointXz4.Z) > - Routines2D.CalculateDistanceToLine(pointXz2.X, pointXz2.Z, pointXz3.X, pointXz3.Z, pointXz4.X, pointXz4.Z)) - flag3 = false; - else - flag4 = false; - } - - if (flag1 & flag2) - { - if (Routines2D.CalculateDistanceToLine(pointXz3.X, pointXz3.Z, pointXz1.X, pointXz1.Z, pointXz2.X, pointXz2.Z) > - Routines2D.CalculateDistanceToLine(pointXz4.X, pointXz4.Z, pointXz1.X, pointXz1.Z, pointXz2.X, pointXz2.Z)) - flag1 = false; - else - flag2 = false; - } + flag3 = DetermineFlagValues1to4(flag3, pointXz1, pointXz3, pointXz4, pointXz2, ref flag4, ref flag1, ref flag2); } bool flag5 = Routines2D.DetermineIfPointsCoincide(pointXz1.X, pointXz1.Z, pointXz3.X, pointXz3.Z, 0.001); @@ -1087,32 +1069,16 @@ bool flag8 = Routines2D.DetermineIfPointsCoincide(pointXz2.X, pointXz2.Z, pointXz4.X, pointXz4.Z, 0.001); if (Math.Abs(Math.Abs(num1) - Math.Abs(num2)) > 5E-12) { - if (flag5 & flag7) - { - if (Routines2D.Compute2DDistance(pointXz1.X, pointXz1.Z, pointXz3.X, pointXz3.Z) > - Routines2D.Compute2DDistance(pointXz2.X, pointXz2.Z, pointXz3.X, pointXz3.Z)) - flag5 = false; - else - flag7 = false; - } - - if (flag6 & flag8) - { - if (Routines2D.Compute2DDistance(pointXz1.X, pointXz1.Z, pointXz4.X, pointXz4.Z) > - Routines2D.Compute2DDistance(pointXz2.X, pointXz2.Z, pointXz4.X, pointXz4.Z)) - flag6 = false; - else - flag8 = false; - } + DetermineFlagValues5to8(pointXz1, pointXz3, pointXz2, pointXz4, ref flag5, ref flag6, ref flag7, ref flag8); } - if (flag3 & flag4 || flag1 & flag2 || flag1 | flag2 && flag3 | flag4) + if (flag3 && flag4 || flag1 && flag2 || flag1 | flag2 && flag3 || flag4) { if (flag1 && !flag2 && !flag5 && !flag7) { GeometryCurve geometryCurve = this.SplitCurve(line1, line2.HeadPoint); line2.HeadPoint = flag4 ? geometryCurve.EndPoint : line1.HeadPoint; - this.CheckAndAddToIntersectedCurveList(line2); + CheckAndAddToIntersectedCurveList(line2); isCurveInserted = true; return true; } @@ -1121,36 +1087,36 @@ { GeometryCurve geometryCurve = this.SplitCurve(line1, line2.EndPoint); line2.EndPoint = flag4 ? geometryCurve.EndPoint : line1.HeadPoint; - this.CheckAndAddToIntersectedCurveList(line2); + CheckAndAddToIntersectedCurveList(line2); isCurveInserted = true; return true; } - if (flag1 & flag2 && !flag5 && !flag8 && !flag7 && !flag6) + if (flag1 && flag2 && !flag5 && !flag8 && !flag7 && !flag6) { if (Routines2D.Compute2DDistance(pointXz1.X, pointXz1.Z, pointXz3.X, pointXz3.Z) < Routines2D.Compute2DDistance(pointXz1.X, pointXz1.Z, pointXz4.X, pointXz4.Z)) - this.SplitCurve(line1, line2.HeadPoint).HeadPoint = line2.EndPoint; + SplitCurve(line1, line2.HeadPoint).HeadPoint = line2.EndPoint; else - this.SplitCurve(line1, line2.EndPoint).HeadPoint = line2.HeadPoint; - this.CheckAndAddToIntersectedCurveList(line2); + SplitCurve(line1, line2.EndPoint).HeadPoint = line2.HeadPoint; + CheckAndAddToIntersectedCurveList(line2); isCurveInserted = true; return true; } - if (flag3 & flag4 && !flag5 && !flag8 && !flag7 && !flag6) + if (flag3 && flag4 && !flag5 && !flag8 && !flag7 && !flag6) { if (Routines2D.Compute2DDistance(pointXz3.X, pointXz3.Z, pointXz1.X, pointXz1.Z) < Routines2D.Compute2DDistance(pointXz3.X, pointXz3.Z, pointXz2.X, pointXz2.Z)) - this.SplitCurve(line2, line1.HeadPoint).HeadPoint = line1.EndPoint; + SplitCurve(line2, line1.HeadPoint).HeadPoint = line1.EndPoint; else - this.SplitCurve(line2, line1.EndPoint).HeadPoint = line1.HeadPoint; - this.CheckAndAddToIntersectedCurveList(line1); + SplitCurve(line2, line1.EndPoint).HeadPoint = line1.HeadPoint; + CheckAndAddToIntersectedCurveList(line1); isCurveInserted = true; return true; } - if (flag2 & flag1 && !(flag7 & flag6) && !(flag8 & flag5)) + if (flag2 && flag1 && !(flag7 && flag6) && !(flag8 && flag5)) { if (flag5) line1.HeadPoint = line2.EndPoint; @@ -1160,8 +1126,8 @@ line1.EndPoint = line2.HeadPoint; else line1.HeadPoint = line2.HeadPoint; - this.CheckAndAddToIntersectedCurveList(line1); - this.CheckAndAddToIntersectedCurveList(line2); + CheckAndAddToIntersectedCurveList(line1); + CheckAndAddToIntersectedCurveList(line2); return true; } @@ -1188,6 +1154,50 @@ return true; } + private static void DetermineFlagValues5to8(Point2D pointXz1, Point2D pointXz3, Point2D pointXz2, Point2D pointXz4, ref bool flag5, ref bool flag6, ref bool flag7, ref bool flag8) + { + if (flag5 && flag7) + { + if (Routines2D.Compute2DDistance(pointXz1.X, pointXz1.Z, pointXz3.X, pointXz3.Z) > + Routines2D.Compute2DDistance(pointXz2.X, pointXz2.Z, pointXz3.X, pointXz3.Z)) + flag5 = false; + else + flag7 = false; + } + + if (flag6 && flag8) + { + if (Routines2D.Compute2DDistance(pointXz1.X, pointXz1.Z, pointXz4.X, pointXz4.Z) > + Routines2D.Compute2DDistance(pointXz2.X, pointXz2.Z, pointXz4.X, pointXz4.Z)) + flag6 = false; + else + flag8 = false; + } + } + + private static bool DetermineFlagValues1to4(bool flag3, Point2D pointXz1, Point2D pointXz3, Point2D pointXz4, Point2D pointXz2, ref bool flag4, ref bool flag1, ref bool flag2) + { + if (flag3 && flag4) + { + if (Routines2D.CalculateDistanceToLine(pointXz1.X, pointXz1.Z, pointXz3.X, pointXz3.Z, pointXz4.X, pointXz4.Z) > + Routines2D.CalculateDistanceToLine(pointXz2.X, pointXz2.Z, pointXz3.X, pointXz3.Z, pointXz4.X, pointXz4.Z)) + flag3 = false; + else + flag4 = false; + } + + if (flag1 && flag2) + { + if (Routines2D.CalculateDistanceToLine(pointXz3.X, pointXz3.Z, pointXz1.X, pointXz1.Z, pointXz2.X, pointXz2.Z) > + Routines2D.CalculateDistanceToLine(pointXz4.X, pointXz4.Z, pointXz1.X, pointXz1.Z, pointXz2.X, pointXz2.Z)) + flag1 = false; + else + flag2 = false; + } + + return flag3; + } + private GeometryCurve SplitCurve(GeometryCurve curve, Point2D aPointOnCurve) { if (aPointOnCurve.LocationEquals(curve.HeadPoint) || aPointOnCurve.LocationEquals(curve.EndPoint)) Index: DamEngine/trunk/src/Deltares.DamEngine.Data/Geotechnics/SoilProfile2D.cs =================================================================== diff -u -r5061 -r5069 --- DamEngine/trunk/src/Deltares.DamEngine.Data/Geotechnics/SoilProfile2D.cs (.../SoilProfile2D.cs) (revision 5061) +++ DamEngine/trunk/src/Deltares.DamEngine.Data/Geotechnics/SoilProfile2D.cs (.../SoilProfile2D.cs) (revision 5069) @@ -19,7 +19,9 @@ // Stichting Deltares and remain full property of Stichting Deltares at all times. // All rights reserved. +using System; using System.Collections.Generic; +using System.Runtime.CompilerServices; using Deltares.DamEngine.Data.Geometry; using Deltares.DamEngine.Data.Standard.Language; using Deltares.DamEngine.Data.Standard.Validation; @@ -171,7 +173,7 @@ SoilLayer2D clonedSurface = surface.Clone(); clonedSoilProfile2D.Surfaces.Add(clonedSurface); } - + clonedSoilProfile2D.Geometry = geometry.Clone(); return clonedSoilProfile2D; } @@ -185,4 +187,172 @@ { return Name; } + + /// + /// Get the soil layer from the old surfaces + /// + /// + /// + /// + public SoilLayer2D GetOriginalLayerFromOldSurfaces(GeometrySurface geometrySurface, IEnumerable oldSurfaces) + { + return GetOriginalLayerFromOldSurfaces(geometrySurface, oldSurfaces, 0.0); + } + + /// + /// Get the soil layer from the old surfaces + /// + /// + /// + /// + /// + private SoilLayer2D GetOriginalLayerFromOldSurfaces(GeometrySurface geometrySurface, IEnumerable oldSurfaces, + double shift) + { + Point2D point = new Point2D(0.0, 0.0); + bool flag1 = false; + foreach (GeometryCurve curve in geometrySurface.OuterLoop.CurveList) + { + point = new Point2D((curve.HeadPoint.X + curve.EndPoint.X) / 2.0, (curve.HeadPoint.Z + curve.EndPoint.Z) / 2.0); + if (IsPointWithinOldSurfaces(point, oldSurfaces, shift, 1E-05) || IsPointWithinOldSurfaces(point, oldSurfaces, shift, -1E-05)) + { + point.Z += 1E-05; + if (Routines2D.CheckIfPointIsInPolygon(geometrySurface.OuterLoop, point.X, point.Z) == PointInPolygon.InsidePolygon) + { + flag1 = true; + break; + } + point.Z -= 2E-05; + if (Routines2D.CheckIfPointIsInPolygon(geometrySurface.OuterLoop, point.X, point.Z) == PointInPolygon.InsidePolygon) + { + flag1 = true; + break; + } + } + } + if (!flag1) + { + GeometryPointString topGeometrySurface = geometrySurface.DetermineTopGeometrySurface(); + topGeometrySurface.SortPointsByXAscending(); + Point2D geometryPoint1 = topGeometrySurface[0]; + geometryPoint1.X -= 1E-05; + geometryPoint1.Z -= 1E-05; + Point2D geometryPoint2 = topGeometrySurface[checked (topGeometrySurface.Count - 1)]; + geometryPoint2.X += 1E-05; + geometryPoint2.Z -= 1E-05; + int num = this.IsPointWithinOldSurfaces(geometryPoint1, oldSurfaces, shift, -1E-05) ? 1 : 0; + bool flag2 = this.IsPointWithinOldSurfaces(geometryPoint2, oldSurfaces, shift, -1E-05); + double d = double.NaN; + if (num != 0 && !flag2) + { + point.X = geometryPoint1.X; + point.Z = geometryPoint1.Z; + flag1 = true; + d = geometryPoint1.X + 1E-05; + } + if (num == 0 & flag2) + { + point.X = geometryPoint2.X; + point.Z = geometryPoint2.Z; + flag1 = true; + d = geometryPoint2.X - 1E-05; + } + if (!double.IsNaN(d)) + { + double xminFromSurfaces = GetXminFromSurfaces(oldSurfaces); + double xmaxFromSurfaces = GetXmaxFromSurfaces(oldSurfaces); + if (d <= xmaxFromSurfaces && d >= xminFromSurfaces) + flag1 = false; + } + } + if (flag1) + { + point.X -= shift; + foreach (SoilLayer2D oldSurface in oldSurfaces) + { + GeometryLoop previousOuterLoop = oldSurface.GeometrySurface.PreviousOuterLoop; + if (previousOuterLoop != null && previousOuterLoop.CurveList.Count > 2 && + Routines2D.CheckIfPointIsInPolygon(previousOuterLoop, point.X, point.Z) == PointInPolygon.InsidePolygon) + { + bool flag3 = true; + foreach (GeometryLoop previousInnerLoop in oldSurface.GeometrySurface.PreviousInnerLoops) + { + if (Routines2D.CheckIfPointIsInPolygon(previousInnerLoop, point.X, point.Z) == PointInPolygon.InsidePolygon) + flag3 = false; + } + if (flag3) + return oldSurface; + } + } + foreach (SoilLayer2D oldSurface in oldSurfaces) + { + GeometryLoop outerLoop = oldSurface.GeometrySurface.OuterLoop; + if (outerLoop != null && outerLoop.CurveList.Count > 2 && Routines2D.CheckIfPointIsInPolygon(outerLoop, point.X, point.Z) == PointInPolygon.InsidePolygon) + { + bool flag4 = true; + foreach (GeometryLoop innerLoop in oldSurface.GeometrySurface.InnerLoops) + { + if (Routines2D.CheckIfPointIsInPolygon(innerLoop, point.X, point.Z) == PointInPolygon.InsidePolygon) + flag4 = false; + } + if (flag4) + return oldSurface; + } + } + } + return null; + } + + private double GetXminFromSurfaces(IEnumerable oldSurfaces) + { + double val1 = double.MaxValue; + foreach (SoilLayer2D oldSurface in oldSurfaces) + val1 = Math.Min(val1, oldSurface.GeometrySurface.OuterLoop.GetMinX()); + return val1; + } + + private double GetXmaxFromSurfaces(IEnumerable oldSurfaces) + { + double val1 = double.MinValue; + foreach (SoilLayer2D oldSurface in oldSurfaces) + val1 = Math.Max(val1, oldSurface.GeometrySurface.OuterLoop.GetMaxX()); + return val1; + } + + private bool IsPointWithinOldSurfaces(Point2D point, IEnumerable oldSurfaces, double shift, double deviation) + { + bool flag = false; + point.X -= shift; + point.Z += deviation; + foreach (SoilLayer2D oldSurface in oldSurfaces) + { + GeometryLoop outerLoop = oldSurface.GeometrySurface.OuterLoop; + if (outerLoop != null && Routines2D.CheckIfPointIsInPolygon(outerLoop, point.X, point.Z) == PointInPolygon.InsidePolygon) + { + flag = true; + foreach (GeometryLoop innerLoop in oldSurface.GeometrySurface.InnerLoops) + { + if (Routines2D.CheckIfPointIsInPolygon(innerLoop, point.X, point.Z) == PointInPolygon.InsidePolygon) + flag = false; + } + } + if (!flag) + { + GeometryLoop previousOuterLoop = oldSurface.GeometrySurface.PreviousOuterLoop; + if (previousOuterLoop != null && Routines2D.CheckIfPointIsInPolygon(previousOuterLoop, point.X, point.Z) == PointInPolygon.InsidePolygon) + { + flag = true; + foreach (GeometryLoop previousInnerLoop in oldSurface.GeometrySurface.PreviousInnerLoops) + { + if (Routines2D.CheckIfPointIsInPolygon(previousInnerLoop, point.X, point.Z) == PointInPolygon.InsidePolygon) + flag = false; + } + } + } + if (flag) + return true; + } + return false; + } + } \ No newline at end of file Index: DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometryPointString.cs =================================================================== diff -u -r4904 -r5069 --- DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometryPointString.cs (.../GeometryPointString.cs) (revision 4904) +++ DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometryPointString.cs (.../GeometryPointString.cs) (revision 5069) @@ -191,6 +191,18 @@ } /// + /// Round all points to the nearest geometry accuracy + /// + public void RoundPointsCoordinates() + { + foreach (GeometryPoint point in (IEnumerable) this.Points) + { + point.X = RoundValue(point.X); + point.Z = RoundValue(point.Z); + } + } + + /// /// Gets the minimum z. /// Note: CalcPoints must be used instead of calcPoints as otherwise unclear behaviour of Linq spoils the result /// Change back to calcPoints and Benchmark4_01l in the WaternetCreatorTests will fail! @@ -317,7 +329,7 @@ /// /// Uses constant extrapolation and linear interpolation. /// - public virtual double GetZatX(double x) + public double GetZatX(double x) { if (calcPoints.Any()) { @@ -364,7 +376,7 @@ /// The x. /// The i. /// - public virtual double GetZatX(double x, ref int i) + public double GetZatX(double x, ref int i) { for (; i < calcPoints.Count - 1; i++) { @@ -794,10 +806,94 @@ /// /// /// - private List IntersectionXzPointsWithGeometryPointList(IList list) + public List IntersectionXzPointsWithGeometryPointList(IList list) { return IntersectWithPointsListCore(list, false); } + + /// + /// Determines the relative position of a point to the geometric line in the XZ-plane. + /// + /// + /// + /// + public RelativeXzPosition GetPosition(GeometryLoop geometryLoop, ExtraPolationMode extraPolationMode = ExtraPolationMode.Beyond) + { + foreach (GeometryPoint point in geometryLoop.Points.Concat(GeometryPointString.GetLoopSegmentMiddlePoints(geometryLoop))) + { + RelativeXzPosition position = PositionXzOfPointRelatedToExtrapolatedLine(point, extraPolationMode); + if (position == RelativeXzPosition.BeyondGeometricLine) + position = RelativeXzPosition.BelowGeometricLine; + if (position != RelativeXzPosition.OnGeometricLine) + return position; + } + return RelativeXzPosition.OnGeometricLine; + } + + private RelativeXzPosition PositionXzOfPointRelatedToExtrapolatedLine(GeometryPoint point, + ExtraPolationMode extraPolationMode = ExtraPolationMode.Beyond) + { + return IsPointConsideredBeyondLine(point, extraPolationMode) ? + RelativeXzPosition.BeyondGeometricLine : DeterminePositionWithRespectToExtrapolatedLine(point, extraPolationMode); + } + + private bool IsPointConsideredBeyondLine(GeometryPoint point, ExtraPolationMode extraPolationMode) + { + if (this.Points.Count == 0) + return true; + return this.IsPointOutsideRangeX(point) && extraPolationMode == ExtraPolationMode.Beyond; + } + + private bool IsPointOutsideRangeX(GeometryPoint point) + { + return point.X < Points[0].X || point.X > Points[checked (Points.Count - 1)].X; + } + + private RelativeXzPosition DeterminePositionWithRespectToExtrapolatedLine(GeometryPoint point, + ExtraPolationMode extraPolationMode) + { + double xToEvaluate = GetXToEvaluate(point, extraPolationMode); + var zAtX = GetZatX(xToEvaluate); + //double zatX = GetZAtX(xToEvaluate); + if (Math.Abs(point.Z - zAtX) < 0.001) + return RelativeXzPosition.OnGeometricLine; + return point.Z <= zAtX ? RelativeXzPosition.BelowGeometricLine : RelativeXzPosition.AboveGeometricLine; + } + + private double GetXToEvaluate(GeometryPoint point, ExtraPolationMode extraPolationMode) + { + double x = point.X; + if (extraPolationMode == ExtraPolationMode.Horizontal) + { + if (x < Points[0].X) + x = Points[0].X; + else if (x > Points[checked (Points.Count - 1)].X) + x = Points[checked (Points.Count - 1)].X; + } + return x; + } + + private static IEnumerable GetLoopSegmentMiddlePoints(GeometryLoop geometryLoop) + { + List loopSegmentMiddlePoints = new List(); + int index = 0; + while (index < checked (geometryLoop.Points.Count - 1)) + { + GeometryPointString.AddLoopSegmentMiddlePoint(geometryLoop.Points[index], geometryLoop.Points[checked (index + 1)], loopSegmentMiddlePoints); + checked { ++index; } + } + if (geometryLoop.Points.Count > 2) + GeometryPointString.AddLoopSegmentMiddlePoint(geometryLoop.Points[0], geometryLoop.Points.Last(), loopSegmentMiddlePoints); + return (IEnumerable) loopSegmentMiddlePoints; + } + + private static void AddLoopSegmentMiddlePoint(GeometryPoint point1, GeometryPoint point2, List loopSegmentMiddlePoints) + { + if (Math.Abs(point1.Z - point2.Z) < 0.001) + loopSegmentMiddlePoints.Insert(0, new GeometryPoint((point1.X + point2.X) / 2.0, point1.Z)); + else + loopSegmentMiddlePoints.Add(new GeometryPoint((point1.X + point2.X) / 2.0, (point1.Z + point2.Z) / 2.0)); + } /// /// Checks if constant extrapolation can be applied, and if so set the Z value. Index: DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometryData.cs =================================================================== diff -u -r5065 -r5069 --- DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometryData.cs (.../GeometryData.cs) (revision 5065) +++ DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometryData.cs (.../GeometryData.cs) (revision 5069) @@ -1314,10 +1314,12 @@ foreach (GeometryCurve curve in this.Curves) { clonedGeometryData.Curves.Add(curve.Clone()); + // Wrong: The points in the new curves must be retrieved by location in the new pointlist in order to properly clone the curves } foreach (GeometryCurve curve in this.NewlyEffectedCurves) { clonedGeometryData.NewlyEffectedCurves.Add(curve.Clone()); + // Wrong: The points in the new curves must be retrieved by location in the new effectedpointlist in order to properly clone the curves } foreach (GeometryLoop loop in this.Loops) { @@ -1327,6 +1329,8 @@ { clonedGeometryData.Surfaces.Add(surface.Clone()); } + // for the clone, use the same geometry generator for now. Should be a new one but as the cloning is not correct, that will make some tests fail + clonedGeometryData.geometryGenerator = geometryGenerator;//new GeometryGenerator(this); return clonedGeometryData; } } \ No newline at end of file Index: DamEngine/trunk/src/Deltares.DamEngine.Data/Geotechnics/SoilProfile2DSurfaceLineHelper.cs =================================================================== diff -u -r5010 -r5069 --- DamEngine/trunk/src/Deltares.DamEngine.Data/Geotechnics/SoilProfile2DSurfaceLineHelper.cs (.../SoilProfile2DSurfaceLineHelper.cs) (revision 5010) +++ DamEngine/trunk/src/Deltares.DamEngine.Data/Geotechnics/SoilProfile2DSurfaceLineHelper.cs (.../SoilProfile2DSurfaceLineHelper.cs) (revision 5069) @@ -20,7 +20,10 @@ // All rights reserved. using System; +using System.Collections.Generic; +using System.Linq; using Deltares.DamEngine.Data.Geometry; +using Deltares.DamEngine.Data.Standard.Language; namespace Deltares.DamEngine.Data.Geotechnics; @@ -66,10 +69,46 @@ /// The default soil. /// The required shift of the geometry so it will fit the surface line /// - public static SoilProfile2D CombineSurfaceLineWithSoilProfile2D(SurfaceLine2 surfaceLine, SoilProfile2D soilProfile2D, + public static SoilProfile2D CombineSurfaceLineWithSoilProfile2D(GeometryPointString surfaceLine, SoilProfile2D soilProfile2D, Soil defaultSoil, double shift) { - return null; + if (soilProfile2D == null || soilProfile2D.Surfaces.Count == 0) + return null; + if (surfaceLine == null || surfaceLine.Points.Count == 0) + return null; + //this.ThrowIfSurfaceLineIsAtLeastPartiallyBeneathGeometry(surfaceLine, soilProfile2D); + SoilProfile2D result = new SoilProfile2D(); + + SoilProfile2D clonedProfile = soilProfile2D.Clone(); + GeometryPointString clonedSurfaceline = surfaceLine.Clone(); + RoundCoordinates(clonedSurfaceline, clonedProfile); + shift = clonedSurfaceline.RoundValue(shift); + if (Math.Abs(shift) >= 0.0) + { + foreach (Point2D point in clonedProfile.Geometry.Points) + point.X += shift; + clonedProfile.Geometry.Rebox(); + } + + if (clonedProfile.Geometry.Right < surfaceLine.GetMinX() || clonedProfile.Geometry.Left > surfaceLine.GetMaxX()) + { + throw new ArgumentException("Surface line is beyond the profile."); + } + // As the original profile may have been moved as cloned profile, a new clone is needed to perform all reset actions with. + // The cloned profile is used to determine the original surfaces and preconsolidations. + SoilProfile2D soilProfile2D2 = clonedProfile.Clone(); + List oldSurfaces = new List(); + oldSurfaces.AddRange((IEnumerable) soilProfile2D2.Surfaces); + AddGeometryIntersectionsToSurfaceLine(clonedProfile.Geometry, ref clonedSurfaceline); + BuildGeometryModel(clonedSurfaceline, clonedProfile, ref result); + RemoveGeometryDataOfSoilProfileAboveSurfaceLine(clonedSurfaceline, ref result); + ReconstructSurfaces(ref result, clonedProfile, oldSurfaces, defaultSoil); + ReconstructPreConsolidations(ref result, clonedProfile, shift); + //ValidationResult[] source = result.Geometry.ValidateGeometry(); + //if (((IEnumerable) source).Any()) + // throw new GeotechnicsUtilitiesException(LocalizationManager.GetTranslatedText((object) this, "CombinedGeometry") + source[0].Text); + result.Geometry.Rebox(); + return result; } /// @@ -91,4 +130,233 @@ return true; } + + private static void RoundCoordinates(GeometryPointString surfaceLine, SoilProfile2D soilProfile2D) + { + foreach (Point2D point in soilProfile2D.Geometry.Points) + { + point.X = point.RoundValue(point.X); + point.Z = point.RoundValue(point.Z); + } + + soilProfile2D.Geometry.Right = soilProfile2D.Geometry.RoundValue(soilProfile2D.Geometry.Right); + soilProfile2D.Geometry.Left = soilProfile2D.Geometry.RoundValue(soilProfile2D.Geometry.Left); + soilProfile2D.Geometry.Bottom = soilProfile2D.Geometry.RoundValue(soilProfile2D.Geometry.Bottom); + surfaceLine.RoundPointsCoordinates(); + } + private static void AddGeometryIntersectionsToSurfaceLine(GeometryData geometry, ref GeometryPointString cloneSurfaceline) + { + List pointList = new List(); + foreach (GeometrySurface surface in geometry.Surfaces) + { + surface.OuterLoop.SyncCalcPoints(); + foreach (Point2D intersectionPoint in cloneSurfaceline.IntersectionXzPointsWithGeometryPointList(surface.OuterLoop.CalcPoints)) + pointList.Add(intersectionPoint); + } + foreach (Point2D point in pointList) + { + if (cloneSurfaceline.GetPointAt(point.X, point.Z) == null) + cloneSurfaceline.CalcPoints.Add(point); + } + cloneSurfaceline.SyncPoints(); + cloneSurfaceline.SortPointsByXAscending(); + } + + private static void BuildGeometryModel(GeometryPointString surfaceLine, SoilProfile2D soilProfile2D, ref SoilProfile2D result) + { + SoilProfile2D clonedProfile2 = (SoilProfile2D) soilProfile2D.Clone(); + result.Geometry = clonedProfile2.Geometry; + if (result.Geometry.Curves == null || result.Geometry.Curves.Count < 3) + return; + result.Geometry.DeleteLooseCurves(); + result.Geometry.Rebox(); + if (surfaceLine == null || surfaceLine.Points.Count <= 1) + return; + AdaptGeometryToSurfaceLineLimits(surfaceLine, ref result); + AddLeadingSurfaceLineToGeometry(surfaceLine, ref result); + result.Geometry.NewlyEffectedPoints.AddRange((IEnumerable) result.Geometry.Points); + result.Geometry.NewlyEffectedCurves.AddRange((IEnumerable) result.Geometry.Curves); + result.Geometry.Surfaces.Clear(); + result.Geometry.RegenerateGeometry(); + result.Geometry.DeleteLoosePoints(); + result.Geometry.DeleteLooseCurves(); + } + + private static void AdaptGeometryToSurfaceLineLimits(GeometryPointString surfaceLine, ref SoilProfile2D result) + { + GeometryBounds geometryBounds = surfaceLine.GetGeometryBounds(); + if (result.Geometry.Left > geometryBounds.Left) + GeometryHelper.ExtendGeometryLeft(result.Geometry, geometryBounds.Left); + else if (result.Geometry.Left < geometryBounds.Left) + GeometryHelper.CutGeometryLeft(result.Geometry, geometryBounds.Left); + if (result.Geometry.Right > geometryBounds.Right) + { + GeometryHelper.CutGeometryRight(result.Geometry, geometryBounds.Right); + } + else + { + if (result.Geometry.Right < geometryBounds.Right) + GeometryHelper.ExtendGeometryRight(result.Geometry, geometryBounds.Right); + } + } + + private static void AddLeadingSurfaceLineToGeometry(GeometryPointString surfaceLine, ref SoilProfile2D result) + { + surfaceLine.SyncCalcPoints(); + IList points = surfaceLine.CalcPoints; + GeometryData geometry = result.Geometry; + Point2D current1 = points[0].Clone(); + Point2D startPoint1 = new Point2D(current1.X, current1.Z); + Point2D current2 = HandleConnectionPointOfSurfaceLineToGeometry(ref result, startPoint1, current1, geometry, true); + int index = 1; + while (index < points.Count) + { + Point2D headPoint = current2; + current2 = points[index].Clone(); + if (index == checked (points.Count - 1)) + { + Point2D startPoint2 = new Point2D(current2.X, current2.Z); + current2 = HandleConnectionPointOfSurfaceLineToGeometry(ref result, startPoint2, current2, geometry, false); + } + else + { + Point2D point3D = new Point2D(current2.X, current2.Z); + Point2D point = result.Geometry.GetPoint(point3D, 0.001); + if (point != null) + current2 = point; + else + geometry.Points.Add(current2); + } + Point2D endPoint = current2; + GeometryCurve newCurve = new GeometryCurve(headPoint, endPoint); + if (!DoesCurveExist(result.Geometry, newCurve)) + result.Geometry.Curves.Add(newCurve); + checked { ++index; } + } + } + + private static bool DoesCurveExist(GeometryData geometry, GeometryCurve newCurve) + { + foreach (GeometryCurve curve in geometry.Curves) + { + if (curve.HeadPoint == newCurve.HeadPoint && curve.EndPoint == newCurve.EndPoint || + curve.HeadPoint == newCurve.EndPoint && curve.EndPoint == newCurve.HeadPoint) + return true; + } + return false; + } + + private static Point2D HandleConnectionPointOfSurfaceLineToGeometry(ref SoilProfile2D result, Point2D startPoint, + Point2D current, GeometryData data, bool forLeftSide) + { + Point2D point = result.Geometry.GetPoint(startPoint, 0.001); + if (point != null) + { + current = point; + } + else + { + List aCurveList = new List(); + result.Geometry.GetCurvesCoincidingInputPoint(startPoint, ref aCurveList); + if (aCurveList.Count == 0) + { + Point2D headPoint = !forLeftSide ? result.Geometry.GetRightPoints().OrderBy((Func) (a => a.Z)).Last() : + result.Geometry.GetLeftPoints().OrderBy((Func) + (a => a.Z)).Last(); + data.Points.Add(current); + GeometryCurve geometryCurve = new GeometryCurve(headPoint, current); + data.Curves.Add(geometryCurve); + } + else if (aCurveList.Count == 1) + SplitCurveAtPoint(data, aCurveList[0], current); + } + return current; + } + + private static void SplitCurveAtPoint(GeometryData geometry, GeometryCurve geometryCurve, Point2D splitPoint) + { + Point2D point = geometry.GetPoint(new Point2D(splitPoint.X, splitPoint.Z), 0.001); + if (point == null) + geometry.Points.Add(splitPoint); + else + splitPoint = point; + GeometryCurve geometryCurve1 = geometryCurve; + Point2D endPoint = geometryCurve1.EndPoint; + geometryCurve1.EndPoint = splitPoint; + GeometryCurve geometryCurve2 = new GeometryCurve(splitPoint, endPoint); + geometry.Curves.Add(geometryCurve2); + } + + private static void RemoveGeometryDataOfSoilProfileAboveSurfaceLine(GeometryPointString surfaceLine, ref SoilProfile2D result) + { + GeometryData geometry = result.Geometry; + RemoveGeometryDataAboveSurfaceLine(surfaceLine, ref geometry); + result.Geometry = geometry; + } + + private static void RemoveGeometryDataAboveSurfaceLine(GeometryPointString surfaceLine, ref GeometryData result) + { + if (surfaceLine == null || surfaceLine.Points.Count <= 1) + return; + foreach (GeometrySurface geometrySurface in result.Surfaces.ToArray()) + { + if (surfaceLine.GetPosition(geometrySurface.OuterLoop) == RelativeXzPosition.AboveGeometricLine) + { + result.Surfaces.Remove(geometrySurface); + foreach (GeometryCurve curve in geometrySurface.OuterLoop.CurveList) + { + if (curve.SurfaceAtLeft == geometrySurface) + curve.SurfaceAtLeft = null; + if (curve.SurfaceAtRight == geometrySurface) + curve.SurfaceAtRight = null; + } + foreach (GeometryLoop geometryLoop in result.Loops.ToArray()) + { + GeometryLoop loop = geometryLoop; + if (result.Surfaces.Count((Func) (s => s.OuterLoop == loop)) <= 0 && + result.Surfaces.Count((Func) (s => s.InnerLoops.Contains(loop))) <= 0) + { + result.Loops.Remove(loop); + } + } + } + } + //result.geom + result.DeleteLooseCurves(); + result.DeleteLoosePoints(); + } + + private static void ReconstructSurfaces(ref SoilProfile2D result, SoilProfile2D soilProfile2D, List oldSurfaces, + Soil defaultSoil) + { + result.Surfaces.Clear(); + foreach (GeometrySurface surface in result.Geometry.Surfaces) + { + SoilLayer2D soilLayer2D1 = new SoilLayer2D(); + soilLayer2D1.GeometrySurface = surface; + soilLayer2D1.Soil = defaultSoil; + SoilLayer2D soilLayer2D2 = soilLayer2D1; + SoilLayer2D layerFromOldSurfaces = soilProfile2D.GetOriginalLayerFromOldSurfaces(surface, (IEnumerable) oldSurfaces); + if (layerFromOldSurfaces != null && layerFromOldSurfaces.Soil != null) + { + soilLayer2D2.Soil = layerFromOldSurfaces.Soil; + soilLayer2D2.IsAquifer = layerFromOldSurfaces.IsAquifer; + soilLayer2D2.WaterpressureInterpolationModel = layerFromOldSurfaces.WaterpressureInterpolationModel; + //soilLayer2D2.Color = layerFromOldSurfaces.Color; + //soilLayer2D2.Description = layerFromOldSurfaces.Description; + } + + result.Surfaces.Add(soilLayer2D2); + } + } + + private static void ReconstructPreConsolidations(ref SoilProfile2D result, SoilProfile2D cloneProfile, double shift) + { + foreach (PreConsolidationStress preconsolidationStress in cloneProfile.PreconsolidationStresses) + { + preconsolidationStress.X += shift; + result.PreconsolidationStresses.Add(preconsolidationStress); + } + } } \ No newline at end of file Index: DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/Point2D.cs =================================================================== diff -u -r5065 -r5069 --- DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/Point2D.cs (.../Point2D.cs) (revision 5065) +++ DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/Point2D.cs (.../Point2D.cs) (revision 5069) @@ -117,8 +117,19 @@ X = aPoint.X; Z = aPoint.Z; } - + /// + /// Rounds the value using the GeometryConstants.Accuracy + /// + /// + /// + public double RoundValue(double value) + { + var res = Math.Round(value / GeometryConstants.Accuracy); + return res * GeometryConstants.Accuracy; + } + + /// /// The name for a point (must be here for the IGeometryObject interface) /// public string Name { get; set; }