Refactor and prepare for closed shapes

geo/src/algorithm/offset/cross_product.rs

@@ -1,7 +1,7 @@
 use crate::CoordFloat;
 use geo_types::Coord;
 
-/// 2D "Cross Product"
+/// The signed magnitude of the 3D "Cross Product" assuming z ordinates are zero
 ///
 /// > Note: `cross_prod` is already defined on `Point`... but that it seems to be
 /// >       some other operation on 3 points
@@ -40,7 +40,7 @@ use geo_types::Coord;
 ///         = |               (a_x·b_y       - a_y·b_x)·k |
 ///         =                  a_x·b_y       - a_y·b_x
 /// ```
-pub(super) fn cross_product<T>(left: Coord<T>, right: Coord<T>) -> T
+pub(super) fn cross_product_2d<T>(left: Coord<T>, right: Coord<T>) -> T
 where
     T: CoordFloat,
 {
@@ -53,7 +53,7 @@ mod test {
     use crate::Coord;
 
     // private imports
-    use super::cross_product;
+    use super::cross_product_2d;
 
     #[test]
     fn test_cross_product() {
@@ -65,9 +65,9 @@ mod test {
         let ac = c - a;
 
         // expect the area of the parallelogram
-        assert_eq!(cross_product(ac, ab), 1f64);
+        assert_eq!(cross_product_2d(ac, ab), 1f64);
         // expect swapping will result in negative
-        assert_eq!(cross_product(ab, ac), -1f64);
+        assert_eq!(cross_product_2d(ab, ac), -1f64);
 
         // Add skew
         let a = Coord { x: 0f64, y: 0f64 };
@@ -78,8 +78,8 @@ mod test {
         let ac = c - a;
 
         // expect the area of the parallelogram
-        assert_eq!(cross_product(ac, ab), 1f64);
+        assert_eq!(cross_product_2d(ac, ab), 1f64);
         // expect swapping will result in negative
-        assert_eq!(cross_product(ab, ac), -1f64);
+        assert_eq!(cross_product_2d(ab, ac), -1f64);
     }
 }

geo/src/algorithm/offset/line_intersection.rs

@@ -1,4 +1,4 @@
-use super::cross_product;
+use super::cross_product::cross_product_2d;
 use crate::{CoordFloat, CoordNum};
 use geo_types::Coord;
 
@@ -60,15 +60,15 @@ where
 /// M·T=ac
 /// ```
 ///
-/// The determinant of the matrix `M` is the reciprocal of the cross product
-/// of `ab` and `cd`.
+/// Inverting the matrix `M` involves taking the reciprocal of the determinant
+/// (the determinant is same as the of the [cross_product()] of `ab` and `cd`)
 ///
 /// ```text
 /// 1/(ab×cd)
 /// ```
 ///
 /// Therefore if `ab×cd = 0` the determinant is undefined and the matrix cannot
-/// be inverted this means the lines are either
+/// be inverted. The lines are either
 ///   a) parallel or
 ///   b) collinear
 ///
@@ -105,32 +105,41 @@ pub(super) fn line_intersection_with_parameter<T>(
     b: &Coord<T>,
     c: &Coord<T>,
     d: &Coord<T>,
-) -> LineIntersectionWithParameterResult<T>
+) -> Option<LineIntersectionWithParameterResult<T>>
 where
     T: CoordFloat,
 {
     let ab = *b - *a;
     let cd = *d - *c;
     let ac = *c - *a;
 
-    let ab_cross_cd = cross_product(ab, cd);
-
-    if ab_cross_cd == num_traits::zero() {
-        // TODO: We can't tolerate this situation as it will cause a divide by
-        //       zero in the next step. Even values close to zero are a problem,
-        //       but I don't know how to deal with that problem jut yet
-
-        // TODO: this is prevented anyway by the only use of this function.
-        todo!("")
-    }
-
-    let t_ab = cross_product(ac, cd) / ab_cross_cd;
-    let t_cd = -cross_product(ab, ac) / ab_cross_cd;
-    let intersection = *a + ab * t_ab;
-    LineIntersectionWithParameterResult {
-        t_ab,
-        t_cd,
-        intersection,
+    // TODO: I'm still confused about how to use Kernel / RobustKernel;
+    //       the following did not work. I need to read more code
+    //       from the rest of this repo to understand.
+    // if Kernel::orient2d(*a, *b, *d) == Orientation::Collinear {
+    //       note that it is sufficient to check that only one of
+    //       c or d are colinear with ab because of how they are
+    //       related by the original line string.
+    // TODO: The following line
+    //       - Does not use the Kernel
+    //       - uses an arbitrary threshold value which needs more thought
+    let ab_cross_cd = cross_product_2d(ab, cd);
+    if <f64 as num_traits::NumCast>::from(ab_cross_cd)
+        .unwrap()
+        .abs()
+        < num_traits::cast(0.0000001f64).unwrap()
+    {
+        // Segments are parallel or colinear
+        None
+    } else {
+        let t_ab = cross_product_2d(ac, cd) / ab_cross_cd;
+        let t_cd = -cross_product_2d(ab, ac) / ab_cross_cd;
+        let intersection = *a + ab * t_ab;
+        Some(LineIntersectionWithParameterResult {
+            t_ab,
+            t_cd,
+            intersection,
+        })
     }
 }
 
@@ -144,19 +153,49 @@ mod test {
         let b = Coord { x: 2f64, y: 2f64 };
         let c = Coord { x: 0f64, y: 1f64 };
         let d = Coord { x: 1f64, y: 0f64 };
-        let LineIntersectionWithParameterResult {
+        if let Some(LineIntersectionWithParameterResult {
             t_ab,
             t_cd,
             intersection,
-        } = line_intersection_with_parameter(&a, &b, &c, &d);
-        assert_eq!(t_ab, 0.25f64);
-        assert_eq!(t_cd, 0.5f64);
-        assert_eq!(
+        }) = line_intersection_with_parameter(&a, &b, &c, &d)
+        {
+            assert_eq!(t_ab, 0.25f64);
+            assert_eq!(t_cd, 0.5f64);
+            assert_eq!(
+                intersection,
+                Coord {
+                    x: 0.5f64,
+                    y: 0.5f64
+                }
+            );
+        } else {
+            assert!(false)
+        }
+    }
+
+    #[test]
+    fn test_intersection_colinear() {
+        let a = Coord { x: 3f64, y: 4f64 };
+        let b = Coord { x: 6f64, y: 8f64 };
+        let c = Coord { x: 7f64, y: 7f64 };
+        let d = Coord { x: 10f64, y: 9f64 };
+        if let Some(LineIntersectionWithParameterResult {
+            t_ab,
+            t_cd,
             intersection,
-            Coord {
-                x: 0.5f64,
-                y: 0.5f64
-            }
-        );
+        }) = line_intersection_with_parameter(&a, &b, &c, &d)
+        {
+            assert_eq!(t_ab, 0.25f64);
+            assert_eq!(t_cd, 0.5f64);
+            assert_eq!(
+                intersection,
+                Coord {
+                    x: 0.5f64,
+                    y: 0.5f64
+                }
+            );
+        } else {
+            assert!(false)
+        }
     }
 }

geo/src/algorithm/offset/mod.rs

@@ -1,7 +1,7 @@
+
 mod cross_product;
+mod slice_itertools;
 mod line_intersection;
-mod offset_trait;
 
-use cross_product::cross_product;
-use line_intersection::{line_intersection_with_parameter, LineIntersectionWithParameterResult};
+mod offset_trait;
 pub use offset_trait::Offset;

geo/src/algorithm/offset/offset_trait.rs

@@ -1,11 +1,14 @@
-use super::{cross_product, line_intersection_with_parameter, LineIntersectionWithParameterResult};
+use super::line_intersection::{line_intersection_with_parameter, LineIntersectionWithParameterResult};
+use super::slice_itertools::pairwise;
+
 use crate::{
     CoordFloat,
     // Kernel,
     // Orientation,
     Line,
     LineString,
     MultiLineString,
+    Polygon
 };
 use geo_types::Coord;
 
@@ -77,17 +80,6 @@ where
     }
 }
 
-/// Iterate over a slice in overlapping pairs
-///
-/// ```ignore
-/// let items = vec![1, 2, 3, 4, 5];
-/// let actual_result: Vec<(i32, i32)> = pairwise(&items[..]).map(|(a, b)| (*a, *b)).collect();
-/// let expected_result = vec![(1, 2), (2, 3), (3, 4), (4, 5)];
-/// assert_eq!(actual_result, expected_result);
-/// ```
-fn pairwise<T>(iterable: &[T]) -> std::iter::Zip<std::slice::Iter<T>, std::slice::Iter<T>> {
-    iterable.iter().zip(iterable[1..].iter())
-}
 
 impl<T> Offset<T> for LineString<T>
 where
@@ -112,52 +104,32 @@ where
         result.push(first_point);
         result.extend(pairwise(&offset_segments[..]).flat_map(
             |(Line { start: a, end: b }, Line { start: c, end: d })| {
-                let ab = *b - *a;
-                let cd = *d - *c;
-                let ab_cross_cd = cross_product(ab, cd);
-                // TODO: I'm still confused about how to use Kernel / RobustKernel;
-                //       the following did not work. I need to read more code
-                //       from the rest of this repo to understand.
-                // if Kernel::orient2d(*a, *b, *d) == Orientation::Collinear {
-                //       note that it is sufficient to check that only one of
-                //       c or d are colinear with ab because of how they are
-                //       related by the original line string.
-                // TODO: The following line
-                //       - Does not use the Kernel
-                //       - uses an arbitrary threshold value which needs more thought
-                if <f64 as num_traits::NumCast>::from(ab_cross_cd)
-                    .unwrap()
-                    .abs()
-                    < num_traits::cast(0.0000001f64).unwrap()
-                {
-                    vec![*b]
-                } else {
-                    // TODO: if we can inline this function we only need to
-                    //       calculate `ab_cross_cd` once
-                    let LineIntersectionWithParameterResult {
+                match line_intersection_with_parameter(a, b, c, d) {
+                    None => vec![*b], // colinear
+                    Some(LineIntersectionWithParameterResult {
                         t_ab,
                         t_cd,
                         intersection,
-                    } = line_intersection_with_parameter(a, b, c, d);
-
-                    let zero = num_traits::zero::<T>();
-                    let one = num_traits::one::<T>();
-
-                    let tip_ab = zero <= t_ab && t_ab <= one;
-                    let fip_ab = !tip_ab;
-                    let pfip_ab = fip_ab && t_ab > zero;
-
-                    let tip_cd = zero <= t_cd && t_cd <= one;
-                    let fip_cd = !tip_cd;
-
-                    if tip_ab && tip_cd {
-                        // TODO: test for mitre limit
-                        vec![intersection]
-                    } else if fip_ab && fip_cd && pfip_ab {
-                        // TODO: test for mitre limit
-                        vec![intersection]
-                    } else {
-                        vec![*b, *c]
+                    }) => {
+                        let zero = num_traits::zero::<T>();
+                        let one = num_traits::one::<T>();
+
+                        let tip_ab = zero <= t_ab && t_ab <= one;
+                        let fip_ab = !tip_ab;
+                        let pfip_ab = fip_ab && t_ab > zero;
+
+                        let tip_cd = zero <= t_cd && t_cd <= one;
+                        let fip_cd = !tip_cd;
+
+                        if tip_ab && tip_cd {
+                            // TODO: test for mitre limit
+                            vec![intersection]
+                        } else if fip_ab && fip_cd && pfip_ab {
+                            // TODO: test for mitre limit
+                            vec![intersection]
+                        } else {
+                            vec![*b, *c]
+                        }
                     }
                 }
             },
@@ -178,23 +150,25 @@ where
     }
 }
 
+
+// impl<T> Offset<T> for Polygon<T>
+// where
+//     T: CoordFloat,
+// {
+//     fn offset(&self, distance: T) -> Self {
+//         // TODO: not finished yet... need to do interiors
+//         // self.interiors()
+//         // TODO: is the winding order configurable?
+//         self.exterior();
+//         todo!("Not finished")
+//     }
+// }
+
 #[cfg(test)]
 mod test {
 
-    // crate dependencies
     use crate::{line_string, Coord, Line, MultiLineString, Offset};
 
-    // private imports
-    use super::pairwise;
-
-    #[test]
-    fn test_pairwise() {
-        let items = vec![1, 2, 3, 4, 5];
-        let actual_result: Vec<(i32, i32)> = pairwise(&items[..]).map(|(a, b)| (*a, *b)).collect();
-        let expected_result = vec![(1, 2), (2, 3), (3, 4), (4, 5)];
-        assert_eq!(actual_result, expected_result);
-    }
-
     #[test]
     fn test_offset_line() {
         let input = Line::new(Coord { x: 1f64, y: 1f64 }, Coord { x: 1f64, y: 2f64 });