Index: DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometryGenerator.cs =================================================================== diff -u --- DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometryGenerator.cs (revision 0) +++ DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometryGenerator.cs (revision 3079) @@ -0,0 +1,933 @@ +// Copyright (C) Stichting Deltares 2019. All rights reserved. +// +// This file is part of the Dam Engine. +// +// The Dam Engine is free software: you can redistribute it and/or modify +// it under the terms of the GNU Affero General Public License as published by +// the Free Software Foundation, either version 3 of the License, or +// (at your option) any later version. +// +// This program is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +// GNU Affero General Public License for more details. +// +// You should have received a copy of the GNU Affero General Public License +// along with this program. If not, see . +// +// All names, logos, and references to "Deltares" are registered trademarks of +// Stichting Deltares and remain full property of Stichting Deltares at all times. +// All rights reserved. + +using System; +using System.Collections.Generic; +using System.Diagnostics; + +namespace Deltares.DamEngine.Data.Geometry +{ + // ---------------------------------------------------------- + // The main geometry regeneration manager + public class GeometryGenerator + { + // aValue1 list that collects all surfaces that are at the right and left of effected curves (above) + // will be used to re-assign after regeneration + + private readonly Dictionary geometryCurveForwardsIsUsed = new Dictionary(); + private readonly Dictionary geometryCurveReversedIsUsed = new Dictionary(); + private readonly GeometryData geometryData; + private readonly Dictionary> geometryLoopDirections = new Dictionary>(); +// private readonly List intersectedCurveList = new List(); + private readonly List newlyDetectedSurfaceList = new List(); + + /// + /// Regenerates the geometry. + /// + public GeometryGenerator(GeometryData aGeometryData) + { + geometryData = aGeometryData; + } + + /// + /// Regenerates the geometry. + /// + public void GenerateGeometry() + { + FilterOutDoublePoints(); + + SetupCurveSurfaceAssociations(); + + bool firstRegeneration = true; + + while (true) + { + // break up all curves at intersections + //intersectedCurveList.Clear(); + //RegenerateAllCurvesIntersection(); Only add if really needed + + newlyDetectedSurfaceList.Clear(); + geometryLoopDirections.Clear(); + + int curvesCount = geometryData.Curves.Count; + + // initialise IsUsed of all curves to false + for (int index1 = 0; index1 < curvesCount; index1++) + { + SetIsUsed(geometryData.Curves[index1], CurveDirection.Forward, false); + SetIsUsed(geometryData.Curves[index1], CurveDirection.Reverse, false); + } + + // detect surfaces... the plaxis algorithm + int result = DetectSurfaces(firstRegeneration); + + // clear the effected & intersected curve list + //intersectedCurveList.Clear(); only when really needed + + if (result < 0) + { + if (firstRegeneration == false) + { + break; + } + } + else + { + break; + } + firstRegeneration = false; + } + } + + /// + /// Adds the curve to the list in given direction. + /// + /// + /// + /// + /// + private bool AddCurve(GeometryCurve aCurve, CurveDirection aDirection, GeometryLoop aEditedLoop) + { + if (FindCurveIndex(aEditedLoop, aCurve, aDirection) > -1) + { + return false; + } + aEditedLoop.CurveList.Add(aCurve); + if (!geometryLoopDirections.ContainsKey(aEditedLoop)) + { + geometryLoopDirections.Add(aEditedLoop, new List()); + } + geometryLoopDirections[aEditedLoop].Add(aDirection); + return true; + } + + // ------------------------------------------------------------------------------------- + // Loop Detection (Surface Creation) Algorithm from Plaxis + + /// + /// Creates the surface. + /// + /// a loop. + /// + private GeometrySurface CreateSurface(GeometryLoop aLoop) + { + if (aLoop.IsLoop() && aLoop.HasArea()) + { + var newSurface = new GeometrySurface(); + //newSurface.SetOuterLoop(aLoop); #Bka: replaced with + newSurface.OuterLoop = aLoop; + geometryData.Surfaces.Add(newSurface); + return newSurface; + } + return null; + } + + /// + /// Finds the index of the curve. + /// + /// a geometry loop. + /// a curve. + /// a direction. + /// + private int FindCurveIndex(GeometryLoop aGeometryLoop, GeometryCurve aCurve, CurveDirection aDirection) + { + if (aGeometryLoop.CurveList.Contains(aCurve) == false) + { + return -1; + } + + int curvesCount = aGeometryLoop.CurveList.Count; + + for (int index = 0; index < curvesCount; index++) + { + // Note Bka: Checking the direction does allow one curve to be added to ONE loop twice! + // This produces some strange surfaces (see LoopDetectionCase5) but that seems to be required + if (aGeometryLoop.CurveList[index] == aCurve && geometryLoopDirections[aGeometryLoop][index] == aDirection) + { + return index; + } + } + return -1; + } + + /// + /// + /// + /// + /// + /// + private void SetIsUsed(GeometryCurve aCurve, CurveDirection aDirection, bool aValue) + { + if (aDirection == CurveDirection.Forward) + { + if (!(geometryCurveForwardsIsUsed.ContainsKey(aCurve))) + { + geometryCurveForwardsIsUsed.Add(aCurve, aValue); + } + else + { + geometryCurveForwardsIsUsed[aCurve] = aValue; + } + } + else + { + if (!(geometryCurveReversedIsUsed.ContainsKey(aCurve))) + { + geometryCurveReversedIsUsed.Add(aCurve, aValue); + } + else + { + geometryCurveReversedIsUsed[aCurve] = aValue; + } + } + } + + /// + /// + /// + /// + /// + /// + private bool GetIsUsed(GeometryCurve aCurve, CurveDirection aDirection) + { + if (aDirection == CurveDirection.Forward) + { + return geometryCurveForwardsIsUsed[aCurve]; + } + else + { + return geometryCurveReversedIsUsed[aCurve]; + } + } + + /// + /// Removes all double points, replaces them in curves and newlyeffectedpoints. + /// Also deletes "zero" curves (beginpoint = endpoint) also from surface loops and newlyeffectedcurrves + /// + private void FilterOutDoublePoints() + { + var doublePoints = new Dictionary(); + + // Make sure all points (as pointers) in curves are in the point list + foreach (var curve in geometryData.Curves) + { + if (!geometryData.Points.Contains(curve.HeadPoint)) + { + geometryData.Points.Add(curve.HeadPoint); + } + + if (!geometryData.Points.Contains(curve.EndPoint)) + { + geometryData.Points.Add(curve.EndPoint); + } + } + + // find double points (by location!) in point list and register them with original + for (int i = 0; i < geometryData.Points.Count; i++) + { + for (int j = 0; j < i; j++) + { + if (i != j && geometryData.Points[i].LocationEquals(geometryData.Points[j])) + { + // register the double point and the original point + doublePoints[geometryData.Points[i]] = geometryData.Points[j]; + break; + } + } + } + + // replace double points in curves with originals + foreach (var doublePoint in doublePoints.Keys) + { + foreach (var curve in geometryData.Curves) + { + if (curve.HeadPoint == doublePoint) + { + curve.HeadPoint = doublePoints[doublePoint]; + } + + if (curve.EndPoint == doublePoint) + { + curve.EndPoint = doublePoints[doublePoint]; + } + } + } + + // remove curves which have the same head as end point + foreach (var curve in geometryData.Curves.ToArray()) + { + if (curve.HeadPoint == curve.EndPoint) + { + geometryData.Curves.Remove(curve); + if (geometryData.NewlyEffectedCurves.Contains(curve)) + { + geometryData.NewlyEffectedCurves.Remove(curve); + } + + foreach (var surface in geometryData.Surfaces) + { + surface.OuterLoop.CurveList.Remove(curve); + foreach (var loop in surface.InnerLoops) + { + loop.CurveList.Remove(curve); + } + } + } + } + + // removing curves from loops in surfaces may have created invalid surfaces, remove those here + foreach (var surface in geometryData.Surfaces.ToArray()) + { + if (!surface.OuterLoop.HasArea()) + { + geometryData.Surfaces.Remove(surface); + } + } + // remove double points from point list + foreach (var point in doublePoints.Keys) + { + geometryData.Points.Remove(point); + + if (geometryData.NewlyEffectedPoints.Contains(point)) + { + geometryData.NewlyEffectedPoints.Remove(point); + if (!geometryData.NewlyEffectedPoints.Contains(doublePoints[point])) + { + geometryData.NewlyEffectedPoints.Add(doublePoints[point]); + } + } + } + } + + private int DetectSurfaces(bool aFirstRegeneration) + { + try + { + // declare some variables + int curvesCount = geometryData.Curves.Count; + + var newLoopList = new List(); + var attachedCurveList = new List(); + + // start the first iteration + for (int index = 0; index < curvesCount * 2; index++) + { + var curveIndex = index / 2; + // look for current curve + var currentCurve = geometryData.Curves[curveIndex]; + + if (currentCurve == null) + { + continue; + } + // check the direction + CurveDirection currentCurveDirection; + if (index % 2 == 0) + { + if (GetIsUsed(currentCurve, CurveDirection.Forward)) + { + continue; + } + currentCurveDirection = CurveDirection.Forward; // get the direction of the current curve + } + else + { + if (GetIsUsed(currentCurve, CurveDirection.Reverse)) + { + continue; + } + currentCurveDirection = CurveDirection.Reverse; // get the direction of the current curve + } + + // create aValue1 new loop + var newLoop = new GeometryLoop(); + + //this.geometryData.Loops.Add(newLoop); + newLoopList.Add(newLoop); + + // initialise LoopBeginCurve + var loopBeginCurve = geometryData.Curves[curveIndex]; + var loopBeginDirection = currentCurveDirection; + + while (true) + { + // set the IsUsed status + SetIsUsed(currentCurve, currentCurveDirection, true); + + // add the current curve to new loop + if (AddCurve(currentCurve, currentCurveDirection, newLoop) == false) + { + // the curve wasn't added bcos the curve-direction pair was already present in loop. + // problem case - break here, else we'd get aValue1 hang! + // Todo: Solve this problem + if (aFirstRegeneration == false) + { + //TODO:Show error message box + break; + } + else + { + return -1; + } + } + + var curveEndPoint = currentCurve.GetEndPoint(currentCurveDirection); + attachedCurveList.Clear(); + var minAngle = 365.0; + var minIndex = 0; + + // find all the curves that are connected to the Current Curve at curveEndPoint + for (int index2 = 0; index2 < curvesCount; index2++) + { + var curve = geometryData.Curves[index2]; + if (curve.LocationEquals(currentCurve)) // lets not get the reverse direction of the current curve here + { + continue; + } + + if (curve.HeadPoint == curveEndPoint) + { + attachedCurveList.Add(new DirectionCurve(curve, CurveDirection.Forward)); + } + else if (curve.EndPoint == curveEndPoint) + { + attachedCurveList.Add(new DirectionCurve(curve, CurveDirection.Reverse)); + } + } + + if (attachedCurveList.Count == 0) // no curves found + { + CurveDirection oppCurrentDirection = currentCurveDirection == CurveDirection.Forward ? CurveDirection.Reverse : CurveDirection.Forward; + // if the current curve is not used in the opposite direction, it is considered in the opposite direction + if (GetIsUsed(currentCurve, oppCurrentDirection) == false) + { + currentCurveDirection = oppCurrentDirection; + continue; + } + else + { + break; + } + } + + // we have aValue1 set of curves, find the one that turns right the most + if (attachedCurveList.Count > 1) + { + minIndex = -1; + + var point1 = currentCurve.GetEndPoint(currentCurveDirection); + var point2 = currentCurve.GetHeadPoint(currentCurveDirection); + + for (int index2 = 0; index2 < attachedCurveList.Count; index2++) + { + var point3 = new Point2D(); + var point4 = new Point2D(); + + point3.X = attachedCurveList[index2].GetHeadPoint().X; + point3.Z = attachedCurveList[index2].GetHeadPoint().Z; + point4.X = attachedCurveList[index2].GetEndPoint().X; + point4.Z = attachedCurveList[index2].GetEndPoint().Z; + + var angle = Routines2D.FindAngle(point1, point2, point3, point4); + + if (angle < minAngle) + { + minAngle = angle; + minIndex = index2; + } + } + } + + DirectionCurve pickedDirectionCurve = attachedCurveList[minIndex]; + + if (pickedDirectionCurve.Curve == loopBeginCurve && pickedDirectionCurve.Direction == loopBeginDirection) + { + break; + } + + // assign the CurrentCurve from the picked one + currentCurve = pickedDirectionCurve.Curve; + currentCurveDirection = pickedDirectionCurve.Direction; + } + } + + // create surfaces! + return CreateSurfaces(newLoopList); + } +#if DEBUG + catch //(Exception ex) + { + //Todo:Show error msg box + //var lm = new LogMessage(LogMessageType.FatalError, this, "DetectSurfaces:" + ex.Message); + //LogManager.Add(lm); + return 0; + } +#else + catch + { + return 0; + } +#endif + } + + private int CreateSurfaces(List aNewLoopList) + { + var newSurfacesGeoDtaObjectList = new List(); + + int loopsCount = aNewLoopList.Count; + int curvesCount; + GeometrySurface newSurface; + var newSurfaceList = new List(); + + for (int index = 0; index < loopsCount; index++) + { + var loop = aNewLoopList[index]; + curvesCount = loop.CurveList.Count; + + if (curvesCount < 2) // dont create aValue1 surface for loops that have less than 2 curves + { + continue; + } + if (curvesCount == 2) // if only 2 curves in loop, make sure they are distinct (non-repeated) + { + if (loop.CurveList[0] == loop.CurveList[1]) + { + continue; + } + } + + if (!loop.IsLoop()) + { + continue; + } + + // if the loop is clockwise, create aValue1 surface + if (loop.IsClockWise() && !CheckIfLoopEnclosesOpenPolyline(loop)) + { + // TEMP: aVector1 mechanism to remember surfaces after Geometry Regeneration + // find the surface that is most repeated in the loop's in curves + newSurface = GetReassignmentSurfaceFromCurves(loop); + + // an existing surface has been found + if (newSurface != null) + { + newSurface = CreateSurface(loop); + } + else // no existing surface found from its comprising curves... create aValue1 brand new surface! + { + newSurface = CreateSurface(loop); + newlyDetectedSurfaceList.Add(newSurface); // populate the newly detected surface list + newSurfacesGeoDtaObjectList.Add(newSurface); + } + + if (newSurface != null) + { + newSurfaceList.Add(newSurface); + AssignSurfaceAtLeftOrRightToCurves(newSurface); + } + } + } + + // clear the left and right surfaces for all curves (some curves will have redundant data) + curvesCount = geometryData.Curves.Count; + for (int index = 0; index < curvesCount; index++) + { + geometryData.Curves[index].SurfaceAtRight = null; + geometryData.Curves[index].SurfaceAtLeft = null; + } + + // for the new surfaces -- assign the left/right surfaces for comprising curves, and find inner loops + var surfacesCount = newSurfaceList.Count; + for (int index = 0; index < surfacesCount; index++) + { + newSurface = newSurfaceList[index]; + AssignSurfaceAtLeftOrRightToCurves(newSurface); + object newSurfaceObject = newSurface /*new object()*/; + CheckAndAddInnerLoops(ref newSurfaceObject); + //newSurface = (GeometrySurface)newSurfaceObject; + } + if (newSurfacesGeoDtaObjectList.Count > 0) + { + var lNewSurfaces = new List(); + GetNewlyDetectedSurfaces(ref lNewSurfaces); + } + return surfacesCount; + } + + /// + /// Checks if loop encloses open polyline. + /// + /// a loop. + /// + private bool CheckIfLoopEnclosesOpenPolyline(GeometryLoop aLoop) + { + int curvesCount = aLoop.CurveList.Count; + if (curvesCount < 3) + { + return true; + } + + for (int index = 0; index < curvesCount; index++) + { + var curve = aLoop.CurveList[index]; + var direction = geometryLoopDirections[aLoop][index]; + var foundOppDirection = false; + + for (int index1 = 0; index1 < curvesCount; index1++) + { + if (index == index1) + { + continue; + } + + if (aLoop.CurveList[index1] == curve && geometryLoopDirections[aLoop][index] != direction) + { + foundOppDirection = true; + break; + } + } + + if (foundOppDirection == false) + { + return false; + } + } + return true; + } + + private GeometrySurface GetReassignmentSurfaceFromCurves(GeometryLoop aLoop) + { + GeometrySurface surface = null; + GeometrySurface reassignmentSurface = null; + int curvesCount = aLoop.CurveList.Count; + + if (!geometryLoopDirections.ContainsKey(aLoop)) + { + return null; + } + + for (int index = 0; index < curvesCount; index++) + { + var curve = aLoop.CurveList[index]; + if (geometryLoopDirections[aLoop][index] == CurveDirection.Forward) + { + if (curve.SurfaceAtRight != null) + { + surface = curve.SurfaceAtRight; + } + } + else + { + if (curve.SurfaceAtLeft != null) + { + surface = curve.SurfaceAtLeft; + } + } + + if (surface == null) + { + continue; + } + + var maxTimesSurfacesFound = 0; + var noTimesSurfaceFound = 0; + for (int index1 = 0; index1 < curvesCount; index1++) + { + if (geometryLoopDirections[aLoop][index] == CurveDirection.Forward) + { + if (curve.SurfaceAtRight == surface) + { + noTimesSurfaceFound++; + } + } + else + { + if (curve.SurfaceAtLeft == surface) + { + noTimesSurfaceFound++; + } + } + } + + if (noTimesSurfaceFound > maxTimesSurfacesFound) + { + maxTimesSurfacesFound = noTimesSurfaceFound; + reassignmentSurface = surface; + } + } + return reassignmentSurface; + } + + /// + /// Assigns the surface at left or right to its curves on the outerloop. This tells at which side of the curve the surface is present. + /// So after this, for each curve in the outerloop of this surface, it is known at which side the surface is located. + /// + /// The surface. + private void AssignSurfaceAtLeftOrRightToCurves(GeometrySurface surface) + { + GeometryLoop loop = surface.OuterLoop; + int curvesCount = loop.CurveList.Count; + bool isClockwise = true; + + try + { + isClockwise = loop.IsClockWise(); + } + catch (GeometryLoop.NotEnoughUniquePointsException e) + { + Debug.WriteLine(e.Message); + } + catch (InvalidOperationException e) + { + Debug.WriteLine(e.Message); + } + + for (int index = 0; index < curvesCount; index++) + { + if (isClockwise) + { + if (geometryLoopDirections[loop][index] == CurveDirection.Forward) + { + loop.CurveList[index].SurfaceAtRight = surface; + } + else + { + loop.CurveList[index].SurfaceAtLeft = surface; + } + } + else + { + if (geometryLoopDirections[loop][index] == CurveDirection.Forward) + { + loop.CurveList[index].SurfaceAtLeft = surface; + } + else + { + loop.CurveList[index].SurfaceAtRight = surface; + } + } + } + } + + private void SetupCurveSurfaceAssociation() + { + // clear the data + int count = geometryData.Curves.Count; + for (int i = 0; i < count; i++) + { + geometryData.Curves[i].SurfaceAtLeft = null; + geometryData.Curves[i].SurfaceAtRight = null; + } + + // reset + count = geometryData.Surfaces.Count; + for (int i = 0; i < count; i++) + { + AssignSurfaceAtLeftOrRightToCurves(geometryData.Surfaces[i]); + } + } + + /// + /// Setups the curve surface associations. + /// + private void SetupCurveSurfaceAssociations() + { + SetUpGeometryLoopDirections(); + // only try to connect curves to surfaces when there are loops (i.e. surfaces) + if (geometryLoopDirections.Count > 0) + { + SetupCurveSurfaceAssociation(); + } + } + + private void SetUpGeometryLoopDirections() + { + geometryLoopDirections.Clear(); + foreach (var loop in geometryData.Loops) + { + if (!geometryLoopDirections.ContainsKey(loop)) + { + SetUpGeometryLoopDirections(loop); + } + } + } + + private void SetUpGeometryLoopDirections(GeometryLoop aLoop) + { + if (aLoop.CurveList.Count > 0) + { + Point2D loopPoint; + + geometryLoopDirections.Add(aLoop, new List()); + // get the first curve + if (aLoop.CurveList[0].EndPoint == aLoop.CurveList[1].HeadPoint || aLoop.CurveList[0].EndPoint == aLoop.CurveList[1].EndPoint) + { + geometryLoopDirections[aLoop].Add(CurveDirection.Forward); + loopPoint = aLoop.CurveList[0].EndPoint; + } + else + { + geometryLoopDirections[aLoop].Add(CurveDirection.Reverse); + loopPoint = aLoop.CurveList[0].HeadPoint; + } + + // the rest of the curves + for (int index1 = 1; index1 < aLoop.CurveList.Count; index1++) + { + if (loopPoint == aLoop.CurveList[index1].HeadPoint) + { + geometryLoopDirections[aLoop].Add(CurveDirection.Forward); + loopPoint = aLoop.CurveList[index1].EndPoint; + } + else + { + geometryLoopDirections[aLoop].Add(CurveDirection.Reverse); + loopPoint = aLoop.CurveList[index1].HeadPoint; + } + } + } + } + + /// + /// Checks and adds inner loop to the new surface. + /// + /// + private void CheckAndAddInnerLoops(ref object aNewSurface) + { + var newSurface = (GeometrySurface)aNewSurface; + int surfaceCount = geometryData.Surfaces.Count; + var newLoop = newSurface.OuterLoop; + int newPointCount = newLoop.CurveList.Count; + List newPolygon = newLoop.CalcPoints; + + newSurface.InnerLoops.Clear(); + + for (int index = 0; index < surfaceCount; index++) + { + if (newSurface == geometryData.Surfaces[index]) + { + continue; + } + + var innerPointCount = 0; + var outerPointCount = 0; + var isOnPointCount = 0; + var hasOnPointCount = 0; + var loop = geometryData.Surfaces[index].OuterLoop; + var polygon = loop.CalcPoints; + var existingLoopPointCount = polygon.Count; + + // check if it is an inner loop + for (int innerIndex = 0; innerIndex < newPointCount; innerIndex++) + { + PointInPolygon location = Routines2D.CheckIfPointIsInPolygon(loop, newPolygon[innerIndex].X, newPolygon[innerIndex].Z); + + if (location == PointInPolygon.InsidePolygon) + { + innerPointCount++; + } + else if (location == PointInPolygon.OnPolygonEdge) + { + isOnPointCount++; + } + } + + // check if it has an inner loop + for (int innerIndex1 = 0; innerIndex1 < existingLoopPointCount; innerIndex1++) + { + PointInPolygon location = Routines2D.CheckIfPointIsInPolygon(newLoop, polygon[innerIndex1].X, polygon[innerIndex1].Z); + + if (location == PointInPolygon.InsidePolygon) + { + outerPointCount++; + } + else if (location == PointInPolygon.OnPolygonEdge) + { + hasOnPointCount++; + } + } + + //Add New Loop as inner loop to the existing Surface + if ((innerPointCount == newPointCount) || ((innerPointCount > 0) && (newPointCount == (innerPointCount + isOnPointCount)))) + { + geometryData.Surfaces[index].AddInnerLoop(newLoop); + } + + //Add Inner Loop to the New Surface + if ((outerPointCount == existingLoopPointCount) || ((outerPointCount > 0) && (existingLoopPointCount == (outerPointCount + hasOnPointCount)))) + { + newSurface.AddInnerLoop(loop); + } + } + } + + /// + /// Gets the newly detected surface from the list. + /// + /// + private void GetNewlyDetectedSurfaces(ref List aNewSurfaceList) + { + aNewSurfaceList.AddRange(newlyDetectedSurfaceList); + } + + #region Nested type: DirectionCurve + + internal struct DirectionCurve + { + private readonly GeometryCurve curve; + private readonly CurveDirection direction; + + internal DirectionCurve(GeometryCurve aCurve, CurveDirection aDirection) + { + curve = aCurve; + direction = aDirection; + } + + internal GeometryCurve Curve + { + get + { + return curve; + } + } + + internal CurveDirection Direction + { + get + { + return direction; + } + } + + internal Point2D GetHeadPoint() + { + return curve.GetHeadPoint(direction); + } + + internal Point2D GetEndPoint() + { + return curve.GetEndPoint(direction); + } + } + + #endregion + } +} Index: DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/LineHelper.cs =================================================================== diff -u -r1974 -r3079 --- DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/LineHelper.cs (.../LineHelper.cs) (revision 1974) +++ DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/LineHelper.cs (.../LineHelper.cs) (revision 3079) @@ -153,8 +153,7 @@ var point3 = new Point2D(p3.X, p3.Z); var point4 = new Point2D(p4.X, p4.Z); - var ip = new Point2D(); - var res = Routines2D.DetermineIf2DLinesIntersectWithExtrapolation(point1, point2, point3, point4, ref ip); + Routines2D.DetermineIf2DLinesIntersectWithExtrapolation(point1, point2, point3, point4, out var ip); if (ip != null) { intersectPoint = new GeometryPoint(ip.X, ip.Z); Index: DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometryPointString.cs =================================================================== diff -u -r1974 -r3079 --- DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometryPointString.cs (.../GeometryPointString.cs) (revision 1974) +++ DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometryPointString.cs (.../GeometryPointString.cs) (revision 3079) @@ -591,6 +591,23 @@ } /// + /// Removes the non ascending points on surface line. + /// + public void RemoveNonAscendingPointsOnSurfaceLine() + { + for (int i = points.Count - 1; i > 0; i--) + { + if (Math.Abs(points[i].X - points[i - 1].X) < GeometryConstants.Accuracy) + { + if (Math.Abs(points[i].Z - points[i - 1].Z) > GeometryConstants.Accuracy) + { + points.Remove(points[i]); + } + } + } + } + + /// /// Gets a part of the geometry point string defined by a begin and end point /// /// Index: DamEngine/trunk/src/Deltares.DamEngine.Data.Tests/Deltares.DamEngine.Data.Tests.csproj =================================================================== diff -u -r2193 -r3079 --- DamEngine/trunk/src/Deltares.DamEngine.Data.Tests/Deltares.DamEngine.Data.Tests.csproj (.../Deltares.DamEngine.Data.Tests.csproj) (revision 2193) +++ DamEngine/trunk/src/Deltares.DamEngine.Data.Tests/Deltares.DamEngine.Data.Tests.csproj (.../Deltares.DamEngine.Data.Tests.csproj) (revision 3079) @@ -51,8 +51,10 @@ + + @@ -66,6 +68,10 @@ {b7a49c1a-1c91-4d72-aba9-9fbac2509d8e} Deltares.DamEngine.Data + + {a41683f4-c3e0-4386-aee7-faacd535cf93} + Deltares.DamEngine.TestHelpers + Index: DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometryData.cs =================================================================== diff -u -r1974 -r3079 --- DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometryData.cs (.../GeometryData.cs) (revision 1974) +++ DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometryData.cs (.../GeometryData.cs) (revision 3079) @@ -33,6 +33,11 @@ /// public class GeometryData : GeometryObject { + private GeometryGenerator geometryGenerator; + private readonly List newlyEffectedCurves = new List(); + private readonly List newlyEffectedPoints = new List(); + private bool isRegeneratingGeometry; + private readonly List curveDataList = new List(); private readonly List loopDataList = new List(); @@ -47,6 +52,14 @@ private bool updatingSurfaceLine; + /// + /// Initializes a new instance of the class. + /// + public GeometryData() + { + geometryGenerator = null; + } + #region properties /// @@ -56,7 +69,7 @@ /// The points. /// [Validate] - public IList Points + public List Points { get { @@ -65,6 +78,20 @@ } /// + /// Gets the newly effected points. + /// + /// + /// The newly effected points. + /// + public List NewlyEffectedPoints + { + get + { + return newlyEffectedPoints; + } + } + + /// /// gets the Curve data list. /// /// @@ -79,6 +106,20 @@ } /// + /// Gets the newly effected curves. + /// + /// + /// The newly effected curves. + /// + public List NewlyEffectedCurves + { + get + { + return newlyEffectedCurves; + } + } + + /// /// gets the Loop data list. /// /// @@ -103,6 +144,35 @@ } } + public void RegenerateGeometry() + { + if (isRegeneratingGeometry) + { + return; + } + + isRegeneratingGeometry = true; + if (geometryGenerator == null) + { + geometryGenerator = new GeometryGenerator(this); + } + lock (this) + { + SynchronizeLoops(); + RemoveDoublesFromNewlyEffectedPointsAndCurves(); + Curves.Clear(); + Points.Clear(); + Points.AddRange(newlyEffectedPoints); + Curves.AddRange(newlyEffectedCurves); + + geometryGenerator.GenerateGeometry(); + newlyEffectedPoints.Clear(); + newlyEffectedCurves.Clear(); + UpdateSurfaceLine(); + } + isRegeneratingGeometry = false; + } + /// /// Gets the minimum geometry points x. /// @@ -225,6 +295,64 @@ } } + public void RemoveDoublesFromNewlyEffectedPointsAndCurves() + { + var pdel = new List(); + var par = newlyEffectedPoints.ToArray(); + for (int i = 0; i < par.Length; i++) + { + for (int j = i; j < par.Length; j++) + { + if (i != j && par[i].LocationEquals(par[j])) + { + if (!pdel.Contains(par[j])) + { + pdel.Add(par[j]); + } + } + } + } + foreach (var point in pdel) + { + newlyEffectedPoints.Remove(point); + } + + var cdel = new List(); + var car = newlyEffectedCurves.ToArray(); + // First remove all "illegal" newlyeffected curves + for (int i = 0; i < car.Length; i++) + { + if (car[i].HeadPoint == null || car[i].EndPoint == null || car[i].HeadPoint.LocationEquals(car[i].EndPoint)) + { + cdel.Add(car[i]); + } + } + foreach (var curve in cdel) + { + newlyEffectedCurves.Remove(curve); + } + // then remove all real doubles + var car2 = newlyEffectedCurves.ToArray(); + cdel.Clear(); + for (int i = 0; i < car2.Length; i++) + { + for (int j = i; j < car2.Length; j++) + { + if (i != j && car2[i].LocationEquals(car2[j])) + { + if (!cdel.Contains(car2[j])) + { + cdel.Add(car2[j]); + } + } + } + } + foreach (var curve in cdel) + { + newlyEffectedCurves.Remove(curve); + } + } + /// /// Gets all points on the Left boundary. /// @@ -295,10 +423,10 @@ pointDataList.Clear(); curveDataList.Clear(); surfaceDataList.Clear(); - -// newlyEffectedPoints.Clear(); -// newlyEffectedCurves.Clear(); + newlyEffectedPoints.Clear(); + newlyEffectedCurves.Clear(); } + /// /// deletes all the Loop from IGeometryLoop. /// @@ -338,7 +466,7 @@ } /// - /// CheckIfIntersectStricktly + /// CheckIfIntersect /// Determines if two lines intersect each other stricktly (so no extrapolated points). /// /// Line 1 GeometryPoint 1 @@ -543,27 +671,26 @@ } /// - /// Updates the line at the top of the geometry + /// Deletes the point and the curves it belongs too. /// - private void UpdateSurfaceLine() + /// The point. + public void DeletePointWithCurves(Point2D point) { - if (updatingSurfaceLine) + var curvesToDelete = new List(); + foreach (var curve in Curves) { - return; + if (curve.ContainsPoint(point)) + { + curvesToDelete.Add(curve); + } } - updatingSurfaceLine = true; - - var bCurves = GetBoundaryCurves(); - if (bCurves.Count == 0) + foreach (var curveToDelete in curvesToDelete) { - surfaceLine.CalcPoints.Clear(); + Curves.Remove(curveToDelete); } - var curvesCopy = GetCurvesCopy(bCurves); - var curves = GetTopCurves(curvesCopy); - CreateSurfaceLinePointString(curves); - updatingSurfaceLine = false; + Points.Remove(point); } /// @@ -789,5 +916,30 @@ return (curve1.HeadPoint == curve2.HeadPoint || curve1.HeadPoint == curve2.EndPoint || curve1.EndPoint == curve2.HeadPoint || curve1.EndPoint == curve2.EndPoint); } + + /// + /// Updates the line at the top of the geometry + /// + private void UpdateSurfaceLine() + { + if (updatingSurfaceLine) + { + return; + } + + updatingSurfaceLine = true; + + var bCurves = GetBoundaryCurves(); + if (bCurves.Count == 0) + { + surfaceLine.CalcPoints.Clear(); + } + var curvesCopy = GetCurvesCopy(bCurves); + var curves = GetTopCurves(curvesCopy); + CreateSurfaceLinePointString(curves); + + updatingSurfaceLine = false; + surfaceLine.SyncPoints(); + } } } \ No newline at end of file Index: DamEngine/trunk/src/Deltares.DamEngine.Calculators/KernelWrappers/MacroStabilityInwards/MacroStabilityInwardsKernelWrapper.cs =================================================================== diff -u -r3041 -r3079 --- DamEngine/trunk/src/Deltares.DamEngine.Calculators/KernelWrappers/MacroStabilityInwards/MacroStabilityInwardsKernelWrapper.cs (.../MacroStabilityInwardsKernelWrapper.cs) (revision 3041) +++ DamEngine/trunk/src/Deltares.DamEngine.Calculators/KernelWrappers/MacroStabilityInwards/MacroStabilityInwardsKernelWrapper.cs (.../MacroStabilityInwardsKernelWrapper.cs) (revision 3079) @@ -156,22 +156,16 @@ var soilProfile2D = subSoilScenario.SoilProfile2D; if (soilProfile2D == null) { - // var soilSurfaceProfile = new SoilSurfaceProfile - // { - // SoilProfile = subSoilScenario.SoilProfile1D, - // SurfaceLine2 = surfaceLine2, - // Name = surfaceLine2.Name + "_" + subSoilScenario.SoilProfile1D.Name, - // DikeEmbankmentMaterial = dikeEmbankmentSoil - // }; - // // Convert the soilsurfacesoilprofile to a SoilProfile2D to be able to edit it properly. - // var soilProfile2DNew = soilSurfaceProfile.ConvertToSoilProfile2D(); - // //soilSurfaceProfile.Dispose(); - // // For some obscure reason, the created soilProfile2D is handled wrong in the UI (see DSB=786). - // // Its curves do not seem to match the surfaces in the events send on selecting and dragging a curve. - // // This causes the strange behaviour. To solve this, a Clone is made of the soilProfile2D and that clone is added instead. - // //soilProfile2D = (SoilProfile2D)soilProfile2D.Clone(); - // //soilProfile2DNew.Dispose(); - // soilProfile2D = soilProfile2DNew; + var soilSurfaceProfile = new SoilSurfaceProfile + { + SoilProfile = subSoilScenario.SoilProfile1D, + SurfaceLine2 = surfaceLine2, + Name = surfaceLine2.Name + "_" + subSoilScenario.SoilProfile1D.Name, + DikeEmbankmentMaterial = dikeEmbankmentSoil + }; + // Convert the soilsurfacesoilprofile to a SoilProfile2D to be able to edit it properly. + var soilProfile2DNew = soilSurfaceProfile.ConvertToSoilProfile2D(); + soilProfile2D = soilProfile2DNew; } } Index: DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometrySurface.cs =================================================================== diff -u -r1974 -r3079 --- DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometrySurface.cs (.../GeometrySurface.cs) (revision 1974) +++ DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometrySurface.cs (.../GeometrySurface.cs) (revision 3079) @@ -76,6 +76,18 @@ } /// + /// Add to if it isn't already + /// in that collection. + /// + public void AddInnerLoop(GeometryLoop aLoop) + { + if (InnerLoops.Contains(aLoop) == false) + { + InnerLoops.Add(aLoop); + } + } + + /// /// determine top of geometry surface outerloop as GeometryPointString /// /// Index: DamEngine/trunk/src/Deltares.DamEngine.Data/Geotechnics/SoilProfile1D.cs =================================================================== diff -u -r1974 -r3079 --- DamEngine/trunk/src/Deltares.DamEngine.Data/Geotechnics/SoilProfile1D.cs (.../SoilProfile1D.cs) (revision 1974) +++ DamEngine/trunk/src/Deltares.DamEngine.Data/Geotechnics/SoilProfile1D.cs (.../SoilProfile1D.cs) (revision 3079) @@ -48,18 +48,18 @@ /// public class SoilProfile1D : SoilProfile { - private const double defaultBottomLayerHeight = 20.0; + private const double DefaultBottomLayerHeight = 20.0; private readonly DelegatedList layers = new DelegatedList(); private double bottomLevel = double.NaN; - private SoilLayer1D infiltrationLayer = null; + private SoilLayer1D infiltrationLayer; /// /// Default constructor /// public SoilProfile1D() { - Name = LocalizationManager.GetTranslatedText(this, "DefaultNameSoilProfile1D"); + base.Name = LocalizationManager.GetTranslatedText(this, "DefaultNameSoilProfile1D"); layers.AddMethod = AddLayer; } @@ -76,16 +76,33 @@ layers.Add(new SoilLayer1D(soil, topLevel)); } + /// + /// Assigns the specified profile. + /// + /// The profile. public void Assign(SoilProfile1D profile) { - this.Assign((SoilProfile)profile); - this.Layers.Clear(); - foreach (SoilLayer1D layer in (IEnumerable)profile.Layers) - this.Layers.Add((SoilLayer1D)layer.Clone()); - this.BottomLevel = profile.BottomLevel; + Assign((SoilProfile)profile); + Layers.Clear(); + foreach (SoilLayer1D layer in profile.Layers) + Layers.Add((SoilLayer1D)layer.Clone()); + BottomLevel = profile.BottomLevel; } /// + /// Creates a new object that is a copy of the current instance. + /// + /// + /// A new object that is a copy of this instance. + /// + public object Clone() + { + var cloneSoilProfile = new SoilProfile1D(); + cloneSoilProfile.Assign(this); + return cloneSoilProfile; + } + + /// /// Gets the soil layer collection for this profile /// public IList Layers @@ -122,7 +139,7 @@ { if (double.IsNaN(bottomLevel) && layers.Count > 0) { - bottomLevel = Layers.Last().TopLevel - defaultBottomLayerHeight; + bottomLevel = Layers.Last().TopLevel - DefaultBottomLayerHeight; } return bottomLevel; @@ -171,34 +188,15 @@ /// public SoilLayer1D GetLayerWithName(string name) { - return this.layers.FirstOrDefault((Func)(layer => + return layers.FirstOrDefault(layer => { if (layer.Name != null) return layer.Name.Equals(name); return false; - })); + }); } /// - /// This method makes sure all names are unique - /// - public void EnsureUniqueLayerNames() - { - if (AreAllLayerNamesUnique()) - { - return; - } - foreach (var layer in Layers) - { - layer.Name = null; - } - foreach (var layer in Layers) - { - layer.Name = GetNewUniqueLayerName(); - } - } - - /// /// Checks weather all layer names are unique. /// /// @@ -245,31 +243,6 @@ } /// - /// The highest aquifer layer, can be null if not aquifer is present - /// - [Validate] - public SoilLayer1D TopAquiferLayer - { - get - { - IList sortedLayers = Layers.OrderByDescending(l => l.TopLevel).ToList(); - SoilLayer1D aquiferLayer = null; - - // Search the highest aquifer layer - for (int layerIndex = 0; layerIndex < sortedLayers.Count; layerIndex++) - { - var layer = sortedLayers[layerIndex]; - if (IsAquiferLayer(layer)) - { - aquiferLayer = layer; - break; - } - } - return aquiferLayer; - } - } - - /// /// Gets the toppest aquifer layer of the deepest cluster of aquifers /// /// @@ -352,62 +325,6 @@ } /// - /// Gets the lowest layer in the highest cluster of in-between aquifers - /// - /// - /// The lowest layer in the highest cluster of in-between aquifers - /// - public SoilLayer1D BottomLayerOfInBetweenAquiferCluster - { - get - { - IList sortedLayers = Layers.OrderByDescending(l => l.TopLevel).ToList(); - SoilLayer1D highestAquiferLayer = null; - SoilLayer1D aquiferLayer = null; - int indexHighestAquifer = -1; - - // Search the index of the highest aquifer layer - for (int layerIndex = 1; layerIndex < sortedLayers.Count; layerIndex++) - { - var previousLayer = sortedLayers[layerIndex - 1]; - var layer = sortedLayers[layerIndex]; - if (IsAquiferLayer(layer) && !IsAquiferLayer(previousLayer)) - { - highestAquiferLayer = layer; - indexHighestAquifer = layerIndex; - break; - } - } - - // If highest aquifer layer is bottom aquifer layer, there is no in between aquiferlayer - if (highestAquiferLayer == BottomAquiferLayer) - { - return null; - } - - // in-between aquifers cluster may consists of more then 1 connected (aquifer) layers. - // Search all layers below the found highest aquifer to find bottom aquifer layer. - if (indexHighestAquifer >= 0) - { - for (int layerIndex = indexHighestAquifer; layerIndex < sortedLayers.Count; layerIndex++) - { - var layer = sortedLayers[layerIndex]; - - if (IsAquiferLayer(layer)) - { - aquiferLayer = layer; - } - else - { - break; - } - } - } - return aquiferLayer; - } - } - - /// /// Gets (calculates) the height for a given layer in the profile /// /// The layer to process @@ -428,7 +345,7 @@ var bottomLayer = Layers.Last(); if (bottomLayer.Height.IsZero()) { - BottomLevel -= defaultBottomLayerHeight; + BottomLevel -= DefaultBottomLayerHeight; } } @@ -669,6 +586,16 @@ } /// + /// Gets the layer at a given z coordinate. + /// + /// The z. + /// the found layer + public SoilLayer1D GetLayerAt(double z) + { + return Layers.FirstOrDefault(soilLayer => soilLayer.BottomLevel < z && soilLayer.TopLevel >= z); + } + + /// /// Gets the layer below the given layer. /// /// The given layer. @@ -719,78 +646,6 @@ } /// - /// Gets the total soil pressure, including the eventual free water on surface - /// - /// The level for wich the total pressure is calculated - /// The phreatic level - /// Unit weight of the water - /// The soil pressure - public double GetTotalPressure(double z, double phreaticLevel, double unitWeightWater) - { - double pressure = 0; - - // see if free water is at surface - double surfaceLevel = Layers[0].TopLevel; - if (z > surfaceLevel) - { - return Math.Max(0, (phreaticLevel - z)*unitWeightWater); - } - - if (phreaticLevel > surfaceLevel) - { - if (z < phreaticLevel) - { - pressure += (phreaticLevel - Math.Max(Layers[0].TopLevel, z))*unitWeightWater; - } - } - - foreach (var layer in Layers) - { - if (layer.Height > 0 && layer.Soil != null) - { - if (z >= layer.TopLevel) - { - continue; - } - if (z <= layer.BottomLevel) - { - if (phreaticLevel <= layer.BottomLevel) - { - pressure += layer.Height*layer.Soil.AbovePhreaticLevel; - } - else if (phreaticLevel >= layer.TopLevel) - { - pressure += layer.Height*layer.Soil.BelowPhreaticLevel; - } - else - { - pressure += (phreaticLevel - layer.BottomLevel)*layer.Soil.BelowPhreaticLevel; - pressure += (layer.TopLevel - phreaticLevel)*layer.Soil.AbovePhreaticLevel; - } - } - else - { - if (phreaticLevel <= z) - { - pressure += (layer.TopLevel - z)*layer.Soil.AbovePhreaticLevel; - } - else if (phreaticLevel >= layer.TopLevel) - { - pressure += (layer.TopLevel - z)*layer.Soil.BelowPhreaticLevel; - } - else - { - pressure += (layer.TopLevel - phreaticLevel)*layer.Soil.AbovePhreaticLevel; - pressure += (phreaticLevel - z)*layer.Soil.BelowPhreaticLevel; - } - } - } - } - - return pressure; - } - - /// /// Returns a that represents this instance. /// /// Index: DamEngine/trunk/src/Deltares.DamEngine.Data/Deltares.DamEngine.Data.csproj =================================================================== diff -u -r2963 -r3079 --- DamEngine/trunk/src/Deltares.DamEngine.Data/Deltares.DamEngine.Data.csproj (.../Deltares.DamEngine.Data.csproj) (revision 2963) +++ DamEngine/trunk/src/Deltares.DamEngine.Data/Deltares.DamEngine.Data.csproj (.../Deltares.DamEngine.Data.csproj) (revision 3079) @@ -110,6 +110,7 @@ + Index: DamEngine/trunk/src/Deltares.DamEngine.Data/Geotechnics/SoilSurfaceProfile.cs =================================================================== diff -u -r1974 -r3079 --- DamEngine/trunk/src/Deltares.DamEngine.Data/Geotechnics/SoilSurfaceProfile.cs (.../SoilSurfaceProfile.cs) (revision 1974) +++ DamEngine/trunk/src/Deltares.DamEngine.Data/Geotechnics/SoilSurfaceProfile.cs (.../SoilSurfaceProfile.cs) (revision 3079) @@ -21,6 +21,7 @@ using System; using System.Collections.Generic; +using System.Linq; using Deltares.DamEngine.Data.Geometry; using Deltares.DamEngine.Data.Standard; using Deltares.DamEngine.Data.Standard.Validation; @@ -36,43 +37,19 @@ public class SoilSurfaceProfile : SoilProfile2D { private Soil dikeEmbankmentMaterial; - private bool dirty; - private bool initial = true; - private SoilProfile1D soilProfile = null; - private SoilProfile1D orgSoilProfile; - private GeometryPointString surfaceLine; + private SoilProfile1D soilProfile; + private SoilProfile1D originalSoilProfile1D; private SurfaceLine2 surfaceLine2; - private bool createSettlementZones; - + /// /// Initializes a new instance of the class. /// public SoilSurfaceProfile() { - Name = "Soil Surface Profile"; + base.Name = "Soil Surface Profile"; } /// - /// Gets or sets the surface line. - /// - /// - /// The surface line. - /// - [Validate] - public GeometryPointString SurfaceLine - { - get - { - return surfaceLine; - } - set - { - surfaceLine = value; - dirty = true; - } - } - - /// /// Gets or sets the surface line2. /// [Validate] @@ -85,15 +62,6 @@ set { surfaceLine2 = value; - if (surfaceLine2 != null) - { - SurfaceLine = surfaceLine2.Geometry; - } - else - { - SurfaceLine = null; - } - dirty = true; } } @@ -108,12 +76,11 @@ { get { - return orgSoilProfile; + return originalSoilProfile1D; } set { - orgSoilProfile = value; - dirty = true; + originalSoilProfile1D = value; } } @@ -132,121 +99,10 @@ set { dikeEmbankmentMaterial = value; - dirty = true; } } /// - /// List of preconsolidation stresses - /// - public override List PreconsolidationStresses - { - get - { - if (orgSoilProfile != null) - { - return orgSoilProfile.PreconsolidationStresses; - } - return null; - } - } - - /// - /// Gets the surfaces. - /// - /// - /// The surfaces. - /// - [Validate] - public override IList Surfaces - { - get - { - if (initial || dirty) - { - UpdateLayers(); - } - return base.Surfaces; - } - } - - /// - /// Gets or sets a value indicating whether to [create settlement zones]. - /// - /// - /// true if [create settlement zones]; otherwise, false. - /// - public bool CreateSettlementZones - { - get - { - return createSettlementZones; - } - set - { - createSettlementZones = value; - } - } - - /// - /// Generates a 1D profile at a given x - /// - /// - /// Generated 1D soil profile - public override SoilProfile1D GetSoilProfile1D(double x) - { - if (initial || dirty) - { - initial = false; - - UpdateLayers(); - } - - // Try to obtain a cached soil profile (performance optimization) - var soilProfile1D = GetCachedSoilProfile1D(x); - if (soilProfile1D != null) - { - return soilProfile1D; - } - - // Otherwise, create and configure a new soil profile - soilProfile1D = new SoilProfile1D(); - - double top = SurfaceLine.GetZatX(x); - soilProfile1D.BottomLevel = soilProfile.BottomLevel; - - if (top > soilProfile.TopLevel) - { - var emptyTopLayer = new SoilLayer1D(null, top) - { - SoilProfile = soilProfile1D - }; - soilProfile1D.Layers.Add(emptyTopLayer); - } - - for (int i = 0; i < soilProfile.Layers.Count; i++) - { - SoilLayer1D layer = soilProfile.Layers[i]; - if (layer.BottomLevel < top) - { - // Perform without publishing events because otherwise the soil profile cache would get lost - var newLayer = new SoilLayer1D(layer.Soil, Math.Min(layer.TopLevel, top)) - { - IsAquifer = layer.IsAquifer, - WaterpressureInterpolationModel = layer.WaterpressureInterpolationModel, - SoilProfile = soilProfile1D - }; - soilProfile1D.Layers.Add(newLayer); - } - } - - // Add the newly created soil profile to the cache - cachedSoilProfiles1D[x] = soilProfile1D; - - return soilProfile1D; - } - - /// /// Returns a that represents this instance (either Name or SoilProfile1D). /// /// @@ -258,9 +114,9 @@ { return Name; } - if (orgSoilProfile != null) + if (originalSoilProfile1D != null) { - return orgSoilProfile.ToString(); + return originalSoilProfile1D.ToString(); } return "Soil Surface Profile with null Soil Profile"; } @@ -271,6 +127,10 @@ /// public SoilProfile2D ConvertToSoilProfile2D() { + // Generate the geometry for the new 2D profile + Create2DGeometryBasedOn1DProfileAndSurfaceLine(); + + // Create the actual new 2D profile var soilProfile2D = new SoilProfile2D { Geometry = Geometry @@ -280,15 +140,15 @@ return soilProfile2D; } - + /// - /// Updates the dike material. + /// Make sure that the top layer of the 1D profile is above the surface line by adding extra layer if needed using the dike material when needed. /// - private void UpdateDikeMaterial() + private void EnsureAvailabilityOfMaterialForNew2DGeometry() { - if (soilProfile != null && SurfaceLine != null) + if (soilProfile != null && surfaceLine2.Geometry != null) { - double topSurfaceLine = SurfaceLine.GetMaxZ(); + double topSurfaceLine = surfaceLine2.Geometry.GetMaxZ(); if (!double.IsNaN(topSurfaceLine) && (soilProfile.Layers.Count == 0 || soilProfile.TopLevel < topSurfaceLine)) { double newTopLevel = topSurfaceLine + 0.5; @@ -310,7 +170,7 @@ // special handling for dummy soil profiles if (soilProfile.Layers.Count == 0 && DikeEmbankmentMaterial == null) { - soilProfile.BottomLevel = SurfaceLine.GetMinZ() - 1; + soilProfile.BottomLevel = surfaceLine2.Geometry.GetMinZ() - 1; } soilProfile.Layers.Insert(0, newLayer); @@ -327,7 +187,7 @@ /// The layer. /// The minimum x. /// The maximum x. - private static void AddLayerGeometry(GeometryData data, SoilLayer1D layer, double minX, double maxX) + private static void Add2DLayerToGeometry(GeometryData data, SoilLayer1D layer, double minX, double maxX) { var top = layer.TopLevel; var bottom = layer.BottomLevel; @@ -353,140 +213,158 @@ data.Curves.Add(new GeometryCurve(topLeft, bottomLeft)); } - private void AddSurfaceLineGeometry(GeometryData data) + private GeometryPointString GetClonedSurfaceLineWithIntersections() { - var current = new Point2D(surfaceLine.Points[0].X, surfaceLine.Points[0].Z); - data.Points.Add(current); - for (int i = 1; i < surfaceLine.Points.Count; ++i) + var locsurfaceline = surfaceLine2.Geometry.Clone(); + var intersectionPoints = new List(); + var bounds = surfaceLine2.Geometry.GetGeometryBounds(); + double minX = bounds.Left; + double maxX = bounds.Right; + foreach (var locGeometryCurve in originalSoilProfile1D.Layers) { - var previous = current; - current = new Point2D(surfaceLine.Points[i].X, surfaceLine.Points[i].Z); - - data.Points.Add(current); - data.Curves.Add(new GeometryCurve(previous, current)); + var line = new Line + { + BeginPoint = new Point2D(minX - 5, locGeometryCurve.TopLevel), + EndPoint = new Point2D(maxX + 5, locGeometryCurve.TopLevel) + }; + foreach (var point2D in locsurfaceline.IntersectionPointsXzWithLineXz(line)) + { + intersectionPoints.Add(new GeometryPoint(point2D.X, point2D.Z)); + } } + locsurfaceline.Points.AddRange(intersectionPoints); + locsurfaceline.SortPointsByXAscending(); + locsurfaceline.RemoveNonAscendingPointsOnSurfaceLine(); + return locsurfaceline; } - private void AddSettlementZones(GeometryData locGeometry) + private GeometryPointString AddSurfaceLineToGeometry(GeometryData data) { - - if (surfaceLine2 != null && surfaceLine2.HasDike()) + var locSurfaceLine = GetClonedSurfaceLineWithIntersections(); + locSurfaceLine.SyncCalcPoints(); + var currentPoint = locSurfaceLine.CalcPoints[0]; + data.Points.Add(currentPoint); + for (int i = 1; i < locSurfaceLine.CalcPoints.Count; ++i) { - var xToeRiver = surfaceLine2.CharacteristicPoints.GetGeometryPoint(CharacteristicPointType.DikeToeAtRiver).X; - AddVerticalAt(locGeometry, xToeRiver); - var xTopRiver = surfaceLine2.CharacteristicPoints.GetGeometryPoint(CharacteristicPointType.DikeTopAtRiver).X; - var xRiver = xToeRiver + (xTopRiver - xToeRiver) / 3.0; - AddVerticalAt(locGeometry, xRiver); + var previousPoint = currentPoint; + currentPoint = locSurfaceLine.CalcPoints[i]; - var xToePolder = surfaceLine2.CharacteristicPoints.GetGeometryPoint(CharacteristicPointType.DikeToeAtPolder).X; - var xTopPolder = surfaceLine2.CharacteristicPoints.GetGeometryPoint(CharacteristicPointType.DikeTopAtPolder).X; - var xPolder = xTopPolder + 2 * (xToePolder - xTopPolder) / 3.0; - AddVerticalAt(locGeometry, xPolder); - AddVerticalAt(locGeometry, xToePolder); + data.Points.Add(currentPoint); + data.Curves.Add(new GeometryCurve(previousPoint, currentPoint)); } - } - private void AddVerticalAt(GeometryData aGeometry, double x) - { - const double zMax = 100000.0; - var topX = new Point2D(x, zMax); - var bottomX = new Point2D(x, -zMax); - aGeometry.Points.Add(topX); - aGeometry.Points.Add(bottomX); - aGeometry.Curves.Add(new GeometryCurve(topX, bottomX)); - //aGeometry.NewlyEffectedCurves.Add(aGeometry.Curves.Last());##Bka + return locSurfaceLine; } private void BuildGeometryModel() { - if (surfaceLine == null + if (surfaceLine2.Geometry == null || soilProfile == null || soilProfile.Layers.Count == 0 - || surfaceLine.GetMaxZ() < soilProfile.BottomLevel - || surfaceLine.Points.Count < 2) + || surfaceLine2.Geometry.GetMaxZ() < soilProfile.BottomLevel + || surfaceLine2.Geometry.Points.Count < 2) { return; } - var bounds = SurfaceLine.GetGeometryBounds(); + var bounds = surfaceLine2.Geometry.GetGeometryBounds(); double minX = bounds.Left; double maxX = bounds.Right; Geometry.Clear(); Geometry.Left = minX; Geometry.Right = maxX; - Geometry.Bottom = Math.Min(soilProfile.BottomLevel, surfaceLine.GetMinZ() - 1); + Geometry.Bottom = Math.Min(soilProfile.BottomLevel, surfaceLine2.Geometry.GetMinZ() - 1); - // add profile geometry + // add 1D profile layers as 2D layers (points and curves) to geometry soilProfile.Layers.Sort(); foreach (var layer in soilProfile.Layers) { - AddLayerGeometry(Geometry, layer, minX, maxX); + Add2DLayerToGeometry(Geometry, layer, minX, maxX); } - // add surface line geometry - AddSurfaceLineGeometry(Geometry); + // add surface line as points and curves to geometry, return the geometry of the surface line with the intersection points + var localSurfaceLineGeometry = AddSurfaceLineToGeometry(Geometry); - // add the extra lines to create the settlementzones if required - if (createSettlementZones) + // Remove all points and curve above the surface line + RemoveGeometryDataAboveSurfaceLine(localSurfaceLineGeometry); + + // Add curves from start and end of surface line to top layer. + Geometry.Curves.Add(new GeometryCurve(Geometry.Points[0], localSurfaceLineGeometry.CalcPoints[0])); + Geometry.Curves.Add(new GeometryCurve(Geometry.Points[1], localSurfaceLineGeometry.CalcPoints.Last())); + + // Now regenerate surfaces based on added points and curves + Geometry.NewlyEffectedPoints.AddRange(Geometry.Points); + Geometry.NewlyEffectedCurves.AddRange(Geometry.Curves); + Geometry.RegenerateGeometry(); + } + + private void RebuildSurfaces(GeometryData locGeometry) + { + Surfaces.Clear(); + if (locGeometry == null || surfaceLine2.Geometry == null) { - AddSettlementZones(Geometry); + return; } -// Geometry.NewlyEffectedPoints.AddRange(Geometry.Points); -// Geometry.NewlyEffectedCurves.AddRange(Geometry.Curves); -// Geometry.RegenerateGeometry(); -// Geometry.DeleteLooseCurves();##Bka + foreach (var geometrySurface in locGeometry.Surfaces) + { + var bounds = geometrySurface.GetGeometryBounds(); + var z = (bounds.Top + bounds.Bottom) * 0.5; + + var soilLayer = soilProfile.GetLayerAt(z); + + if (soilLayer != null) + { + Surfaces.Add(new SoilLayer2D + { + GeometrySurface = geometrySurface, + IsAquifer = soilLayer.IsAquifer, + WaterpressureInterpolationModel = soilLayer.WaterpressureInterpolationModel, + Soil = soilLayer.Soil + }); + } + } } - private void RebuildSurfaces(GeometryData locGeometry) + /// + /// Removes the geometry data above surface line. + /// + /// The surface line. + private void RemoveGeometryDataAboveSurfaceLine(GeometryPointString surfaceLine) { -// Surfaces.Clear(); -// if (locGeometry == null || surfaceLine == null) -// { -// return; -// } -// var gu = new GeotechnicsUtilities(); -// gu.RemoveGeometryDataAboveSurfaceLine(ref locGeometry, surfaceLine); -// foreach (var geometrySurface in locGeometry.Surfaces) -// { -// var bounds = geometrySurface.GetGeometryBounds(); -// var z = (bounds.Top + bounds.Bottom) * 0.5; -// -// var soilLayer = soilProfile.GetLayerAt(z); -// -// if (soilLayer != null) -// { -// Surfaces.Add(new SoilLayer2D -// { -// GeometrySurface = geometrySurface, -// IsAquifer = soilLayer.IsAquifer, -// WaterpressureInterpolationModel = soilLayer.WaterpressureInterpolationModel, -// Soil = soilLayer.Soil, -// SoilProfile1D = this -// }); -// } -// } ##Bka + // remove all points and attached curves that are above the surfaceline. + if (surfaceLine != null && Geometry.Curves.Count > 1) + { + var pointsToDelete = new List(); + foreach (var point in Geometry.Points.ToArray()) + { + // See if point is above surface line and if so remove it with its belonging curves + if (surfaceLine.GetZatX(point.X) < point.Z - GeometryConstants.Accuracy) + { + pointsToDelete.Add(point); + } + } + // finally remove points with their belonging curves + foreach (var pointToDelete in pointsToDelete) + { + Geometry.DeletePointWithCurves(pointToDelete); + } + } } /// /// Updates the layers. /// - private void UpdateLayers() + private void Create2DGeometryBasedOn1DProfileAndSurfaceLine() { -// dirty = false; -// initial = false; -// -// // Clear all cached soil profiles -// cachedSoilProfiles1D.Clear(); -// -// if (orgSoilProfile != null) -// { -// soilProfile = (SoilProfile1D)orgSoilProfile.Clone(); -// -// UpdateDikeMaterial(); -// BuildGeometryModel(); -// RebuildSurfaces(Geometry); -// } ##Bka + if (originalSoilProfile1D != null) + { + soilProfile = (SoilProfile1D)originalSoilProfile1D.Clone(); + + EnsureAvailabilityOfMaterialForNew2DGeometry(); + BuildGeometryModel(); + RebuildSurfaces(Geometry); + } } } } \ No newline at end of file Index: DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/Routines2D.cs =================================================================== diff -u -r1974 -r3079 --- DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/Routines2D.cs (.../Routines2D.cs) (revision 1974) +++ DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/Routines2D.cs (.../Routines2D.cs) (revision 3079) @@ -21,6 +21,8 @@ using System; using System.Collections.Generic; +using System.Linq; +using Deltares.DamEngine.Data.Standard; namespace Deltares.DamEngine.Data.Geometry { @@ -63,6 +65,17 @@ } /// + /// Enum for the clockwise options + /// + public enum Clockwise + { + IsClockwise = 0, + AntiClockwise = 1, + NotEnoughUniquePoints = 2, + PointsOnLine = 3 + } + + /// /// Static class Routines2D /// public static class Routines2D @@ -74,7 +87,7 @@ public double Z; } - private const double cEpsilon = 1.0e-4; + private const double CEpsilon = 1.0e-4; /// /// Checks if the 2D Lines strictly intersect. Strictly means that the intersection point must be part of both lines so @@ -99,12 +112,51 @@ double point3X, double point3Z, double point4X, double point4Z, out Point2D intersectionPoint) { return DetermineIf2DLinesIntersectStrickly(point1X, point1Z, point2X, point2Z, point3X, point3Z, point4X, point4Z, - out intersectionPoint, cEpsilon); + out intersectionPoint, CEpsilon); } + /// + /// Determines the point (aResultX, aResultY) on the line (aLine1X, aLine1Y, aLine2X, aLine2Y) with the + /// shortest possible distance to a point (aPointX, aPointY). + /// + /// + /// + /// + /// + /// + /// + /// + /// + public static void GetPointOnLineClosestTo(double aPointX, double aPointY, double aLine1X, double aLine1Y, double aLine2X, double aLine2Y, out double aResultX, out double aResultY) + { + + var lAb = Compute2DDistance(aLine1X, aLine1Y, aLine2X, aLine2Y); + var lAc = Compute2DDistance(aLine1X, aLine1Y, aPointX, aPointY); + var lBc = Compute2DDistance(aLine2X, aLine2Y, aPointX, aPointY); + var lAClose = (Math.Pow(lAc, 2) - Math.Pow(lBc, 2) + Math.Pow(lAb, 2)) / (2 * lAb); + if (lAClose <= 0) + { + aResultX = aLine1X; + aResultY = aLine1Y; + } + else + { + if (lAClose >= lAb) + { + aResultX = aLine2X; + aResultY = aLine2Y; + } + else + { + aResultX = aLine1X + ((lAClose / lAb) * (aLine2X - aLine1X)); + aResultY = aLine1Y + ((lAClose / lAb) * (aLine2Y - aLine1Y)); + } + } + } + /// - /// Doeses the point exist in line. + /// Does the point exist in line. /// /// The line point1 x. /// The line point1 z. @@ -164,6 +216,41 @@ } /// + /// To check if the points of the polygon are clock-wise + /// + /// + /// + public static Clockwise IsClockWise(IEnumerable aPolygon) + { + var distinctPoints = aPolygon.Distinct().ToArray(); + if (distinctPoints.Length < 3) + { + return Clockwise.NotEnoughUniquePoints; + } + + double sumClockwise = 0; + double clockwise; + for (int ii = 0; ii < distinctPoints.Length - 1; ii++) + { + clockwise = (distinctPoints[ii + 1].X - distinctPoints[ii].X) * + (distinctPoints[ii + 1].Z + distinctPoints[ii].Z); + sumClockwise += clockwise; + } + clockwise = (distinctPoints[0].X - distinctPoints[distinctPoints.Length - 1].X) * + (distinctPoints[0].Z + distinctPoints[distinctPoints.Length - 1].Z); + + sumClockwise += clockwise; + var result = Math.Sign(sumClockwise); + + if (result == 0) + { + return Clockwise.PointsOnLine; + } + + return result > 0.0 ? Clockwise.IsClockwise : Clockwise.AntiClockwise; + } + + /// /// Find if points Ax,ay are between a boundary polygon. /// /// The polygon. @@ -321,8 +408,38 @@ return result; } - /// + /// Determines angle between 2 lines. + /// Clockwise Input: point1, commonpoint, point2, commonpoint, normalvect(ex <0,0,-1>), validate + /// + /// + /// + /// + /// + /// + public static double FindAngle(Point2D line1Point1, Point2D line1Point2, Point2D line2Point1, + Point2D line2Point2) + { + var x1 = line1Point2.X - line1Point1.X; + var z1 = line1Point2.Z - line1Point1.Z; + + var x2 = line2Point2.X - line2Point1.X; + var z2 = line2Point2.Z - line2Point1.Z; + + var angle1 = GeneralMathRoutines.RadianToDegree(Math.Atan2(z1, x1)); + var angle2 = GeneralMathRoutines.RadianToDegree(Math.Atan2(z2, x2)); + + var angle = angle2 - angle1; + + if (angle < 0) + { + angle += 360; + } + + return angle; + } + + /// /// Checks if values are equal /// /// @@ -382,7 +499,16 @@ } - public static LineIntersection DetermineIf2DLinesIntersectWithExtrapolation(Point2D aPoint1, Point2D aPoint2, Point2D aPoint3, Point2D aPoint4, ref Point2D aIntersectionPoint) + /// + /// Determines the if 2D lines intersect using extrapolation when needed. + /// + /// a point1. + /// a point2. + /// a point3. + /// a point4. + /// a intersection point. + /// + public static LineIntersection DetermineIf2DLinesIntersectWithExtrapolation(Point2D aPoint1, Point2D aPoint2, Point2D aPoint3, Point2D aPoint4, out Point2D aIntersectionPoint) { var lLineConstant1 = CalculateNormalLineConstants(aPoint1.X, aPoint1.Z, aPoint2.X, aPoint2.Z); var lLineConstant2 = CalculateNormalLineConstants(aPoint3.X, aPoint3.Z, aPoint4.X, aPoint4.Z); @@ -412,7 +538,7 @@ var lB2 = aLine2Constant.Y; var lC2 = aLine2Constant.Z; - if (AreLinesParallel(lA1, lA2, lB1, lB2, cEpsilon)) + if (AreLinesParallel(lA1, lA2, lB1, lB2, CEpsilon)) { aIntersectionPoint = null; return LineIntersection.Parallel; @@ -594,7 +720,7 @@ { double q = Math.Sqrt(pointX * pointX + pointY * pointY); - if (q != 0) + if (!q.AlmostEquals(0, 1e-8)) { pointX = pointX / q; pointY = pointY / q; @@ -610,7 +736,7 @@ /// The point2 z. /// The tolerance. /// - private static bool DetermineIfPointsCoincide(double point1X, double point1Z, double point2X, double point2Z, double tolerance) + public static bool DetermineIfPointsCoincide(double point1X, double point1Z, double point2X, double point2Z, double tolerance) { if ((Math.Abs(point1X - point2X)) < tolerance && Math.Abs(point1Z - point2Z) < tolerance) { Index: DamEngine/trunk/src/Deltares.DamEngine.Data.Tests/General/GeometryDataTests.cs =================================================================== diff -u --- DamEngine/trunk/src/Deltares.DamEngine.Data.Tests/General/GeometryDataTests.cs (revision 0) +++ DamEngine/trunk/src/Deltares.DamEngine.Data.Tests/General/GeometryDataTests.cs (revision 3079) @@ -0,0 +1,281 @@ +// Copyright (C) Stichting Deltares 2019. All rights reserved. +// +// This file is part of the Dam Engine. +// +// The Dam Engine is free software: you can redistribute it and/or modify +// it under the terms of the GNU Affero General Public License as published by +// the Free Software Foundation, either version 3 of the License, or +// (at your option) any later version. +// +// This program is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +// GNU Affero General Public License for more details. +// +// You should have received a copy of the GNU Affero General Public License +// along with this program. If not, see . +// +// All names, logos, and references to "Deltares" are registered trademarks of +// Stichting Deltares and remain full property of Stichting Deltares at all times. +// All rights reserved. + +using System; +using System.Collections.Generic; +using System.Linq; +using System.Text; +using System.Threading.Tasks; +using Deltares.DamEngine.Data.Geometry; +using NUnit.Framework; + +namespace Deltares.DamEngine.Data.Tests.General +{ + [TestFixture] + public class GeometryDataTests + { + + [Test] + public void TestRemoveDoublesFromNewlyEffectedPointsAndCurves() + { + var point1 = new Point2D(1, 0); // to be kept + var point2 = new Point2D(1, 0); // to be removed as == p1 + var point3 = new Point2D(5, 0); // to be kept + var point4 = new Point2D(5, 0); // to be removed as == p3 + var point5 = new Point2D(1, 0); // to be removed as == p1 + var point6 = new Point2D(5, 0); // to be removed as == p3 + + var curve1 = new GeometryCurve(point1, point2); // to be removed as no lenght + var curve2 = new GeometryCurve(point2, point3); // to be kept + var curve3 = new GeometryCurve(point1, point2); // to be removed as no lenght + var curve4 = new GeometryCurve(point4, point1); // to be removed as == c2 (in reverse order) + var curve5 = new GeometryCurve(point3, point5); // to be removed as == c2 (in reverse order) + var curve6 = new GeometryCurve(point1, point1); // to be removed as no lenght + var curve7 = new GeometryCurve(point1, point1); // to be removed as no lenght + var curve8 = new GeometryCurve(null, null); // to be removed as no points + var curve9 = new GeometryCurve(null, point2); // to be removed as no headpoint + var curve10 = new GeometryCurve(point3, null); // to be removed as no endpoint + + var geometrymodel = new GeometryData(); + geometrymodel.NewlyEffectedPoints.Add(point1); + geometrymodel.NewlyEffectedPoints.Add(point2); + geometrymodel.NewlyEffectedPoints.Add(point3); + geometrymodel.NewlyEffectedPoints.Add(point4); + geometrymodel.NewlyEffectedPoints.Add(point5); + geometrymodel.NewlyEffectedPoints.Add(point6); + + geometrymodel.NewlyEffectedCurves.Add(curve1); + geometrymodel.NewlyEffectedCurves.Add(curve2); + geometrymodel.NewlyEffectedCurves.Add(curve3); + geometrymodel.NewlyEffectedCurves.Add(curve4); + geometrymodel.NewlyEffectedCurves.Add(curve5); + geometrymodel.NewlyEffectedCurves.Add(curve6); + geometrymodel.NewlyEffectedCurves.Add(curve7); + geometrymodel.NewlyEffectedCurves.Add(curve8); + geometrymodel.NewlyEffectedCurves.Add(curve9); + geometrymodel.NewlyEffectedCurves.Add(curve10); + + geometrymodel.RemoveDoublesFromNewlyEffectedPointsAndCurves(); + + // only two points are unique by location + Assert.AreEqual(2, geometrymodel.NewlyEffectedPoints.Count); + // only two curves hold both head and endpoint with unique locations + Assert.AreEqual(1, geometrymodel.NewlyEffectedCurves.Count); + } + + [Test] + public void GeometryHelperTestSurfaceTwoVertLayers() + { + var gData = CreateGeometrySurface2(); + var line = gData.SurfaceLine; + Assert.AreEqual(3, line.Points.Count); + Assert.AreEqual(1.0, line.Points.OrderBy(p => p.X).First().X); + Assert.AreEqual(10.0, line.Points.OrderByDescending(p => p.X).First().X); + } + + [Test] + public void GeometryDataTestSurfaceTwoHorLayers() + { + var gData = CreateGeometrySurface(); + var line = gData.SurfaceLine; + Assert.AreEqual(2, line.CalcPoints.Count); + Assert.AreEqual(1.0, line.CalcPoints.OrderBy(p => p.X).First().X); + Assert.AreEqual(10.0, line.CalcPoints.OrderByDescending(p => p.X).First().X); + } + + [Test] + public void TestSimpleGeneration() + { + var geom = new GeometryData + { + Left = 0, + Bottom = 0, + Right = 100 + }; + var p1 = new Point2D(0, 5); + geom.Points.Add(p1); + var p2 = new Point2D(100, 5); + geom.Points.Add(p2); + var p3 = new Point2D(0, 0); + geom.Points.Add(p3); + var p4 = new Point2D(100, 0); + geom.Points.Add(p4); + geom.Curves.Add(new GeometryCurve(p1, p2)); + geom.Curves.Add(new GeometryCurve(p3, p4)); + geom.Curves.Add(new GeometryCurve(p1, p3)); + geom.Curves.Add(new GeometryCurve(p2, p4)); + geom.NewlyEffectedCurves.Add(geom.Curves[0]); + geom.RegenerateGeometry(); + + Assert.AreEqual(4, geom.Points.Count); + Assert.AreEqual(4, geom.Curves.Count); + Assert.AreEqual(1, geom.Surfaces.Count); + Assert.AreEqual(0, geom.Surfaces[0].InnerLoops.Count); + Assert.AreEqual(4, geom.Surfaces[0].OuterLoop.Count); + Assert.AreEqual(4, geom.Surfaces[0].OuterLoop.CurveList.Count); + Assert.AreEqual(4, geom.Surfaces[0].OuterLoop.Points.Count); + Assert.IsTrue(geom.Surfaces[0].OuterLoop.HasArea()); + Assert.IsTrue(geom.Surfaces[0].OuterLoop.IsLoop()); + Assert.IsTrue(geom.Surfaces[0].OuterLoop.IsPointInLoopArea(new Point2D(25, 3))); + Assert.IsFalse(geom.Surfaces[0].OuterLoop.IsPointInLoopArea(new Point2D(25, 5.1))); + } + + private GeometryData CreateGeometrySurface() + { + var geometryModel = new GeometryData(); + + /* The following model looks as follows + * + * |----------| + * | | + * |----------| + * | | + * |----------| + * + */ + + var point1 = new Point2D(1, 0); + var point2 = new Point2D(10, 0); + var point3 = new Point2D(10, 10); + var point4 = new Point2D(1, 10); + var point5 = new Point2D(1, 5); + var point6 = new Point2D(10, 5); + + var curve1 = new GeometryCurve(point1, point2); + var curve2 = new GeometryCurve(point2, point3); + var curve3 = new GeometryCurve(point3, point4); + var curve4 = new GeometryCurve(point4, point1); + var curve5 = new GeometryCurve(point5, point6); + + geometryModel.Points.AddRange(new[] + { + point1, point2, point3, point4, point5, point6 + }); + geometryModel.Curves.AddRange(new[] + { + curve1, curve2, curve3, curve4, curve5 + }); + + geometryModel.NewlyEffectedPoints.AddRange(geometryModel.Points); + geometryModel.NewlyEffectedCurves.AddRange(geometryModel.Curves); + geometryModel.RegenerateGeometry(); + return geometryModel; + } + + private GeometryData CreateGeometrySurface2() + { + var geometryModel = new GeometryData(); + + /* The following model looks as follows + * + * |-----|----| + * | | | + * | | | + * | | | + * |-----|----| + * + */ + + var point1 = new Point2D(1, 0); + var point2 = new Point2D(1, 10); + var point3 = new Point2D(5, 10); + var point4 = new Point2D(5, 0); + var point5 = new Point2D(10, 10); + var point6 = new Point2D(10, 0); + + var curve1 = new GeometryCurve(point1, point2); + var curve2 = new GeometryCurve(point2, point3); + var curve3 = new GeometryCurve(point3, point4); + var curve4 = new GeometryCurve(point4, point1); + var curve5 = new GeometryCurve(point3, point5); + var curve6 = new GeometryCurve(point5, point6); + var curve7 = new GeometryCurve(point6, point4); + + geometryModel.Points.AddRange(new[] + { + point1, point2, point3, point4, point5, point6 + }); + geometryModel.Curves.AddRange(new[] + { + curve1, curve2, curve3, curve4, curve5, curve5, curve6, curve7 + }); + + geometryModel.NewlyEffectedPoints.AddRange(geometryModel.Points); + geometryModel.NewlyEffectedCurves.AddRange(geometryModel.Curves); + geometryModel.RegenerateGeometry(); + return geometryModel; + } + + private GeometryData CreateDonutGeometry() + { + var geometryModel = new GeometryData + { + Left = -10, + Bottom = -10, + Right = 20 + }; + + /* The following model looks as follows + * + * |----------| + * | ------ | + * | | | | + * | | | | + * | ------ | + * |----------| + * + */ + + var point1 = new Point2D(0, 0); + var point2 = new Point2D(0, 10); + var point3 = new Point2D(10, 10); + var point4 = new Point2D(10, 0); + var point5 = new Point2D(3, 3); + var point6 = new Point2D(3, 7); + var point7 = new Point2D(7, 7); + var point8 = new Point2D(7, 3); + + var curve1 = new GeometryCurve(point1, point2); + var curve2 = new GeometryCurve(point2, point3); + var curve3 = new GeometryCurve(point3, point4); + var curve4 = new GeometryCurve(point4, point1); + var curve5 = new GeometryCurve(point5, point6); + var curve6 = new GeometryCurve(point6, point7); + var curve7 = new GeometryCurve(point7, point8); + var curve8 = new GeometryCurve(point8, point5); + + geometryModel.Points.AddRange(new[] + { + point1, point2, point3, point4, point5, point6, point7, point8 + }); + geometryModel.Curves.AddRange(new[] + { + curve1, curve2, curve3, curve4, curve5, curve5, curve6, curve7, curve8 + }); + + geometryModel.NewlyEffectedPoints.AddRange(geometryModel.Points); + geometryModel.NewlyEffectedCurves.AddRange(geometryModel.Curves); + geometryModel.RegenerateGeometry(); + return geometryModel; + } + + } +} Index: DamEngine/trunk/src/Deltares.DamEngine.Data/Geotechnics/SoilProfile2D.cs =================================================================== diff -u -r1974 -r3079 --- DamEngine/trunk/src/Deltares.DamEngine.Data/Geotechnics/SoilProfile2D.cs (.../SoilProfile2D.cs) (revision 1974) +++ DamEngine/trunk/src/Deltares.DamEngine.Data/Geotechnics/SoilProfile2D.cs (.../SoilProfile2D.cs) (revision 3079) @@ -19,6 +19,7 @@ // Stichting Deltares and remain full property of Stichting Deltares at all times. // All rights reserved. +using System; using System.Collections.Generic; using Deltares.DamEngine.Data.Geometry; using Deltares.DamEngine.Data.Standard.Language; @@ -31,7 +32,6 @@ /// public class SoilProfile2D : SoilProfile { - protected readonly Dictionary cachedSoilProfiles1D = new Dictionary(); protected GeometryData geometry = new GeometryData(); protected readonly List surfaces = new List(); @@ -84,15 +84,9 @@ public virtual SoilProfile1D GetSoilProfile1D(double x) { const double diff = 0.001; - SoilProfile1D soilProfile = GetCachedSoilProfile1D(x); - - if (soilProfile != null) + + var soilProfile = new SoilProfile1D { - return soilProfile; - } - - soilProfile = new SoilProfile1D - { Name = "Generated at x = " + x + " from " + Name }; var detector = new LayerDetector(Surfaces); @@ -120,9 +114,6 @@ } } - - cachedSoilProfiles1D[x] = soilProfile; - return soilProfile; } @@ -170,19 +161,5 @@ { return Name; } - - /// - /// Gets the cached soil profile1 d. - /// - /// The x. - /// - protected SoilProfile1D GetCachedSoilProfile1D(double x) - { - if (cachedSoilProfiles1D.ContainsKey(x)) - { - return cachedSoilProfiles1D[x]; - } - return null; - } } } \ No newline at end of file Index: DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometryLoop.cs =================================================================== diff -u -r1974 -r3079 --- DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometryLoop.cs (.../GeometryLoop.cs) (revision 1974) +++ DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometryLoop.cs (.../GeometryLoop.cs) (revision 3079) @@ -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.Linq; namespace Deltares.DamEngine.Data.Geometry { @@ -88,6 +90,52 @@ } /// + /// Determines whether this instance has area. + /// + /// + public bool HasArea() + { + if (CurveList.Count < 3) + { + return false; + } + + var points = new List(); + points.AddRange(Points); + points.Add(points[0]); + + // Calculate area of polygon using Shoelace algorithm: + var sum = 0.0; + for (int i = 1; i < points.Count; i++) + { + sum += points[i - 1].X * points[i].Z - points[i - 1].Z * points[i].X; + } + + return Math.Abs(sum) > 1e-6; + } + + /// + /// Determines whether [is clock wise]. + /// + /// + /// + /// Cannot determine if loop is clockwise if checked location forms a straight line with its neighboring vectors. + public bool IsClockWise() + { + var polyGon = GetLocalPoint2DList(); + var isClockWise = Routines2D.IsClockWise(polyGon); + if (isClockWise == Clockwise.NotEnoughUniquePoints) + { + throw new NotEnoughUniquePointsException(); + } + if (isClockWise == Clockwise.PointsOnLine) + { + throw new InvalidOperationException("Cannot determine if loop is clockwise if checked location forms a straight line with its neighboring vectors."); + } + return isClockWise == Clockwise.IsClockwise; + } + + /// /// See if a point lies in a closed surface /// /// @@ -98,6 +146,15 @@ } /// + /// Gets the local 2D Points List. + /// + /// + private List GetLocalPoint2DList() + { + return CalcPoints;//Points.Select(loopPoint => new Point2D(loopPoint.X, loopPoint.Z)).ToList(); + } + + /// /// Populates the loop's GeometryPoint list. /// private void PopulateLoopPointList() @@ -150,6 +207,18 @@ calcPoints.RemoveAt(calcPoints.Count - 1); } } - } + } + + /// + /// Helper class NotEnoughUniquePointsException + /// + public class NotEnoughUniquePointsException : InvalidOperationException + { + /// + /// Initializes a new instance of the class. + /// + public NotEnoughUniquePointsException() + : base("At least 3 unique points are required to determine if the loop is running clockwise.") { } + } } } \ No newline at end of file Index: DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometryCurve.cs =================================================================== diff -u -r1974 -r3079 --- DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometryCurve.cs (.../GeometryCurve.cs) (revision 1974) +++ DamEngine/trunk/src/Deltares.DamEngine.Data/Geometry/GeometryCurve.cs (.../GeometryCurve.cs) (revision 3079) @@ -125,6 +125,17 @@ } /// + /// Locations the equals. + /// + /// a curve. + /// + public bool LocationEquals(GeometryCurve aCurve) + { + return (headPoint.LocationEquals(aCurve.HeadPoint) && endPoint.LocationEquals(aCurve.EndPoint)) || + (headPoint.LocationEquals(aCurve.EndPoint) && endPoint.LocationEquals(aCurve.HeadPoint)); + } + + /// /// Swaps and . /// public void Reverse() @@ -186,7 +197,7 @@ /// The point. /// The tolerance. /// - public bool ContainsPoint(Point2D point, double tolerance) + public bool ContainsPoint(Point2D point, double tolerance = GeometryConstants.Accuracy * 0.5) { return HeadPoint != null && EndPoint != null && Routines2D.DoesPointExistInLine(HeadPoint.X, HeadPoint.Z, EndPoint.X, EndPoint.Z, point.X, point.Z, tolerance); @@ -202,6 +213,18 @@ public override string ToString() { return string.Empty; - } + } + + /// + /// Gets or sets the surface at the Left. + /// + public virtual GeometrySurface SurfaceAtLeft { get; set; } + + + /// + /// Gets or sets the surface at the Right. + /// + public virtual GeometrySurface SurfaceAtRight { get; set; } + } } \ No newline at end of file Index: DamEngine/trunk/src/Deltares.DamEngine.Data.Tests/Geotechnics/SoilSurfaceProfileTests.cs =================================================================== diff -u --- DamEngine/trunk/src/Deltares.DamEngine.Data.Tests/Geotechnics/SoilSurfaceProfileTests.cs (revision 0) +++ DamEngine/trunk/src/Deltares.DamEngine.Data.Tests/Geotechnics/SoilSurfaceProfileTests.cs (revision 3079) @@ -0,0 +1,86 @@ +// Copyright (C) Stichting Deltares 2019. All rights reserved. +// +// This file is part of the Dam Engine. +// +// The Dam Engine is free software: you can redistribute it and/or modify +// it under the terms of the GNU Affero General Public License as published by +// the Free Software Foundation, either version 3 of the License, or +// (at your option) any later version. +// +// This program is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +// GNU Affero General Public License for more details. +// +// You should have received a copy of the GNU Affero General Public License +// along with this program. If not, see . +// +// All names, logos, and references to "Deltares" are registered trademarks of +// Stichting Deltares and remain full property of Stichting Deltares at all times. +// All rights reserved.using System; + +using System.Collections.Generic; +using System.Linq; +using System.Text; +using System.Threading.Tasks; +using Deltares.DamEngine.Data.Geometry; +using Deltares.DamEngine.Data.Geotechnics; +using Deltares.DamEngine.TestHelpers.Factories; +using NUnit.Framework; + +namespace Deltares.DamEngine.Data.Tests.Geotechnics +{ + [TestFixture] + public class SoilSurfaceProfileTests + { + + [Test] + public void TestCreate2DProfileWithSimpleSurfaceLineAboveProfile() + { + // Create a 1D profiles with toplayer at -2 + var prof1D = CreateTwoLayerProfile(-2); + // surfaceLine has 7 points, varies in height from -1 to 5 m + var surfaceLine = FactoryForSurfaceLines.CreateSimpleSurfaceLineForExitPointTest(); + + // Create profile + var soilSurfaceProfile = new SoilSurfaceProfile + { + SoilProfile = prof1D, + SurfaceLine2 = surfaceLine, + DikeEmbankmentMaterial = new Soil + { + Name = "Dike material" + } + }; + var profile2D = soilSurfaceProfile.ConvertToSoilProfile2D(); + Assert.AreEqual(13, profile2D.Geometry.Points.Count); + Assert.AreEqual(15, profile2D.Geometry.Curves.Count); + Assert.AreEqual(3, profile2D.Surfaces.Count); + } + + private static SoilProfile1D CreateTwoLayerProfile(double offset = 0) + { + var soilProfile = new SoilProfile1D(); + + var layer = new SoilLayer1D(); + layer.TopLevel = 0.0 + offset; + layer.Soil = new Soil("HW-OBO", 12.0, 10.0); + layer.Soil.DryUnitWeight = 0.01; + layer.IsAquifer = false; + layer.WaterpressureInterpolationModel = WaterpressureInterpolationModel.Automatic; + soilProfile.Layers.Add(layer); + + layer = new SoilLayer1D(); + layer.TopLevel = -5.0 + offset; + layer.Soil = new Soil("Alg-zand (0-30)", 22.0, 20.0); + layer.Soil.DryUnitWeight = 0.01; + layer.IsAquifer = true; + layer.WaterpressureInterpolationModel = WaterpressureInterpolationModel.Hydrostatic; + soilProfile.Layers.Add(layer); + + soilProfile.BottomLevel = -10.0 + offset; + + return soilProfile; + } + } +} Index: DamEngine/trunk/src/Deltares.DamEngine.Data/General/Location.cs =================================================================== diff -u -r2908 -r3079 --- DamEngine/trunk/src/Deltares.DamEngine.Data/General/Location.cs (.../Location.cs) (revision 2908) +++ DamEngine/trunk/src/Deltares.DamEngine.Data/General/Location.cs (.../Location.cs) (revision 3079) @@ -1139,11 +1139,11 @@ /// public Soil GetDikeEmbankmentSoil() { - if ((stabilityOptions != null && String.IsNullOrEmpty(stabilityOptions.SoilDatabaseName)) || String.IsNullOrEmpty(DikeEmbankmentMaterial)) + if (SoilList != null && SoilList.Soils.Count > 0 && !String.IsNullOrEmpty(DikeEmbankmentMaterial)) { - return null; + return SoilList.GetSoilByName(DikeEmbankmentMaterial); } - return SoilList.GetSoilByName(DikeEmbankmentMaterial); + return null; } /// Index: DamEngine/trunk/src/Deltares.DamEngine.Data/Geotechnics/CharacteristicPointSet.cs =================================================================== diff -u -r1974 -r3079 --- DamEngine/trunk/src/Deltares.DamEngine.Data/Geotechnics/CharacteristicPointSet.cs (.../CharacteristicPointSet.cs) (revision 1974) +++ DamEngine/trunk/src/Deltares.DamEngine.Data/Geotechnics/CharacteristicPointSet.cs (.../CharacteristicPointSet.cs) (revision 3079) @@ -427,9 +427,12 @@ // Only add new point if no point on same location was found PerformListInsertWithEvents(Geometry.Points, item.GeometryPoint, geometryIndex); } - } + } - UpdateCharacteristicPointHeight(item); + if (geometry != null) + { + UpdateCharacteristicPointHeight(item); + } UpdateTypeCache(item); }