faster upper bound for line-rect intersection

Author: cookiedan42

geo-benches/src/intersection.rs

@@ -2,6 +2,46 @@ use criterion::{Criterion, criterion_group, criterion_main};
 use geo::intersects::Intersects;
 use geo::{MultiPolygon, Triangle};
 
+fn line_rect_intersection(c: &mut Criterion) {
+    use geo::{Line, Rect, coord};
+
+    c.bench_function("Line intersects Rect (end intersect)", |bencher| {
+        let line = Line::new(coord! {x:0., y:0.}, coord! {x:10., y:10.});
+        let rect = Rect::new(coord! {x:0., y:0.}, coord! {x:10., y:10.});
+
+        bencher.iter(|| {
+            assert!(criterion::black_box(&line).intersects(criterion::black_box(&rect)));
+        });
+    });
+
+    c.bench_function("Line intersects Rect cross", |bencher| {
+        let line = Line::new(coord! {x:-1., y:5.}, coord! {x:11., y:16.});
+        let rect = Rect::new(coord! {x:0., y:0.}, coord! {x:10., y:10.});
+
+        bencher.iter(|| {
+            assert!(criterion::black_box(&line).intersects(criterion::black_box(&rect)));
+        });
+    });
+
+    c.bench_function("Line disjoint Rect bbox disjoint", |bencher| {
+        let line = Line::new(coord! {x:0., y:0.}, coord! {x:10., y:10.});
+        let rect = Rect::new(coord! {x:11., y:11.}, coord! {x:12., y:12.});
+
+        bencher.iter(|| {
+            assert!(!criterion::black_box(&line).intersects(criterion::black_box(&rect)));
+        });
+    });
+
+    c.bench_function("Line disjoint Rect bbox overlap", |bencher| {
+        let line = Line::new(coord! {x:0., y:0.}, coord! {x:10., y:10.});
+        let rect = Rect::new(coord! {x:6., y:4.}, coord! {x:10., y:0.});
+
+        bencher.iter(|| {
+            assert!(!criterion::black_box(&line).intersects(criterion::black_box(&rect)));
+        });
+    });
+}
+
 fn multi_polygon_intersection(c: &mut Criterion) {
     let plot_polygons: MultiPolygon = geo_test_fixtures::nl_plots_wgs84();
     let zone_polygons: MultiPolygon = geo_test_fixtures::nl_zones();
@@ -253,8 +293,10 @@ criterion_group! {
 }
 
 criterion_group! { bench_linestring_poly,linestring_polygon_intersection}
+criterion_group! { bench_line_rect_intersection,line_rect_intersection}
 
 criterion_main!(
+    bench_line_rect_intersection,
     bench_linestring_poly,
     bench_multi_polygons,
     bench_point_rect,

geo/src/algorithm/intersects/line_string.rs

@@ -63,7 +63,19 @@ where
     }
 }
 
-intersects_line_string_impl!(area: Rect<T>);
+impl<T> Intersects<Rect<T>> for LineString<T>
+where
+    T: GeoNum,
+{
+    fn intersects(&self, rhs: &Rect<T>) -> bool {
+        if has_disjoint_bboxes(self, rhs) {
+            return false;
+        }
+        // splitting into `Line` intersects `Rect`
+        self.lines_iter().any(|ln| ln.intersects(rhs))
+    }
+}
+
 intersects_line_string_impl!(area: Triangle<T>);
 
 // Blanket implementation from LineString<T>

geo/src/algorithm/intersects/rect.rs

@@ -23,17 +23,53 @@ where
     T: GeoNum,
 {
     fn intersects(&self, rhs: &Line<T>) -> bool {
-        let lb = self.min();
-        let rt = self.max();
-        let lt = Coord::from((lb.x, rt.y));
-        let rb = Coord::from((rt.x, lb.y));
-        // If either rhs.{start,end} lies inside Rect, then true
-        self.intersects(&rhs.start)
-            || self.intersects(&rhs.end)
-            || Line::new(lt, rt).intersects(rhs)
-            || Line::new(rt, rb).intersects(rhs)
-            || Line::new(lb, rb).intersects(rhs)
-            || Line::new(lt, lb).intersects(rhs)
+        // adds overhead to worst case
+        // but short circuits if a terminal intersects the rectangle
+
+        // if self.intersects(&rhs.start)|| self.intersects(&rhs.end) {
+        //         return true;
+        //     }
+
+        if !self.intersects(&rhs.bounding_rect()) {
+            return false;
+        }
+
+        /*
+        o3  o2
+        |  /
+        | /
+        o0--o1
+        If `rhs` line extended to infinity crosses any of these three lines,
+        and bounding boxes intersect,
+        then `rhs` (Line) intersects `self` (Rect).
+        */
+
+        let c0 = self.min();
+        let c1 = coord! {x: self.max().x, y: self.min().y};
+        let o0 = T::Ker::orient2d(rhs.start, rhs.end, c0);
+        let o1 = T::Ker::orient2d(rhs.start, rhs.end, c1);
+        if o0 != o1 {
+            return true;
+        }
+
+        let c2 = self.max();
+        let o2 = T::Ker::orient2d(rhs.start, rhs.end, c2);
+        if o0 != o2 {
+            return true;
+        }
+
+        let c3 = coord! {x: self.min().x, y: self.max().y};
+        let o3 = T::Ker::orient2d(rhs.start, rhs.end, c3);
+        if o0 != o3 {
+            return true;
+        }
+
+        // At this point we know all the orientations are equal and that the bounding boxes overlap.
+        // The only ways there could be an intersection is if
+        // 1. `self` (Rect) has degenerated to a line, and `rhs` (Line) is on that line.
+        // 2. `self` (Rect) has degenerated to a point, and lies on `rhs` (Line).
+        // 3. `rhs` (Line) is degenerated to a point.
+        o0 == Orientation::Collinear
     }
 }
 
@@ -117,10 +153,52 @@ where
     }
 }
 
+#[cfg(test)]
+mod test_line {
+    use super::*;
+    use crate::wkt;
+
+    #[test]
+    fn test_overlap_bbox_no_overlap() {
+        let rect = wkt! {RECT(6 4, 10 0)};
+        let line = wkt! {LINE(0 0, 10 10)};
+
+        assert!(!rect.intersects(&line));
+    }
+
+    #[test]
+    fn test_degen_line() {
+        let rect = wkt! {RECT(0 0, 10 10)};
+        let line = wkt! {LINE(0 0, 0 0)};
+
+        assert!(rect.intersects(&line));
+    }
+
+    #[test]
+    fn test_degen_rect() {
+        let rect_pt = wkt! {RECT(0 0, 0 10)};
+        let rect_line1 = wkt! {RECT(0 0, 0 10)};
+        let rect_line2 = wkt! {RECT(0 0, 10 0)};
+        let line = wkt! {LINE(0 0, 10 0)};
+
+        assert!(rect_pt.intersects(&line));
+        assert!(rect_line1.intersects(&line));
+        assert!(rect_line2.intersects(&line));
+    }
+}
+
 #[cfg(test)]
 mod test_triangle {
     use super::*;
 
+    #[test]
+    fn test_rhs_degenerate_line() {
+        let s: Rect<f64> = Rect::new((0, 0), (0, 0)).convert();
+        let l = Line::new(coord! {x: 0.0, y: 0.0}, coord! { x: 0.0, y: 0.0 });
+
+        assert!(s.intersects(&l));
+    }
+
     #[test]
     fn test_disjoint() {
         let rect: Rect<f64> = Rect::new((0, 0), (10, 10)).convert();