// Copyright (C) Stichting Deltares 2024. 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 Deltares.DamEngine.Data.Geometry; using Deltares.DamEngine.Data.Standard.Language; namespace Deltares.DamEngine.Data.Geotechnics; /// /// Helper class for SoilProfile2D and SurfaceLine2 /// public static class SoilProfile2DSurfaceLineHelper { /// /// Determines if the surface line is inside the soil profile. /// /// The surface-line. /// The soil profile 2D /// public static bool IsSurfaceLineInsideSoilProfile2D(SurfaceLine2 surfaceLine, SoilProfile2D soilProfile2D) { const int precisionDecimals = 3; if (soilProfile2D == null) { return false; } if (Math.Round(surfaceLine.CharacteristicPoints.GetGeometryPoint(CharacteristicPointType.SurfaceLevelOutside).X, precisionDecimals) < Math.Round(soilProfile2D.Geometry.GetGeometryBounds().Left, precisionDecimals)) { return false; } if (Math.Round(surfaceLine.CharacteristicPoints.GetGeometryPoint(CharacteristicPointType.SurfaceLevelInside).X, precisionDecimals) > Math.Round(soilProfile2D.Geometry.GetGeometryBounds().Right, precisionDecimals)) { return false; } return true; } /// /// Combines the surface line with SoilProfile2D to create a new profile2D. /// /// The surface line. /// The SoilProfile2D. /// The default soil. /// The required shift of the geometry so it will fit the surface line /// public static SoilProfile2D CombineSurfaceLineWithSoilProfile2D(GeometryPointString surfaceLine, SoilProfile2D soilProfile2D, Soil defaultSoil, double shift) { 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; } /// /// Check if all the points of the surface line are above the bottom of the soil profile. /// /// /// /// public static bool IsSurfaceLineAboveBottomSoilProfile2D(SurfaceLine2 surfaceLine, SoilProfile2D soilProfile2D) { double bottom = soilProfile2D.Geometry.Bottom; foreach (GeometryPoint point in surfaceLine.Geometry.Points) { if (point.Z < bottom) { return false; } } 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.Curves.Clear(); result.Geometry.Points.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) { var soilLayer2D1 = new SoilLayer2D { GeometrySurface = surface, 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; } if (layerFromOldSurfaces == null) { SoilLayer2D oldLayer = GetSoilIfSurfaceIsLeftOrRightOfOldSurfaces(surface, oldSurfaces); if (oldLayer != null) { soilLayer2D2.IsAquifer = oldLayer.IsAquifer; soilLayer2D2.WaterpressureInterpolationModel = oldLayer.WaterpressureInterpolationModel; soilLayer2D2.SoilName = oldLayer.SoilName; soilLayer2D2.Soil = oldLayer.Soil; } else { soilLayer2D2.IsAquifer = false; soilLayer2D2.WaterpressureInterpolationModel = WaterpressureInterpolationModel.Automatic; soilLayer2D2.SoilName = defaultSoil.Name; soilLayer2D2.Soil = defaultSoil; } } 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); } } private static SoilLayer2D GetSoilIfSurfaceIsLeftOrRightOfOldSurfaces(GeometrySurface surface, List oldSurfaces) { double[] xMin = surface.OuterLoop.CalcPoints.Select(p => p.X).OrderBy(x => x).Distinct().ToArray(); double[] zMin = surface.OuterLoop.CalcPoints.Select(p => p.Z).OrderBy(z => z).Distinct().ToArray(); var leftPoint = new Point2D { X = xMin[0] - 0.1, Z = zMin[0] + 0.1 }; SoilLayer2D soilLayer = FindLayerFromPointInOldSurfaces(oldSurfaces, leftPoint); if (soilLayer != null) { return soilLayer; } var rightPoint = new Point2D { X = xMin[^1] + 0.1, Z = zMin[0] + 0.1 }; soilLayer = FindLayerFromPointInOldSurfaces(oldSurfaces, rightPoint); return soilLayer; } private static SoilLayer2D FindLayerFromPointInOldSurfaces(List oldSurfaces, Point2D searchPoint) { for (var i = 0; i < oldSurfaces.Count; i++) { bool find = oldSurfaces[i].GeometrySurface.OuterLoop.IsPointInLoopArea(searchPoint); if (find) { return oldSurfaces[i]; } } return null; } }