From a22faae428f7ff51a0b2eaea902e5b5ed1a7c8a5 Mon Sep 17 00:00:00 2001
From: Arne Dußin
Date: Wed, 11 Nov 2020 21:36:44 +0100
Subject: Add polygon that can check if a point is inside

..except for one super super edgy edge case, but I wanted to get this
algorithm out into a commit before I ruin it completely (probably)
---
 src/math/line_segment.rs | 116 +++++++++++++++++++++++++++++++++++++++++++++++
 1 file changed, 116 insertions(+)
 create mode 100644 src/math/line_segment.rs

(limited to 'src/math/line_segment.rs')

diff --git a/src/math/line_segment.rs b/src/math/line_segment.rs
new file mode 100644
index 0000000..d7bcdb1
--- /dev/null
+++ b/src/math/line_segment.rs
@@ -0,0 +1,116 @@
+use super::Vec2;
+use alga::general::{ClosedMul, ClosedSub};
+use nalgebra::Scalar;
+use num_traits::Zero;
+
+pub struct LineSegment<T: Scalar + Copy> {
+    pub start: Vec2<T>,
+    pub end: Vec2<T>,
+}
+
+impl<T: Scalar + Copy> LineSegment<T> {
+    pub fn new(start: Vec2<T>, end: Vec2<T>) -> Self {
+        Self { start, end }
+    }
+
+    /* Helper function to check if this line contains a point. This function is very efficient, but
+     * must only be used for points that are collinear with the segment. This is checked only by
+     * assertion, so make sure that everything is ok in release mode.
+     */
+    pub(crate) fn contains_collinear(&self, point: Vec2<T>) -> bool
+    where
+        T: PartialOrd,
+    {
+        point.x <= super::partial_max(self.start.x, self.end.x)
+            && point.x >= super::partial_max(self.start.x, self.end.x)
+            && point.y <= super::partial_min(self.start.y, self.end.y)
+            && point.y >= super::partial_min(self.start.y, self.end.y)
+    }
+
+    /// Checks if two segments intersect. This is much more efficient than trying to find the exact
+    /// point of intersection, so it should be used if it is not strictly necessary to know it.
+    pub fn intersect(ls1: &LineSegment<T>, ls2: &LineSegment<T>) -> bool
+    where
+        T: Scalar + Copy + ClosedSub + ClosedMul + PartialOrd + Zero,
+    {
+        /* This algorithm works by checking the triplet orientation of the first lines points with the
+         * first point of the second line. After that it does the same for the second point of the
+         * second line. If the orientations are different, that must mean the second line starts
+         * "before" the first line and ends "after" it. It does the same, but with the roles of first
+         * and second line reversed. If both of these conditions are met, it follows that the lines
+         * must have crossed.
+         *
+         * Edge case: If an orientation comes out as collinear, the point of the other line that was
+         * checked may be on the other line or after/before it (if you extend the segment). If it is on
+         * the other line, this also means, the lines cross, since one line "stands" on the other.
+         * If however none of the collinear cases are like this, the lines cannot touch, because line
+         * segment a point was collinear with was not long enough.
+         */
+
+        // Cache the necessary orientations.
+        let ls1_ls2start_orientation = triplet_orientation(ls1.start, ls1.end, ls2.start);
+        let ls1_ls2end_orientation = triplet_orientation(ls1.start, ls1.end, ls2.end);
+        let ls2_ls1start_orientation = triplet_orientation(ls2.start, ls2.end, ls1.start);
+        let ls2_ls1end_orientation = triplet_orientation(ls2.start, ls2.end, ls1.end);
+
+        // Check for the first case that wase described (general case).
+        if ls1_ls2start_orientation != ls1_ls2end_orientation
+            && ls2_ls1start_orientation != ls2_ls1end_orientation
+        {
+            return true;
+        }
+
+        // Check if the start of ls2 lies on ls1
+        if ls1_ls2start_orientation == TripletOrientation::Collinear
+            && ls1.contains_collinear(ls2.start)
+        {
+            return true;
+        }
+        // Check if the end of ls2 lies on ls1
+        if ls1_ls2end_orientation == TripletOrientation::Collinear
+            && ls1.contains_collinear(ls2.end)
+        {
+            return true;
+        }
+
+        // Check if the start of ls1 lies on ls2
+        if ls2_ls1start_orientation == TripletOrientation::Collinear
+            && ls2.contains_collinear(ls1.start)
+        {
+            return true;
+        }
+        // Check if the end of ls1 lies on ls2
+        if ls2_ls1end_orientation == TripletOrientation::Collinear
+            && ls2.contains_collinear(ls1.end)
+        {
+            return true;
+        }
+
+        false
+    }
+}
+
+#[derive(PartialEq, Eq)]
+pub(crate) enum TripletOrientation {
+    Clockwise,
+    Counterclockwise,
+    Collinear,
+}
+
+/// Helper function to determine which direction one would turn to traverse from the first point
+/// over the second to the third point. The third option is collinear, in which case the three points
+/// are on the same line.
+pub(crate) fn triplet_orientation<T>(a: Vec2<T>, b: Vec2<T>, c: Vec2<T>) -> TripletOrientation
+where
+    T: Scalar + Copy + ClosedSub + ClosedMul + PartialOrd + Zero,
+{
+    /* Check the slopes of the vector from a to b and b to c. If the slope of ab is greater than
+     * that of bc, the rotation is clockwise. If ab is smaller than bc it's counterclockwise. If
+     * they are the same it follows that the three points are collinear.
+     */
+    match (b.y - a.y) * (c.x - b.x) - (b.x - a.x) * (c.y - b.y) {
+        q if q > T::zero() => TripletOrientation::Clockwise,
+        q if q < T::zero() => TripletOrientation::Counterclockwise,
+        _ => TripletOrientation::Collinear,
+    }
+}
-- 
cgit v1.2.3-70-g09d2


From 8ae528e9a2b219c5ac2f0b967c0cb91e79921e21 Mon Sep 17 00:00:00 2001
From: Arne Dußin
Date: Thu, 12 Nov 2020 23:37:43 +0100
Subject: Fix the polygon point containment algorithm

The algorithm before was really not working for a lot of edge cases and
very difficult to adapt. This version is definitely not the be-all and
end-all, but it should work for most (well, hopefully all) cases. After
refactoring and hopefully simplifying and straightening out the logic a
little more, it should be verifiable.
---
 src/math/line_segment.rs |   1 +
 src/math/polygon.rs      | 111 ++++++++++++++++++++++++++++++++++-------------
 2 files changed, 83 insertions(+), 29 deletions(-)

(limited to 'src/math/line_segment.rs')

diff --git a/src/math/line_segment.rs b/src/math/line_segment.rs
index d7bcdb1..338a681 100644
--- a/src/math/line_segment.rs
+++ b/src/math/line_segment.rs
@@ -3,6 +3,7 @@ use alga::general::{ClosedMul, ClosedSub};
 use nalgebra::Scalar;
 use num_traits::Zero;
 
+#[derive(Debug)]
 pub struct LineSegment<T: Scalar + Copy> {
     pub start: Vec2<T>,
     pub end: Vec2<T>,
diff --git a/src/math/polygon.rs b/src/math/polygon.rs
index e70cd2c..2af669a 100644
--- a/src/math/polygon.rs
+++ b/src/math/polygon.rs
@@ -23,54 +23,106 @@ impl<T: Scalar + Copy> Polygon<T> {
             return false;
         }
 
-        /* What's happening here:
-         *
-         * This algorithm creates a horizontal line for the point that should be checked. It goes
-         * to infinity (since infinite lines are not possible, it just goes to the maximum x value
-         * of all corner points). Then, the times the line crosses a polygon edge is counted. If
-         * the times it hit a polygon edge is uneven, it must have been inside, otherwise outside.
-         *
-         * There is an edge case however (pun not intended); When the line from the point is collinear
-         * to a polygon edge. it might hit this edge once and no other (it's at the top y or bottom y
-         * of the polygon), but may still not be inside the polygon. In this case, it's checked
-         * whether the point is in between the corner points of this polygon edge or not. If it is
-         * in between, it is still inside the polygon, otherwise outside.
-         */
-
         // Get the maximum x for the horizontal lines for "Infinity".
         let mut max_x = self.corners[0].x;
-        for corner in self.corners.iter().skip(1) {
+        for corner in self.corners.iter() {
+            /* Handle the edge case where the point is on a polygon corner, so we don't have to
+             * worry about it later.
+             */
+            if point == *corner {
+                return true;
+            }
+
             max_x = super::partial_max(corner.x, max_x);
         }
-        println!("Max x is {:?}", max_x);
+
         // The horizontally "infinite" point of the test point.
         let line_to_infinity = LineSegment::new(point, Vec2::new(max_x, point.y));
 
-        // Check for intersection for all the polygon edges.
-        let mut num_intersects = 0;
+        // The number of times the point line to infinity pierced through the polygon hull.
+        let mut num_pierced = 0;
+
+        /* Find edges that should be ignored. These are all edges, where the straight line passes
+         * right through a corner of the edge. It might still be counted as piercing the polygon,
+         * if the next corners lie below and the previous corner lies above the line or vice
+         * versa, but not if both points lie above it or below. The corners are exempted from the
+         * run through the edges. The edges that should be ignored are saved by their starting
+         * corner's index.
+         */
+        let mut ignore_edge = vec![false; self.corners.len()];
         for i in 0..self.corners.len() {
-            // Wrap around for the edge that returns to the start.
-            let j = (i + 1) % self.corners.len();
-            let current_edge = LineSegment::new(self.corners[i], self.corners[j]);
+            if point.y == self.corners[i].y {
+                // The previous corner index
+                let prev_i = (i + self.corners.len() - 1) % self.corners.len();
+
+                // Make sure the edges containing this corner will be ignored.
+                ignore_edge[i] = true;
+                ignore_edge[prev_i] = true;
+
+                /* Make sure we don't count a pierce double because of consecutive collinear
+                 * points. (Although the polygon really should be sanitised, but you know)
+                 */
+                if self.corners[prev_i].y == point.y {
+                    continue;
+                }
+
+                /* Ignore coming collinear points, since we are only interested in when the
+                 * y-coordinate finally changes
+                 */
+                let mut j = (i + 1) % self.corners.len();
+                while i != j {
+                    if self.corners[j].y != point.y {
+                        if (self.corners[prev_i].y < point.y && self.corners[j].y > point.y)
+                            || (self.corners[prev_i].y > point.y && self.corners[j].y < point.y)
+                        {
+                            /* Only count as pierced when it would be right of the point, otherwise
+                             * it would be like checking to the left of the line that is drawn to
+                             * infinity, which is incorrect
+                             */
+                            if self.corners[i].x >= point.x {
+                                //println!("Pierced through corner {:?}", self.corners[i]);
+                                num_pierced += 1;
+                            }
+                        }
+
+                        break;
+                    }
+
+                    j = (j + 1) % self.corners.len();
+                }
+            }
+        }
+
+        // Find the points where the x-line may actually collide with a proper line.. ^^
+        for (i, edge_start) in self.corners.iter().enumerate() {
+            if ignore_edge[i] {
+                continue;
+            }
+
+            let next = (i + 1) % self.corners.len();
+            let current_edge = LineSegment::new(*edge_start, self.corners[next]);
 
             if LineSegment::intersect(&current_edge, &line_to_infinity) {
+                num_pierced += 1;
+                //println!("Pierced through edge: {:?}", &current_edge);
+
                 /* Check for the edge case, where the point might lie right on an edge of the
                  * polygon.
                  */
-                if super::triplet_orientation(current_edge.start, current_edge.end, point)
+                if super::triplet_orientation(*edge_start, current_edge.end, point)
                     == TripletOrientation::Collinear
                 {
-                    // The point is right on an edge
-                    if current_edge.contains_collinear(point) {
-                        return true;
-                    }
-                } else {
-                    num_intersects += 1;
+                    // The point should be right on an edge.
+                    // XXX: Should this not just be an assertion? It seems that this always has to
+                    // be the case when reaching this part.
+                    return current_edge.contains_collinear(point);
                 }
             }
         }
 
-        num_intersects % 2 == 1
+        //println!("Num pierced is: {}", num_pierced);
+
+        num_pierced % 2 == 1
     }
 
     /// Join this polygon with another, ensuring the area of the two stays the same, minus the
@@ -108,5 +160,6 @@ mod test {
         assert!(polygon.contains(Vec2::new(0.5, 1.5)));
         assert!(!polygon.contains(Vec2::new(-2., 1.9)));
         assert!(!polygon.contains(Vec2::new(0., 3.)));
+        assert!(polygon.contains(Vec2::new(1., 3.)));
     }
 }
-- 
cgit v1.2.3-70-g09d2


From 63f292010e51ad8622cccc4d4b6da887e4eaec47 Mon Sep 17 00:00:00 2001
From: Arne Dußin
Date: Sun, 15 Nov 2020 20:38:14 +0100
Subject: Add intersection point algorithm

---
 src/math/line_segment.rs |  43 ++++++++-
 src/math/polygon.rs      | 228 +++++++++++++++++++++++------------------------
 2 files changed, 153 insertions(+), 118 deletions(-)

(limited to 'src/math/line_segment.rs')

diff --git a/src/math/line_segment.rs b/src/math/line_segment.rs
index 338a681..7dacf54 100644
--- a/src/math/line_segment.rs
+++ b/src/math/line_segment.rs
@@ -1,6 +1,6 @@
-use super::Vec2;
-use alga::general::{ClosedMul, ClosedSub};
-use nalgebra::Scalar;
+use super::{Rect, Vec2};
+use alga::general::{ClosedDiv, ClosedMul, ClosedSub};
+use nalgebra::{RealField, Scalar};
 use num_traits::Zero;
 
 #[derive(Debug)]
@@ -89,6 +89,43 @@ impl<T: Scalar + Copy> LineSegment<T> {
 
         false
     }
+
+    pub fn intersection(line_a: &LineSegment<T>, line_b: &LineSegment<T>) -> Option<Vec2<T>>
+    where
+        T: ClosedSub + ClosedMul + ClosedDiv + Zero + RealField,
+    {
+        // Calculate the differences of each line value from starting point to end point
+        // coordinate.
+        let dax = line_a.end.x - line_a.start.x;
+        let dbx = line_b.end.x - line_b.start.x;
+        let day = line_a.end.y - line_a.start.y;
+        let dby = line_b.end.y - line_b.start.y;
+
+        // Calculate determinant to see, if the lines are parallel or not
+        // TODO: Collinearity check?
+        let d = (dax * dby) - (day * dbx);
+        if d == T::zero() {
+            None
+        } else {
+            let ad = (line_a.start.x * line_a.end.y) - (line_a.start.y * line_a.end.x);
+            let bd = (line_b.start.x * line_b.end.y) - (line_b.start.y * line_b.end.x);
+
+            let out = Vec2::new(((ad * dbx) - (dax * bd)) / d, ((ad * dby) - (day * bd)) / d);
+
+            /* Since the line intersection, not the line segment intersection is calculated, a
+             * bounding check is necessary, to see if the intersection still lies inside both of
+             * the segments. We know it's on the lines, so checking with the lines bounding box is
+             * faster than checking where on the line exactly it would be.
+             */
+            if Rect::bounding_rect(line_a.start, line_a.end).contains(out)
+                && Rect::bounding_rect(line_b.start, line_b.end).contains(out)
+            {
+                Some(out)
+            } else {
+                None
+            }
+        }
+    }
 }
 
 #[derive(PartialEq, Eq)]
diff --git a/src/math/polygon.rs b/src/math/polygon.rs
index 2af669a..e761136 100644
--- a/src/math/polygon.rs
+++ b/src/math/polygon.rs
@@ -1,8 +1,9 @@
-use super::{LineSegment, TripletOrientation, Vec2};
-use alga::general::{ClosedMul, ClosedSub};
-use nalgebra::Scalar;
+use super::Vec2;
+use nalgebra::{ClosedDiv, ClosedMul, ClosedSub, Scalar};
 use num_traits::Zero;
+use std::ops::Neg;
 
+#[derive(Debug)]
 pub struct Polygon<T: Scalar + Copy> {
     pub corners: Vec<Vec2<T>>,
 }
@@ -14,125 +15,90 @@ impl<T: Scalar + Copy> Polygon<T> {
 
     /// Check whether a point is inside a polygon or not. If a point lies on an edge, it also
     /// counts as inside the polygon.
-    pub fn contains(&self, point: Vec2<T>) -> bool
+    pub fn contains_point(&self, p: Vec2<T>) -> bool
     where
-        T: Scalar + Copy + ClosedSub + ClosedMul + PartialOrd + Zero,
+        T: Zero + ClosedSub + ClosedMul + ClosedDiv + Neg<Output = T> + PartialOrd,
     {
-        // Must be at least a triangle to contain anything
-        if self.corners.len() < 3 {
-            return false;
-        }
-
-        // Get the maximum x for the horizontal lines for "Infinity".
-        let mut max_x = self.corners[0].x;
-        for corner in self.corners.iter() {
-            /* Handle the edge case where the point is on a polygon corner, so we don't have to
-             * worry about it later.
-             */
-            if point == *corner {
-                return true;
+        let n = self.corners.len();
+
+        let a = self.corners[n - 1];
+        let mut b = self.corners[n - 2];
+        let mut ax;
+        let mut ay = a.y - p.y;
+        let mut bx = b.x - p.x;
+        let mut by = b.y - p.y;
+
+        let mut lup = by > ay;
+        for i in 0..n {
+            // ax = bx;
+            ay = by;
+            b = self.corners[i];
+            bx = b.x - p.x;
+            by = b.y - p.y;
+
+            if ay == by {
+                continue;
             }
-
-            max_x = super::partial_max(corner.x, max_x);
+            lup = by > ay;
         }
 
-        // The horizontally "infinite" point of the test point.
-        let line_to_infinity = LineSegment::new(point, Vec2::new(max_x, point.y));
-
-        // The number of times the point line to infinity pierced through the polygon hull.
-        let mut num_pierced = 0;
-
-        /* Find edges that should be ignored. These are all edges, where the straight line passes
-         * right through a corner of the edge. It might still be counted as piercing the polygon,
-         * if the next corners lie below and the previous corner lies above the line or vice
-         * versa, but not if both points lie above it or below. The corners are exempted from the
-         * run through the edges. The edges that should be ignored are saved by their starting
-         * corner's index.
-         */
-        let mut ignore_edge = vec![false; self.corners.len()];
-        for i in 0..self.corners.len() {
-            if point.y == self.corners[i].y {
-                // The previous corner index
-                let prev_i = (i + self.corners.len() - 1) % self.corners.len();
-
-                // Make sure the edges containing this corner will be ignored.
-                ignore_edge[i] = true;
-                ignore_edge[prev_i] = true;
-
-                /* Make sure we don't count a pierce double because of consecutive collinear
-                 * points. (Although the polygon really should be sanitised, but you know)
-                 */
-                if self.corners[prev_i].y == point.y {
-                    continue;
-                }
-
-                /* Ignore coming collinear points, since we are only interested in when the
-                 * y-coordinate finally changes
-                 */
-                let mut j = (i + 1) % self.corners.len();
-                while i != j {
-                    if self.corners[j].y != point.y {
-                        if (self.corners[prev_i].y < point.y && self.corners[j].y > point.y)
-                            || (self.corners[prev_i].y > point.y && self.corners[j].y < point.y)
-                        {
-                            /* Only count as pierced when it would be right of the point, otherwise
-                             * it would be like checking to the left of the line that is drawn to
-                             * infinity, which is incorrect
-                             */
-                            if self.corners[i].x >= point.x {
-                                //println!("Pierced through corner {:?}", self.corners[i]);
-                                num_pierced += 1;
-                            }
-                        }
-
-                        break;
-                    }
-
-                    j = (j + 1) % self.corners.len();
-                }
+        let mut depth = 0;
+        for i in 0..n {
+            ax = bx;
+            ay = by;
+            let b = &self.corners[i];
+            bx = b.x - p.x;
+            by = b.y - p.y;
+
+            if ay < T::zero() && by < T::zero() {
+                // both "up" or both "down"
+                continue;
+            }
+            if ay > T::zero() && by > T::zero() {
+                // both "up" or both "down"
+                continue;
+            }
+            if ax < T::zero() && bx < T::zero() {
+                // both points on the left
+                continue;
             }
-        }
 
-        // Find the points where the x-line may actually collide with a proper line.. ^^
-        for (i, edge_start) in self.corners.iter().enumerate() {
-            if ignore_edge[i] {
+            if ay == by && (if ax < bx { ax } else { bx }) <= T::zero() {
+                return true;
+            }
+            if ay == by {
                 continue;
             }
 
-            let next = (i + 1) % self.corners.len();
-            let current_edge = LineSegment::new(*edge_start, self.corners[next]);
-
-            if LineSegment::intersect(&current_edge, &line_to_infinity) {
-                num_pierced += 1;
-                //println!("Pierced through edge: {:?}", &current_edge);
-
-                /* Check for the edge case, where the point might lie right on an edge of the
-                 * polygon.
-                 */
-                if super::triplet_orientation(*edge_start, current_edge.end, point)
-                    == TripletOrientation::Collinear
-                {
-                    // The point should be right on an edge.
-                    // XXX: Should this not just be an assertion? It seems that this always has to
-                    // be the case when reaching this part.
-                    return current_edge.contains_collinear(point);
-                }
+            let lx = ax + (((bx - ax) * -ay) / (by - ay));
+            if lx == T::zero() {
+                // point on edge
+                return true;
+            }
+            if lx > T::zero() {
+                depth += 1;
+            }
+            if ay == T::zero() && lup && by > ay {
+                // hit vertex, both up
+                depth -= 1;
+            }
+            if ay == T::zero() && !lup && by < ay {
+                // hit vertex, both down
+                depth -= 1;
             }
-        }
 
-        //println!("Num pierced is: {}", num_pierced);
+            lup = by > ay;
+        }
 
-        num_pierced % 2 == 1
+        (depth & 1) == 1
     }
 
-    /// Join this polygon with another, ensuring the area of the two stays the same, minus the
-    /// overlap.
-    ///
-    /// # Panics
-    /// If the two polygons do not intersect, it is impossible to unite them into a single polygon,
-    /// in which case this function is bound to fail.
-    // TODO: Create a function that only tries to join the polygons.
-    pub fn unite(&mut self, mut _other: Polygon<T>) {
+    /// Join this polygon with another, ensuring the area of the two stays the same, but the
+    /// overlap is not doubled, but instead joined into one.
+    /// Returns the Polygons themselves, if there is no overlap
+    // TODO: At the moment, counts holes in the middle as a different polygon. These should be
+    // attached to the polygon they are inside of.
+    pub fn unite(self, _other: Polygon<T>) -> Vec<Polygon<T>> {
         unimplemented!()
     }
 }
@@ -153,13 +119,45 @@ mod test {
             Vec2::new(1., 1.),
         ]);
 
-        assert!(!polygon.contains(Vec2::new(1., -2.)));
-        assert!(!polygon.contains(Vec2::new(-1., 0.5)));
-        assert!(polygon.contains(Vec2::new(0., 0.5)));
-        assert!(polygon.contains(Vec2::new(0.5, 1.)));
-        assert!(polygon.contains(Vec2::new(0.5, 1.5)));
-        assert!(!polygon.contains(Vec2::new(-2., 1.9)));
-        assert!(!polygon.contains(Vec2::new(0., 3.)));
-        assert!(polygon.contains(Vec2::new(1., 3.)));
+        assert!(!polygon.contains_point(Vec2::new(1., -2.)));
+        assert!(!polygon.contains_point(Vec2::new(-1., 0.5)));
+        assert!(polygon.contains_point(Vec2::new(0., 0.5)));
+        assert!(polygon.contains_point(Vec2::new(0.5, 1.)));
+        assert!(polygon.contains_point(Vec2::new(0.5, 1.5)));
+        assert!(!polygon.contains_point(Vec2::new(-2., 1.9)));
+        assert!(!polygon.contains_point(Vec2::new(0., 3.)));
+        assert!(polygon.contains_point(Vec2::new(1., 3.)));
+    }
+
+    #[test]
+    fn polygon_union() {
+        let first = Polygon::new(vec![
+            Vec2::new(-2., 1.),
+            Vec2::new(-0.5, 2.5),
+            Vec2::new(2., 2.),
+            Vec2::new(0.5, 1.5),
+            Vec2::new(1., 0.),
+            Vec2::new(-0.5, 1.),
+        ]);
+
+        let second = Polygon::new(vec![
+            Vec2::new(0., 0.),
+            Vec2::new(-2., 2.),
+            Vec2::new(3., 2.),
+            Vec2::new(1.5, 0.),
+        ]);
+
+        let union = first.unite(second);
+        assert_eq!(union.len(), 1);
+        let union = &union[0];
+
+        println!("Union of the two polygons: {:?}", union);
+
+        assert_eq!(union.corners.len(), 10);
+        assert!(union
+            .corners
+            .iter()
+            .find(|&p| p.x == 0. && p.y == 0.)
+            .is_some());
     }
 }
-- 
cgit v1.2.3-70-g09d2


From cc253346c52425bea2ca9919de87e7cfa5ecea97 Mon Sep 17 00:00:00 2001
From: Arne Dußin
Date: Tue, 17 Nov 2020 20:59:50 +0100
Subject: Add polygon graph

---
 src/math/line_segment.rs  |   8 ++
 src/math/mod.rs           |   2 +
 src/math/polygon_graph.rs | 358 ++++++++++++++++++++++++++++++++++++++++++++++
 src/math/vec2.rs          |   2 +-
 4 files changed, 369 insertions(+), 1 deletion(-)
 create mode 100644 src/math/polygon_graph.rs

(limited to 'src/math/line_segment.rs')

diff --git a/src/math/line_segment.rs b/src/math/line_segment.rs
index 7dacf54..1c2c01c 100644
--- a/src/math/line_segment.rs
+++ b/src/math/line_segment.rs
@@ -126,6 +126,14 @@ impl<T: Scalar + Copy> LineSegment<T> {
             }
         }
     }
+
+    /// Checks if this line segment and the other line segment are the same, ignoring the direction
+    /// in which the lines are going, in other words, which of the vectors the line starts at and which
+    /// vector the line ends at is irrelevant.
+    pub fn eq_ignore_dir(&self, other: &LineSegment<T>) -> bool {
+        (self.start == other.start && self.end == other.end)
+            || (self.end == other.start && self.start == other.end)
+    }
 }
 
 #[derive(PartialEq, Eq)]
diff --git a/src/math/mod.rs b/src/math/mod.rs
index bfc2ab4..afbc4a3 100644
--- a/src/math/mod.rs
+++ b/src/math/mod.rs
@@ -1,10 +1,12 @@
 pub mod line_segment;
 pub mod polygon;
+pub mod polygon_graph;
 pub mod rect;
 pub mod vec2;
 
 pub use self::line_segment::*;
 pub use self::polygon::*;
+pub use self::polygon_graph::*;
 pub use self::rect::*;
 pub use self::vec2::*;
 
diff --git a/src/math/polygon_graph.rs b/src/math/polygon_graph.rs
new file mode 100644
index 0000000..3acd46f
--- /dev/null
+++ b/src/math/polygon_graph.rs
@@ -0,0 +1,358 @@
+use super::{LineSegment, Polygon, Vec2};
+use nalgebra::Scalar;
+use std::cmp::{Ordering, PartialOrd};
+
+struct Node<T: Scalar + Copy> {
+    pub vec: Vec2<T>,
+    pub adjacent: Vec<Vec2<T>>,
+}
+
+/// An undirected graph, that is optimised for polygon edge operations. Since edges of a polygon
+/// are an undirected graph, this structure also holds both directions. This makes it rather space
+/// inefficient, but makes edge operations rather swift. ß
+pub struct PolygonGraph<T: Scalar + Copy + PartialOrd> {
+    /// The nodes of the graph, together with their adjacency list.
+    nodes: Vec<Node<T>>,
+}
+
+// Helper functions to find nodes/vecs in sorted fields, so It doesn't always have to be written
+// out.
+#[inline]
+fn find_vec2<T: Scalar + Copy + PartialOrd>(
+    field: &[Vec2<T>],
+    lookup: &Vec2<T>,
+) -> Result<usize, usize> {
+    field.binary_search_by(|n| n.partial_cmp(lookup).unwrap_or(Ordering::Greater))
+}
+#[inline]
+fn find_node<T: Scalar + Copy + PartialOrd>(
+    field: &[Node<T>],
+    lookup: &Vec2<T>,
+) -> Result<usize, usize> {
+    field.binary_search_by(|n| n.vec.partial_cmp(lookup).unwrap_or(Ordering::Greater))
+}
+
+impl<T: Scalar + Copy + PartialOrd> PolygonGraph<T> {
+    /// Create a new, empty polygon graph.
+    pub fn new() -> Self {
+        Self { nodes: Vec::new() }
+    }
+
+    /// Count the number of edges in the graph. Internally, for each two connected points there are
+    /// two edges, but this returns the amount of polygon edges.
+    pub fn num_edges(&self) -> usize {
+        let mut num_edges = 0;
+        for node in &self.nodes {
+            for _ in &node.adjacent {
+                num_edges += 1;
+            }
+        }
+
+        assert!(num_edges % 2 == 0);
+        num_edges / 2
+    }
+
+    /// Count the number of nodes in this graph. If this graph consists of multiple polygons, this
+    /// can be different than the amount of corners, since corners with the same position are only
+    /// counted once.
+    pub fn num_nodes(&self) -> usize {
+        self.nodes.len()
+    }
+
+    /// Checks if there is an edge between the two given vectors. Is commutative in respect to the
+    /// two arguments.
+    pub fn has_edge(&self, from: &Vec2<T>, to: &Vec2<T>) -> bool {
+        // Binary search the starting and then the end node.
+        if let Ok(from) = find_node(&self.nodes, from) {
+            if let Ok(_) = find_vec2(&self.nodes[from].adjacent, to) {
+                true
+            } else {
+                false
+            }
+        } else {
+            false
+        }
+    }
+
+    // Helper function to add the edge into the internal graph representation for one side only.
+    // Since to the outside the graph should always be undirected, this must be private.
+    fn add_edge_onesided(&mut self, from: &Vec2<T>, to: &Vec2<T>) -> bool {
+        match find_node(&self.nodes, from) {
+            Ok(pos) => match find_vec2(&self.nodes[pos].adjacent, to) {
+                Ok(_) => return false,
+                Err(i) => self.nodes[pos].adjacent.insert(i, *to),
+            },
+            Err(pos) => self.nodes.insert(
+                pos,
+                Node {
+                    vec: *from,
+                    adjacent: vec![*to],
+                },
+            ),
+        }
+
+        true
+    }
+
+    /// Add an edge between the given vectors. If the edge already exists, it does nothing and
+    /// returns false, otherwise it returns true after addition.
+    pub fn add_edge(&mut self, from: &Vec2<T>, to: &Vec2<T>) -> bool {
+        if !self.add_edge_onesided(from, to) {
+            return false;
+        }
+
+        let back_edge_succ = self.add_edge_onesided(to, from);
+        assert!(back_edge_succ);
+
+        true
+    }
+
+    // Helper function to remove the edge in the internal graph representation for one side only.
+    // Since to the outside the graph should always be undirected, this must be private.
+    fn remove_edge_onesided(&mut self, from: &Vec2<T>, to: &Vec2<T>) -> bool {
+        if let Ok(from) = find_node(&self.nodes, from) {
+            if let Ok(to) = find_vec2(&self.nodes[from].adjacent, to) {
+                // Remove the edge from the vector.
+                self.nodes[from].adjacent.remove(to);
+
+                // If the node has no adjacent nodes anymore, remove it entirely.
+                if self.nodes[from].adjacent.is_empty() {
+                    self.nodes.remove(from);
+                }
+
+                true
+            } else {
+                false
+            }
+        } else {
+            false
+        }
+    }
+
+    /// Remove an edge between the given vectors. If there is no edge between them, it does nothing
+    /// and returns false, otherwise it returns true after deletion.
+    pub fn remove_edge(&mut self, from: &Vec2<T>, to: &Vec2<T>) -> bool {
+        if !self.remove_edge_onesided(from, to) {
+            return false;
+        }
+
+        let back_edge_succ = self.remove_edge_onesided(to, from);
+        assert!(back_edge_succ);
+
+        true
+    }
+
+    /// Constructs a new PolygonGraph from the provided polygon. Adds a node for every corner and
+    /// an edge to all connected corners (which should be exactly two, for regular polygons)
+    pub fn from_polygon(polygon: &Polygon<T>) -> Self {
+        let mut graph = PolygonGraph {
+            nodes: Vec::with_capacity(polygon.corners.len()),
+        };
+
+        graph.add_all(polygon);
+        graph
+    }
+
+    /// Add all edges of the provided polygon into the graph. Requires roughly double as much space
+    /// as the normal polygon.
+    pub fn add_all(&mut self, polygon: &Polygon<T>) {
+        for i in 0..polygon.corners.len() {
+            self.add_edge(
+                &polygon.corners[i],
+                &polygon.corners[(i + 1) % polygon.corners.len()],
+            );
+        }
+    }
+
+    /// Calculates all points where the graph edges intersect with one another. It then adds them
+    /// into the adjacency list such that the intersection point lies between the nodes of the
+    /// lines.
+    pub fn intersect_self(&mut self) {
+        // Find all intersections with all other edges.
+        for node in &self.nodes {
+            for to in &node.adjacent {
+                let current_segment = LineSegment::new(node.vec, *to);
+
+                for node in &self.nodes {
+                    for to in &node.adjacent {
+                        let compare_segment = LineSegment::new(node.vec, *to);
+                        if current_segment.eq_ignore_dir(&compare_segment) {
+                            continue;
+                        }
+                    }
+                }
+            }
+        }
+    }
+
+    /// Finds the minimal polygon that could hold this graph together, i.e. could contain the
+    /// entire graph, but with the minimal amount of space. It may however still contain extra
+    /// corner points, meaning an extra edge for three collinear points on this edge, that can be
+    /// cleaned up.
+    pub fn bounding_polygon(&self) -> Polygon<T> {
+        unimplemented!()
+    }
+}
+
+impl<T: Scalar + Copy + PartialOrd> Default for PolygonGraph<T> {
+    fn default() -> Self {
+        Self::new()
+    }
+}
+
+#[cfg(test)]
+mod test {
+    use super::*;
+
+    #[test]
+    fn from_polygon() {
+        let a = Vec2::new(0., 0.);
+        let b = Vec2::new(0., 1.);
+        let c = Vec2::new(0.5, 1.);
+
+        let triangle = Polygon::new(vec![a, b, c]);
+
+        let graph = PolygonGraph::from_polygon(&triangle);
+        assert_eq!(graph.num_edges(), 3);
+
+        assert!(graph.has_edge(&a, &b));
+        assert!(graph.has_edge(&b, &a));
+
+        assert!(graph.has_edge(&a, &c));
+        assert!(graph.has_edge(&c, &a));
+
+        assert!(graph.has_edge(&b, &c));
+        assert!(graph.has_edge(&c, &b));
+    }
+
+    #[test]
+    fn add_all() {
+        let top_left = Vec2::new(0., 0.);
+        let top_right = Vec2::new(1., 0.);
+        let bot_left = Vec2::new(0., 1.);
+        let bot_right = Vec2::new(1., 1.);
+
+        let triangle = Polygon::new(vec![top_left, bot_right, top_right]);
+
+        let square = Polygon::new(vec![bot_left, bot_right, top_right, top_left]);
+
+        let mut graph = PolygonGraph::new();
+        graph.add_all(&triangle);
+        graph.add_all(&square);
+
+        assert_eq!(graph.num_edges(), 5);
+        assert_eq!(graph.num_nodes(), 4);
+
+        assert!(graph.has_edge(&top_left, &top_right));
+        assert!(graph.has_edge(&top_right, &top_left));
+
+        assert!(graph.has_edge(&top_left, &bot_left));
+        assert!(graph.has_edge(&bot_left, &top_left));
+
+        assert!(graph.has_edge(&bot_left, &bot_right));
+        assert!(graph.has_edge(&bot_right, &bot_left));
+
+        assert!(graph.has_edge(&bot_right, &top_right));
+        assert!(graph.has_edge(&top_right, &bot_right));
+
+        assert!(graph.has_edge(&top_left, &bot_right));
+        assert!(graph.has_edge(&bot_right, &top_left));
+    }
+
+    #[test]
+    fn intersect_self() {
+        let first = Polygon::new(vec![
+            Vec2::new(0., 0.),
+            Vec2::new(0., 2.),
+            Vec2::new(2., 2.),
+            Vec2::new(3., 1.),
+            Vec2::new(2., 0.),
+        ]);
+
+        let second = Polygon::new(vec![
+            Vec2::new(2.5, -0.5),
+            Vec2::new(0., 2.),
+            Vec2::new(2., 2.),
+            Vec2::new(2., 0.5),
+            Vec2::new(2.5, 0.),
+        ]);
+
+        let mut graph = PolygonGraph::from_polygon(&first);
+        graph.add_all(&second);
+
+        graph.intersect_self();
+
+        assert_eq!(graph.num_nodes(), 9);
+        assert_eq!(graph.num_edges(), 12);
+
+        assert!(graph.has_edge(&Vec2::new(2., 0.), &Vec2::new(2.25, 0.25)));
+        assert!(graph.has_edge(&Vec2::new(3., 1.), &Vec2::new(2.25, 0.25)));
+        assert!(!graph.has_edge(&Vec2::new(2., 0.), &Vec2::new(3., 1.)));
+        assert!(graph.has_edge(&Vec2::new(0., 2.), &Vec2::new(2., 2.)));
+        assert!(graph.has_edge(&Vec2::new(0., 2.), &Vec2::new(2., 0.)));
+        assert!(!graph.has_edge(&Vec2::new(0., 2.), &Vec2::new(2.5, -0.5)));
+    }
+
+    #[test]
+    fn bounding_polygon() {
+        let first = Polygon::new(vec![
+            Vec2::new(0., 0.),
+            Vec2::new(0., 2.),
+            Vec2::new(2., 2.),
+            Vec2::new(3., 1.),
+            Vec2::new(2., 0.),
+        ]);
+
+        let second = Polygon::new(vec![
+            Vec2::new(2.5, -0.5),
+            Vec2::new(0., 2.),
+            Vec2::new(2., 2.),
+            Vec2::new(2., 0.5),
+            Vec2::new(2.5, 0.),
+        ]);
+
+        let mut graph = PolygonGraph::from_polygon(&first);
+        graph.add_all(&second);
+
+        let bounding = graph.bounding_polygon();
+
+        let num_corners = 8;
+        assert_eq!(bounding.corners.len(), num_corners);
+
+        // Run around the polygon to see if it was constructed correctly.
+        let start_i = bounding
+            .corners
+            .iter()
+            .position(|&x| x == Vec2::new(0., 0.))
+            .expect("Starting vector does not exist in polygon.");
+
+        assert_eq!(
+            bounding.corners[(start_i + 1) % num_corners],
+            Vec2::new(0., 2.)
+        );
+        assert_eq!(
+            bounding.corners[(start_i + 2) % num_corners],
+            Vec2::new(2., 2.)
+        );
+        assert_eq!(
+            bounding.corners[(start_i + 3) % num_corners],
+            Vec2::new(3., 1.)
+        );
+        assert_eq!(
+            bounding.corners[(start_i + 4) % num_corners],
+            Vec2::new(2.25, 0.25)
+        );
+        assert_eq!(
+            bounding.corners[(start_i + 5) % num_corners],
+            Vec2::new(2.5, 0.0)
+        );
+        assert_eq!(
+            bounding.corners[(start_i + 6) % num_corners],
+            Vec2::new(2.5, -0.5)
+        );
+        assert_eq!(
+            bounding.corners[(start_i + 7) % num_corners],
+            Vec2::new(2., 0.)
+        );
+    }
+}
diff --git a/src/math/vec2.rs b/src/math/vec2.rs
index 2b33f7b..cd38889 100644
--- a/src/math/vec2.rs
+++ b/src/math/vec2.rs
@@ -8,7 +8,7 @@ use std::convert::{From, Into};
 use std::ops::{Add, AddAssign, Div, Mul, MulAssign, Neg, Sub, SubAssign};
 use std::{fmt, mem};
 
-#[derive(Clone, Copy, Debug, Default, PartialEq, Serialize, Deserialize, Eq)]
+#[derive(Clone, Copy, Debug, Default, PartialEq, Serialize, Deserialize, Eq, Hash)]
 pub struct Vec2<T: Scalar + Copy> {
     pub x: T,
     pub y: T,
-- 
cgit v1.2.3-70-g09d2


From af4f4b91d013384a4af82f48d768ac041662eb20 Mon Sep 17 00:00:00 2001
From: Arne Dußin
Date: Wed, 18 Nov 2020 13:01:10 +0100
Subject: Add self intersection for polygon graphs

---
 src/math/line_segment.rs  | 69 ++++++++++++++++++++++++++++++----
 src/math/polygon_graph.rs | 95 ++++++++++++++++++++++++++++++++++++++++-------
 2 files changed, 143 insertions(+), 21 deletions(-)

(limited to 'src/math/line_segment.rs')

diff --git a/src/math/line_segment.rs b/src/math/line_segment.rs
index 1c2c01c..2b67783 100644
--- a/src/math/line_segment.rs
+++ b/src/math/line_segment.rs
@@ -2,8 +2,9 @@ use super::{Rect, Vec2};
 use alga::general::{ClosedDiv, ClosedMul, ClosedSub};
 use nalgebra::{RealField, Scalar};
 use num_traits::Zero;
+use std::cmp::Ordering;
 
-#[derive(Debug)]
+#[derive(Debug, Clone)]
 pub struct LineSegment<T: Scalar + Copy> {
     pub start: Vec2<T>,
     pub end: Vec2<T>,
@@ -23,8 +24,8 @@ impl<T: Scalar + Copy> LineSegment<T> {
         T: PartialOrd,
     {
         point.x <= super::partial_max(self.start.x, self.end.x)
-            && point.x >= super::partial_max(self.start.x, self.end.x)
-            && point.y <= super::partial_min(self.start.y, self.end.y)
+            && point.x >= super::partial_min(self.start.x, self.end.x)
+            && point.y <= super::partial_max(self.start.y, self.end.y)
             && point.y >= super::partial_min(self.start.y, self.end.y)
     }
 
@@ -96,10 +97,10 @@ impl<T: Scalar + Copy> LineSegment<T> {
     {
         // Calculate the differences of each line value from starting point to end point
         // coordinate.
-        let dax = line_a.end.x - line_a.start.x;
-        let dbx = line_b.end.x - line_b.start.x;
-        let day = line_a.end.y - line_a.start.y;
-        let dby = line_b.end.y - line_b.start.y;
+        let dax = line_a.start.x - line_a.end.x;
+        let dbx = line_b.start.x - line_b.end.x;
+        let day = line_a.start.y - line_a.end.y;
+        let dby = line_b.start.y - line_b.end.y;
 
         // Calculate determinant to see, if the lines are parallel or not
         // TODO: Collinearity check?
@@ -127,6 +128,42 @@ impl<T: Scalar + Copy> LineSegment<T> {
         }
     }
 
+    /// Find all segments, into which this LineSegment would be splitted, when the points provided
+    /// would cut the segment. The points must be on the line, otherwise this does not make sense.
+    /// Also, no segments of length zero (start point = end point) will be created.
+    pub fn segments(&self, split_points: &[Vec2<T>]) -> Vec<Vec2<T>>
+    where
+        T: ClosedSub + ClosedMul + PartialOrd + Zero,
+    {
+        // Make sure all segments are collinear with the segment and actually on it
+        assert_eq!(
+            split_points
+                .iter()
+                .find(|&p| triplet_orientation(self.start, self.end, *p)
+                    != TripletOrientation::Collinear),
+            None
+        );
+        assert_eq!(
+            split_points.iter().find(|&p| !self.contains_collinear(*p)),
+            None
+        );
+
+        let mut segments = Vec::with_capacity(split_points.len() + 2);
+        segments.push(self.start);
+        segments.push(self.end);
+        for point in split_points {
+            segments.push(*point);
+        }
+
+        segments.sort_unstable_by(|a, b| a.partial_cmp(b).unwrap_or(Ordering::Greater));
+        segments.dedup();
+        if self.start > self.end {
+            segments.reverse();
+        }
+
+        segments
+    }
+
     /// Checks if this line segment and the other line segment are the same, ignoring the direction
     /// in which the lines are going, in other words, which of the vectors the line starts at and which
     /// vector the line ends at is irrelevant.
@@ -160,3 +197,21 @@ where
         _ => TripletOrientation::Collinear,
     }
 }
+
+#[cfg(test)]
+mod test {
+    use super::*;
+
+    #[test]
+    fn contains_collinear() {
+        let segment = LineSegment::new(Vec2::new(0., 0.), Vec2::new(2., 2.));
+
+        assert!(segment.contains_collinear(Vec2::new(0., 0.)));
+        assert!(segment.contains_collinear(Vec2::new(1., 1.)));
+        assert!(segment.contains_collinear(Vec2::new(2., 2.)));
+
+        assert!(!segment.contains_collinear(Vec2::new(-1., -1.)));
+        assert!(!segment.contains_collinear(Vec2::new(3., 3.)));
+        assert!(!segment.contains_collinear(Vec2::new(-1., 0.)));
+    }
+}
diff --git a/src/math/polygon_graph.rs b/src/math/polygon_graph.rs
index 3acd46f..43c3bc4 100644
--- a/src/math/polygon_graph.rs
+++ b/src/math/polygon_graph.rs
@@ -1,20 +1,26 @@
 use super::{LineSegment, Polygon, Vec2};
-use nalgebra::Scalar;
+use nalgebra::{RealField, Scalar};
 use std::cmp::{Ordering, PartialOrd};
 
+#[derive(Debug)]
 struct Node<T: Scalar + Copy> {
     pub vec: Vec2<T>,
     pub adjacent: Vec<Vec2<T>>,
 }
 
+struct EdgeIterator<'a, T: Scalar + Copy> {
+    nodes: &'a Vec<Node<T>>,
+    pos: (usize, usize),
+}
+
 /// An undirected graph, that is optimised for polygon edge operations. Since edges of a polygon
 /// are an undirected graph, this structure also holds both directions. This makes it rather space
 /// inefficient, but makes edge operations rather swift. ß
+#[derive(Debug)]
 pub struct PolygonGraph<T: Scalar + Copy + PartialOrd> {
     /// The nodes of the graph, together with their adjacency list.
     nodes: Vec<Node<T>>,
 }
-
 // Helper functions to find nodes/vecs in sorted fields, so It doesn't always have to be written
 // out.
 #[inline]
@@ -32,6 +38,35 @@ fn find_node<T: Scalar + Copy + PartialOrd>(
     field.binary_search_by(|n| n.vec.partial_cmp(lookup).unwrap_or(Ordering::Greater))
 }
 
+impl<'a, T: Scalar + Copy> EdgeIterator<'a, T> {
+    pub fn new(nodes: &'a Vec<Node<T>>) -> Self {
+        Self { nodes, pos: (0, 0) }
+    }
+}
+
+impl<'a, T: Scalar + Copy> Iterator for EdgeIterator<'a, T> {
+    type Item = LineSegment<T>;
+
+    fn next(&mut self) -> Option<Self::Item> {
+        // Try to find the element in the nodes vector, if it exists.
+        if let Some(node) = self.nodes.get(self.pos.0) {
+            let end = node.adjacent[self.pos.1];
+
+            // Advance the iterator to the next possible element
+            if self.pos.1 + 1 < node.adjacent.len() {
+                self.pos.1 += 1;
+            } else {
+                self.pos.1 = 0;
+                self.pos.0 += 1;
+            }
+
+            Some(LineSegment::new(node.vec, end))
+        } else {
+            None
+        }
+    }
+}
+
 impl<T: Scalar + Copy + PartialOrd> PolygonGraph<T> {
     /// Create a new, empty polygon graph.
     pub fn new() -> Self {
@@ -167,22 +202,51 @@ impl<T: Scalar + Copy + PartialOrd> PolygonGraph<T> {
     /// Calculates all points where the graph edges intersect with one another. It then adds them
     /// into the adjacency list such that the intersection point lies between the nodes of the
     /// lines.
-    pub fn intersect_self(&mut self) {
+    pub fn intersect_self(&mut self)
+    where
+        T: RealField,
+    {
         // Find all intersections with all other edges.
-        for node in &self.nodes {
-            for to in &node.adjacent {
-                let current_segment = LineSegment::new(node.vec, *to);
-
-                for node in &self.nodes {
-                    for to in &node.adjacent {
-                        let compare_segment = LineSegment::new(node.vec, *to);
-                        if current_segment.eq_ignore_dir(&compare_segment) {
-                            continue;
-                        }
-                    }
+        let mut to_delete: Vec<LineSegment<T>> = Vec::new();
+        let mut to_add: Vec<(Vec2<T>, Vec2<T>)> = Vec::new();
+        for segment in EdgeIterator::new(&self.nodes) {
+            /* Save all intersections of this line with any other line, and the line that it's
+             * intersecting with.
+             */
+            let mut intersections: Vec<Vec2<T>> = Vec::new();
+            for compare_segment in EdgeIterator::new(&self.nodes) {
+                if segment.eq_ignore_dir(&compare_segment) {
+                    continue;
                 }
+
+                if let Some(intersection) = LineSegment::intersection(&segment, &compare_segment) {
+                    println!(
+                        "Found intersection between {:?} and {:?}: {:?}",
+                        &segment, &compare_segment, &intersection
+                    );
+                    intersections.push(intersection);
+                }
+            }
+
+            if intersections.is_empty() {
+                continue;
+            }
+
+            to_delete.push(segment.clone());
+
+            // Safe, since at least the line segment itself is represented.
+            let segments = segment.segments(&intersections);
+            for i in 1..segments.len() {
+                to_add.push((segments[i - 1], segments[i]));
             }
         }
+
+        for segment in to_delete {
+            self.remove_edge(&segment.start, &segment.end);
+        }
+        for (start, end) in to_add {
+            self.add_edge(&start, &end);
+        }
     }
 
     /// Finds the minimal polygon that could hold this graph together, i.e. could contain the
@@ -282,6 +346,9 @@ mod test {
 
         graph.intersect_self();
 
+        println!("Intersected graph:");
+        println!("{:#?}", &graph);
+
         assert_eq!(graph.num_nodes(), 9);
         assert_eq!(graph.num_edges(), 12);
 
-- 
cgit v1.2.3-70-g09d2


From 8d99501918393ea0c046d721e6fa6015fe43b70f Mon Sep 17 00:00:00 2001
From: Arne Dußin
Date: Wed, 18 Nov 2020 21:54:41 +0100
Subject: Implement bounding box function

---
 src/math/line_segment.rs  | 74 +++++++++++++++++++++++++++++++++++++++++++++--
 src/math/polygon.rs       | 20 ++++++++-----
 src/math/polygon_graph.rs | 65 +++++++++++++++++++++++++++++++++--------
 3 files changed, 138 insertions(+), 21 deletions(-)

(limited to 'src/math/line_segment.rs')

diff --git a/src/math/line_segment.rs b/src/math/line_segment.rs
index 2b67783..e6ca70f 100644
--- a/src/math/line_segment.rs
+++ b/src/math/line_segment.rs
@@ -192,15 +192,57 @@ where
      * they are the same it follows that the three points are collinear.
      */
     match (b.y - a.y) * (c.x - b.x) - (b.x - a.x) * (c.y - b.y) {
-        q if q > T::zero() => TripletOrientation::Clockwise,
-        q if q < T::zero() => TripletOrientation::Counterclockwise,
+        q if q > T::zero() => TripletOrientation::Counterclockwise,
+        q if q < T::zero() => TripletOrientation::Clockwise,
         _ => TripletOrientation::Collinear,
     }
 }
 
+/// Calculate the angle between three points. Guaranteed to have 0 <= angle < 2Pi. The angle is
+/// calculated as the angle between a line from a to b and then from b to c, increasing
+/// counterclockwise.
+///
+/// # Panics
+/// If the length from a to b or the length from b to c is zero.
+pub(crate) fn triplet_angle<T>(a: Vec2<T>, b: Vec2<T>, c: Vec2<T>) -> T
+where
+    T: Scalar + Copy + ClosedSub + RealField + Zero,
+{
+    assert!(a != b);
+    assert!(b != c);
+
+    // Handle the extreme 0 and 180 degree cases
+    let orientation = triplet_orientation(a, b, c);
+    if orientation == TripletOrientation::Collinear {
+        if LineSegment::new(a, b).contains_collinear(c)
+            || LineSegment::new(b, c).contains_collinear(a)
+        {
+            return T::zero();
+        } else {
+            return T::pi();
+        }
+    }
+
+    // Calculate the vectors out of the points
+    let ba = a - b;
+    let bc = c - b;
+
+    // Calculate the angle between 0 and Pi.
+    let angle = ((ba * bc) / (ba.length() * bc.length())).acos();
+
+    // Make angle into a full circle angle by looking at the orientation of the triplet.
+    let angle = match orientation {
+        TripletOrientation::Counterclockwise => T::pi() + (T::pi() - angle),
+        _ => angle,
+    };
+
+    angle
+}
+
 #[cfg(test)]
 mod test {
     use super::*;
+    use std::f64::consts::PI;
 
     #[test]
     fn contains_collinear() {
@@ -214,4 +256,32 @@ mod test {
         assert!(!segment.contains_collinear(Vec2::new(3., 3.)));
         assert!(!segment.contains_collinear(Vec2::new(-1., 0.)));
     }
+
+    #[test]
+    fn triplet_angle() {
+        assert_eq!(
+            super::triplet_angle(Vec2::new(0., 0.), Vec2::new(0., -1.), Vec2::new(-1., -1.)),
+            1.5 * PI
+        );
+        assert_eq!(
+            super::triplet_angle(Vec2::new(2., 2.), Vec2::new(0., 0.), Vec2::new(-1., -1.)),
+            PI
+        );
+        assert_eq!(
+            super::triplet_angle(Vec2::new(3., 1.), Vec2::new(0., 0.), Vec2::new(6., 2.)),
+            0.
+        );
+        assert_eq!(
+            super::triplet_angle(Vec2::new(3., 1.), Vec2::new(0., 0.), Vec2::new(-6., -2.)),
+            PI
+        );
+        assert_eq!(
+            super::triplet_angle(Vec2::new(0., 1.), Vec2::new(0., 2.), Vec2::new(-3., 2.)),
+            0.5 * PI
+        );
+        assert_eq!(
+            super::triplet_angle(Vec2::new(3., 1.), Vec2::new(2., 2.), Vec2::new(0., 2.)),
+            0.75 * PI
+        );
+    }
 }
diff --git a/src/math/polygon.rs b/src/math/polygon.rs
index e761136..5711049 100644
--- a/src/math/polygon.rs
+++ b/src/math/polygon.rs
@@ -1,9 +1,10 @@
-use super::Vec2;
-use nalgebra::{ClosedDiv, ClosedMul, ClosedSub, Scalar};
+use super::{PolygonGraph, Vec2};
+use nalgebra::{ClosedDiv, ClosedMul, ClosedSub, RealField, Scalar};
 use num_traits::Zero;
 use std::ops::Neg;
 
 #[derive(Debug)]
+// TODO: Support polygons with holes
 pub struct Polygon<T: Scalar + Copy> {
     pub corners: Vec<Vec2<T>>,
 }
@@ -96,10 +97,15 @@ impl<T: Scalar + Copy> Polygon<T> {
     /// Join this polygon with another, ensuring the area of the two stays the same, but the
     /// overlap is not doubled, but instead joined into one.
     /// Returns the Polygons themselves, if there is no overlap
-    // TODO: At the moment, counts holes in the middle as a different polygon. These should be
-    // attached to the polygon they are inside of.
-    pub fn unite(self, _other: Polygon<T>) -> Vec<Polygon<T>> {
-        unimplemented!()
+    pub fn unite(self, other: Polygon<T>) -> Vec<Polygon<T>>
+    where
+        T: RealField,
+    {
+        let mut graph = PolygonGraph::from_polygon(&self);
+        graph.add_all(&other);
+
+        // TODO: Make bounding box support multiple polygons
+        vec![graph.bounding_polygon()]
     }
 }
 
@@ -153,7 +159,7 @@ mod test {
 
         println!("Union of the two polygons: {:?}", union);
 
-        assert_eq!(union.corners.len(), 10);
+        assert_eq!(union.corners.len(), 11);
         assert!(union
             .corners
             .iter()
diff --git a/src/math/polygon_graph.rs b/src/math/polygon_graph.rs
index 43c3bc4..7721cbf 100644
--- a/src/math/polygon_graph.rs
+++ b/src/math/polygon_graph.rs
@@ -220,10 +220,6 @@ impl<T: Scalar + Copy + PartialOrd> PolygonGraph<T> {
                 }
 
                 if let Some(intersection) = LineSegment::intersection(&segment, &compare_segment) {
-                    println!(
-                        "Found intersection between {:?} and {:?}: {:?}",
-                        &segment, &compare_segment, &intersection
-                    );
                     intersections.push(intersection);
                 }
             }
@@ -253,8 +249,52 @@ impl<T: Scalar + Copy + PartialOrd> PolygonGraph<T> {
     /// entire graph, but with the minimal amount of space. It may however still contain extra
     /// corner points, meaning an extra edge for three collinear points on this edge, that can be
     /// cleaned up.
-    pub fn bounding_polygon(&self) -> Polygon<T> {
-        unimplemented!()
+    pub fn bounding_polygon(mut self) -> Polygon<T>
+    where
+        T: RealField,
+    {
+        assert!(self.num_nodes() >= 3);
+        self.intersect_self();
+
+        /* Start with the top-left node. Since the nodes are always sorted by y over x from top to
+         * bottom and left to right, this is the very first element in the vector. This is also a
+         * corner, because for such a node to be enveloped, there would have to be a node further
+         * to the top, in which case that node would have been selected.
+         */
+        let mut current_node = &self.nodes[0];
+        // Pretend we're coming from the top to start in the right direction.
+        let mut last_vec = current_node.vec - Vec2::new(T::zero(), T::one());
+        let mut bounding_polygon = Polygon::new(vec![current_node.vec]);
+        loop {
+            /* Find the next point by choosing the one with the greatest angle. This means we are
+             * "bending" to the leftmost edge at each step. Since we are going around the polygon
+             * in a clockwise direction, this yields the hull around the polygon.
+             * NOTE: Going left is just as viable, but we would have to handle the case where the
+             * algorithm would try to go back because the edge back has 0 degrees, which would be
+             * always preferred. Going right makes going back the absolute worst option.
+             */
+            let next_corner = current_node
+                .adjacent
+                .iter()
+                .max_by(|&a, &b| {
+                    super::triplet_angle(last_vec, current_node.vec, *a)
+                        .partial_cmp(&super::triplet_angle(last_vec, current_node.vec, *b))
+                        .unwrap_or(Ordering::Equal)
+                })
+                .expect("Adjacency list is empty. The polygon has an open edge (is broken)");
+
+            // When we have come back to the start, the traversal is completed
+            if *next_corner == bounding_polygon.corners[0] {
+                break;
+            }
+
+            bounding_polygon.corners.push(*next_corner);
+            last_vec = current_node.vec;
+            current_node = &self.nodes[find_node(&self.nodes, &next_corner)
+                .expect("Failure to find node that should be inside list.")];
+        }
+
+        bounding_polygon
     }
 }
 
@@ -356,6 +396,7 @@ mod test {
         assert!(graph.has_edge(&Vec2::new(3., 1.), &Vec2::new(2.25, 0.25)));
         assert!(!graph.has_edge(&Vec2::new(2., 0.), &Vec2::new(3., 1.)));
         assert!(graph.has_edge(&Vec2::new(0., 2.), &Vec2::new(2., 2.)));
+        assert!(graph.has_edge(&Vec2::new(2., 2.), &Vec2::new(0., 2.)));
         assert!(graph.has_edge(&Vec2::new(0., 2.), &Vec2::new(2., 0.)));
         assert!(!graph.has_edge(&Vec2::new(0., 2.), &Vec2::new(2.5, -0.5)));
     }
@@ -395,15 +436,15 @@ mod test {
 
         assert_eq!(
             bounding.corners[(start_i + 1) % num_corners],
-            Vec2::new(0., 2.)
+            Vec2::new(2., 0.)
         );
         assert_eq!(
             bounding.corners[(start_i + 2) % num_corners],
-            Vec2::new(2., 2.)
+            Vec2::new(2.5, -0.5)
         );
         assert_eq!(
             bounding.corners[(start_i + 3) % num_corners],
-            Vec2::new(3., 1.)
+            Vec2::new(2.5, 0.0)
         );
         assert_eq!(
             bounding.corners[(start_i + 4) % num_corners],
@@ -411,15 +452,15 @@ mod test {
         );
         assert_eq!(
             bounding.corners[(start_i + 5) % num_corners],
-            Vec2::new(2.5, 0.0)
+            Vec2::new(3., 1.)
         );
         assert_eq!(
             bounding.corners[(start_i + 6) % num_corners],
-            Vec2::new(2.5, -0.5)
+            Vec2::new(2., 2.)
         );
         assert_eq!(
             bounding.corners[(start_i + 7) % num_corners],
-            Vec2::new(2., 0.)
+            Vec2::new(0., 2.)
         );
     }
 }
-- 
cgit v1.2.3-70-g09d2