From e1880c65494dadbe2422851e71e6f5df229e7b9c Mon Sep 17 00:00:00 2001
From: Theophile Terraz <theophile.terraz@inrae.fr>
Date: Thu, 28 Sep 2023 17:09:40 +0200
Subject: [PATCH] work on purge

---
 src/Model/Geometry/PointXYZ.py   | 37 +++++++++++++++++++++++++++++++-
 src/Model/Geometry/ProfileXYZ.py | 14 ++++++++++++
 src/Model/Geometry/Reach.py      |  2 +-
 3 files changed, 51 insertions(+), 2 deletions(-)

diff --git a/src/Model/Geometry/PointXYZ.py b/src/Model/Geometry/PointXYZ.py
index a8f19c8b..22d2f674 100644
--- a/src/Model/Geometry/PointXYZ.py
+++ b/src/Model/Geometry/PointXYZ.py
@@ -16,7 +16,7 @@
 
 # -*- coding: utf-8 -*-
 
-from math import dist
+from math import dist, sqrt, acos
 import numpy as np
 
 from Model.Tools.PamhyrDB import SQLSubModel
@@ -212,3 +212,38 @@ class PointXYZ(Point, SQLSubModel):
             Euclidean distance between the two points
         """
         return dist((p1.x, p1.y, p1.z), (p2.x, p2.y, p2.z))
+
+    @staticmethod
+    def angle(p1, p2, p3):
+        """angle p1 p2 p3.
+
+        Args:
+            p1: A XYZ Point
+            p2: A XYZ Point
+            p3: A XYZ Point
+
+        Returns:
+            angle p1 p2 p3
+        """
+        a = PointXYZ.distance(p1, p2)
+        b = PointXYZ.distance(p2, p3)
+        c = PointXYZ.distance(p1, p3)
+        return acos((a**2 + b**2 - c**2) / (2 * a * b))
+
+    @staticmethod
+    def triangle_height(p1, p2, p3):
+        """triangle p1 p2 p3 height starting at p2.
+
+        Args:
+            p1: A XYZ Point
+            p2: A XYZ Point
+            p3: A XYZ Point
+
+        Returns:
+            triangle height starting at p2
+        """
+        a = PointXYZ.distance(p1, p2)
+        b = PointXYZ.distance(p2, p3)
+        c = PointXYZ.distance(p1, p3)
+        h2 = a**2 - ((c**2 - b**2 + a**2) / 2*c)**2
+        return sqrt(h2)
diff --git a/src/Model/Geometry/ProfileXYZ.py b/src/Model/Geometry/ProfileXYZ.py
index 23b577e8..c5066ead 100644
--- a/src/Model/Geometry/ProfileXYZ.py
+++ b/src/Model/Geometry/ProfileXYZ.py
@@ -450,6 +450,20 @@ class ProfileXYZ(Profile, SQLSubModel):
 
         return points
 
+    @timer
+    def purge2(self, np) -> list:
+
+        points = self.points
+        if (self.number_points <= np): return points
+        while len(points) > np:
+            l2 = []
+            for i in range(len(points)-2)
+                l2.append(PointXYZ.angle(points(i),points(i+1),points(i+2)))
+            imax = np.argmax(l2)+1
+            points.pop(imax)
+
+        return points
+
 
 
 
diff --git a/src/Model/Geometry/Reach.py b/src/Model/Geometry/Reach.py
index 0306ae90..dd7e8138 100644
--- a/src/Model/Geometry/Reach.py
+++ b/src/Model/Geometry/Reach.py
@@ -802,7 +802,7 @@ class Reach(SQLSubModel):
     def purge(self, np: int) -> list:
         profiles = []
         for profile in self.profiles():
-            points = profile.purge(np)
+            points = profile.purge2(np)
             profiles.append(ProfileXYZ())
             profiles[-1].import_points(points)