Large geometry intersection check perf speedup

geo/Cargo.toml

@@ -97,3 +97,7 @@ harness = false
 [[bench]]
 name = "rand_line_crossings"
 harness = false
+
+[[bench]]
+name = "intersects"
+harness = false

geo/benches/intersects.rs

@@ -0,0 +1,49 @@
+#[macro_use]
+extern crate criterion;
+extern crate geo;
+
+use criterion::Criterion;
+use geo::{intersects::Intersects, Coordinate, Line, LineString, Polygon};
+
+
+fn criterion_benchmark(c: &mut Criterion) {
+    let big: Vec<[f64;2 ]> = include!("../src/algorithm/test_fixtures/big.rs");
+    let norway: Vec<[f64;2 ]> = include!("../src/algorithm/test_fixtures/norway_main.rs");
+    c.bench_function("Large Polygon-Line intersection", |bencher| {
+        let polygon = {
+            let exterior = LineString::<f64>::from(big.clone());
+            Polygon::new(exterior, vec![])
+        };
+        let line = Line::new(
+            Coordinate { x: 49.92187499999999, y: 70.95969716686398 },
+            Coordinate { x: 55.1953125, y: 62.59334083012024 }
+        );
+
+        bencher.iter(|| {
+            criterion::black_box(
+                criterion::black_box(&polygon).intersects(criterion::black_box(&line)),
+            );
+        });
+    });
+
+    c.bench_function("Polygon-Polygon intersection", |bencher| {
+
+        let polygon_big = {
+            let exterior = LineString::<f64>::from(big.clone());
+            Polygon::new(exterior, vec![])
+        };
+        let polygon_norway = {
+            let exterior = LineString::<f64>::from(norway.clone());
+            Polygon::new(exterior, vec![])
+        };
+
+        bencher.iter(|| {
+            criterion::black_box(
+                criterion::black_box(&polygon_big).intersects(criterion::black_box(&polygon_norway)),
+            );
+        });
+    });
+}
+
+criterion_group!(benches, criterion_benchmark);
+criterion_main!(benches);

geo/src/algorithm/intersects/polygon.rs

@@ -1,10 +1,12 @@
 use super::{has_disjoint_bboxes, Intersects};
+use crate::geometry::*;
 use crate::utils::{coord_pos_relative_to_ring, CoordPos};
 use crate::BoundingRect;
-use crate::{
-    CoordNum, Coordinate, GeoNum, Line, LineString, MultiLineString, MultiPolygon, Point, Polygon,
-    Rect,
-};
+use crate::{CoordNum, GeoNum};
+use rstar::{RTree, RTreeObject};
+// the largest total number of segments geometries can have before the algorithm switches
+// to using an R*-tree for queries
+const MAX_NAIVE_SEGMENTS: usize = 10000;
 
 impl<T> Intersects<Coordinate<T>> for Polygon<T>
 where
@@ -24,6 +26,7 @@ symmetric_intersects_impl!(Polygon<T>, Point<T>);
 impl<T> Intersects<Line<T>> for Polygon<T>
 where
     T: GeoNum,
+    Line<T>: RTreeObject,
 {
     fn intersects(&self, line: &Line<T>) -> bool {
         self.exterior().intersects(line)
@@ -32,13 +35,42 @@ where
             || self.intersects(&line.end)
     }
 }
+
+impl<T> Intersects<LineString<T>> for Polygon<T>
+where
+    T: GeoNum,
+    Line<T>: RTreeObject,
+{
+    fn intersects(&self, linestring: &LineString<T>) -> bool {
+        if (self.exterior().0.len() + self.interiors().iter().map(|ls| ls.0.len()).sum::<usize>())
+            * linestring.0.len()
+            > MAX_NAIVE_SEGMENTS
+        {
+            let lines_a: Vec<_> = self
+                .exterior()
+                .lines()
+                .chain(self.interiors().iter().flat_map(|ls| ls.lines()))
+                .collect();
+            let tree_a = RTree::bulk_load(lines_a);
+            let lines_b: Vec<_> = linestring.lines().collect();
+            let tree_b = RTree::bulk_load(lines_b);
+            let mut candidates = tree_a.intersection_candidates_with_other_tree(&tree_b);
+            candidates.any(|line_pair| {
+                dbg!("{:?}", line_pair);
+                line_pair.0.intersects(line_pair.1)
+            })
+        } else {
+            linestring.lines().any(|line| self.intersects(&line))
+        }
+    }
+}
 symmetric_intersects_impl!(Line<T>, Polygon<T>);
-symmetric_intersects_impl!(Polygon<T>, LineString<T>);
 symmetric_intersects_impl!(Polygon<T>, MultiLineString<T>);
 
 impl<T> Intersects<Rect<T>> for Polygon<T>
 where
     T: GeoNum,
+    Line<T>: RTreeObject,
 {
     fn intersects(&self, rect: &Rect<T>) -> bool {
         self.intersects(&rect.to_polygon())
@@ -49,17 +81,43 @@ symmetric_intersects_impl!(Rect<T>, Polygon<T>);
 impl<T> Intersects<Polygon<T>> for Polygon<T>
 where
     T: GeoNum,
+    Line<T>: RTreeObject,
 {
     fn intersects(&self, polygon: &Polygon<T>) -> bool {
         if has_disjoint_bboxes(self, polygon) {
             return false;
         }
+        // switch to querying trees above some threshold x: polygons' combined segment count is higher than x
+        if (self.exterior().0.len() + self.interiors().iter().map(|ls| ls.0.len()).sum::<usize>())
+            * (polygon.exterior().0.len()
+                + polygon
+                    .interiors()
+                    .iter()
+                    .map(|ls| ls.0.len())
+                    .sum::<usize>())
+            > MAX_NAIVE_SEGMENTS
+        {
+            let lines_a: Vec<_> = self
+                .exterior()
+                .lines()
+                .chain(self.interiors().iter().flat_map(|ls| ls.lines()))
+                .collect();
+            let tree_a = RTree::bulk_load(lines_a);
 
-        // self intersects (or contains) any line in polygon
-        self.intersects(polygon.exterior()) ||
+            let lines_b: Vec<_> = polygon
+                .exterior()
+                .lines()
+                .chain(polygon.interiors().iter().flat_map(|ls| ls.lines()))
+                .collect();
+            let tree_b = RTree::bulk_load(lines_b);
+            let mut candidates = tree_a.intersection_candidates_with_other_tree(&tree_b);
+            candidates.any(|line_pair| line_pair.0.intersects(line_pair.1))
+        } else {
+            self.intersects(polygon.exterior()) ||
             polygon.interiors().iter().any(|inner_line_string| self.intersects(inner_line_string)) ||
             // self is contained inside polygon
             polygon.intersects(self.exterior())
+        }
     }
 }