// 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;
}
}