Speed up Relate for PreparedGeometry

geo/CHANGES.md

@@ -79,6 +79,10 @@
   - <https://github.com/georust/geo/pull/1314>
 - Bump `geo` MSRV to 1.80 and update CI
   - <https://github.com/georust/geo/pull/1311>
+- BREAKING: Speed up `Relate` for `PreparedGeometry` - this did require
+  changing some trait constraints, but they are unlikely to affect you in
+  practice unless you have your own Relate implementation.
+  - <https://github.com/georust/geo/pull/1317>
 
 ## 0.29.3 - 2024.12.03
 

geo/src/algorithm/relate/geomgraph/edge_end_bundle_star.rs

@@ -83,7 +83,7 @@ impl<F: GeoFloat> LabeledEdgeEndBundleStar<F> {
         debug!("edge_end_bundle_star: {:?}", self);
     }
 
-    fn propagate_side_labels(&mut self, geom_index: usize, geometry_graph: &GeometryGraph<F>) {
+    fn propagate_side_labels(&mut self, geom_index: usize, _geometry_graph: &GeometryGraph<F>) {
         let mut start_position = None;
 
         for edge_ends in self.edge_end_bundles_iter() {
@@ -108,11 +108,11 @@ impl<F: GeoFloat> LabeledEdgeEndBundleStar<F> {
                 let left_position = label.position(geom_index, Direction::Left);
                 let right_position = label.position(geom_index, Direction::Right);
 
-                if let Some(right_position) = right_position {
+                if let Some(_right_position) = right_position {
                     #[cfg(debug_assertions)]
-                    if right_position != current_position {
+                    if _right_position != current_position {
                         use crate::algorithm::Validation;
-                        if geometry_graph.geometry().is_valid() {
+                        if _geometry_graph.geometry().is_valid() {
                             debug_assert!(false, "topology position conflict with coordinate — this can happen with invalid geometries. coordinate: {:?}, right_location: {:?}, current_location: {:?}", edge_ends.coordinate(), right_position, current_position);
                         } else {
                             warn!("topology position conflict with coordinate — this can happen with invalid geometries. coordinate: {:?}, right_location: {:?}, current_location: {:?}", edge_ends.coordinate(), right_position, current_position);

geo/src/algorithm/relate/geomgraph/index/prepared_geometry.rs

@@ -1,12 +1,13 @@
 use super::Segment;
 use crate::geometry::*;
 use crate::relate::geomgraph::{GeometryGraph, RobustLineIntersector};
-use crate::GeometryCow;
+use crate::{BoundingRect, GeometryCow, HasDimensions};
 use crate::{GeoFloat, Relate};
 
 use std::cell::RefCell;
 use std::rc::Rc;
 
+use crate::dimensions::Dimensions;
 use rstar::{RTree, RTreeNum};
 
 /// A `PreparedGeometry` can be more efficient than a plain Geometry when performing
@@ -27,17 +28,19 @@ use rstar::{RTree, RTreeNum};
 ///
 /// ```
 pub struct PreparedGeometry<'a, F: GeoFloat + RTreeNum = f64> {
-    geometry_graph: GeometryGraph<'a, F>,
+    pub(crate) geometry_graph: GeometryGraph<'a, F>,
+    pub(crate) bounding_rect: Option<Rect<F>>,
 }
 
 mod conversions {
     use crate::geometry_cow::GeometryCow;
     use crate::relate::geomgraph::{GeometryGraph, RobustLineIntersector};
-    use crate::{GeoFloat, PreparedGeometry};
+    use crate::{BoundingRect, GeoFloat, PreparedGeometry};
     use geo_types::{
         Geometry, GeometryCollection, Line, LineString, MultiLineString, MultiPoint, MultiPolygon,
         Point, Polygon, Rect, Triangle,
     };
+    use rstar::Envelope;
     use std::rc::Rc;
 
     impl<F: GeoFloat> From<Point<F>> for PreparedGeometry<'_, F> {
@@ -156,13 +159,29 @@ mod conversions {
     impl<'a, F: GeoFloat> From<GeometryCow<'a, F>> for PreparedGeometry<'a, F> {
         fn from(geometry: GeometryCow<'a, F>) -> Self {
             let mut geometry_graph = GeometryGraph::new(0, geometry);
-            geometry_graph.set_tree(Rc::new(geometry_graph.build_tree()));
+            let r_tree = geometry_graph.build_tree();
+
+            let envelope = r_tree.root().envelope();
+
+            // geo and rstar have different conventions for how to represet an empty bounding boxes
+            let bounding_rect: Option<Rect<F>> = if envelope == rstar::AABB::new_empty() {
+                None
+            } else {
+                Some(Rect::new(envelope.lower(), envelope.upper()))
+            };
+            // They should be equal - but we can save the computation in the `--release` case
+            // by using the bounding_rect from the RTree
+            debug_assert_eq!(bounding_rect, geometry_graph.geometry().bounding_rect());
+
+            geometry_graph.set_tree(Rc::new(r_tree));
 
             // TODO: don't pass in line intersector here - in theory we'll want pluggable line intersectors
             // and the type (Robust) shouldn't be hard coded here.
             geometry_graph.compute_self_nodes(Box::new(RobustLineIntersector::new()));
-
-            Self { geometry_graph }
+            Self {
+                geometry_graph,
+                bounding_rect,
+            }
         }
     }
 }
@@ -176,6 +195,28 @@ where
     }
 }
 
+impl<F: GeoFloat> BoundingRect<F> for PreparedGeometry<'_, F> {
+    type Output = Option<Rect<F>>;
+
+    fn bounding_rect(&self) -> Option<Rect<F>> {
+        self.bounding_rect
+    }
+}
+
+impl<F: GeoFloat> HasDimensions for PreparedGeometry<'_, F> {
+    fn is_empty(&self) -> bool {
+        self.geometry_graph.geometry().is_empty()
+    }
+
+    fn dimensions(&self) -> Dimensions {
+        self.geometry_graph.geometry().dimensions()
+    }
+
+    fn boundary_dimensions(&self) -> Dimensions {
+        self.geometry_graph.geometry().boundary_dimensions()
+    }
+}
+
 impl<F: GeoFloat> Relate<F> for PreparedGeometry<'_, F> {
     /// Efficiently builds a [`GeometryGraph`] which can then be used for topological
     /// computations.

geo/src/algorithm/relate/geomgraph/intersection_matrix.rs

@@ -1,4 +1,6 @@
-use crate::{coordinate_position::CoordPos, dimensions::Dimensions, GeoNum, GeometryCow};
+use crate::{
+    coordinate_position::CoordPos, dimensions::Dimensions, GeoNum, GeometryCow, HasDimensions,
+};
 
 use crate::geometry_cow::GeometryCow::Point;
 use crate::relate::geomgraph::intersection_matrix::dimension_matcher::DimensionMatcher;
@@ -125,12 +127,11 @@ impl IntersectionMatrix {
 
     /// If the Geometries are disjoint, we need to enter their dimension and boundary dimension in
     /// the `Outside` rows in the IM
-    pub(crate) fn compute_disjoint<F: GeoNum>(
+    pub(crate) fn compute_disjoint(
         &mut self,
-        geometry_a: &GeometryCow<F>,
-        geometry_b: &GeometryCow<F>,
+        geometry_a: &(impl HasDimensions + ?Sized),
+        geometry_b: &(impl HasDimensions + ?Sized),
     ) {
-        use crate::algorithm::dimensions::HasDimensions;
         {
             let dimensions = geometry_a.dimensions();
             if dimensions != Dimensions::Empty {

geo/src/algorithm/relate/mod.rs

@@ -5,7 +5,7 @@ use relate_operation::RelateOperation;
 use crate::geometry::*;
 pub use crate::relate::geomgraph::index::PreparedGeometry;
 pub use crate::relate::geomgraph::GeometryGraph;
-use crate::{GeoFloat, GeometryCow};
+use crate::{BoundingRect, GeoFloat, GeometryCow, HasDimensions};
 
 mod edge_end_builder;
 mod geomgraph;
@@ -54,13 +54,20 @@ mod relate_operation;
 /// ```
 ///
 /// Note: `Relate` must not be called on geometries containing `NaN` coordinates.
-pub trait Relate<F: GeoFloat> {
-    /// Construct a [`GeometryGraph`]
-    fn geometry_graph(&self, arg_index: usize) -> GeometryGraph<F>;
+pub trait Relate<F: GeoFloat>: BoundingRect<F> + HasDimensions {
+    /// Returns a noded topology graph for the geometry.
+    ///
+    /// # Params
+    ///
+    /// `idx`: 0 or 1, designating A or B (respectively) in the role this geometry plays
+    ///        in the relation. e.g. in `a.relate(b)`
+    fn geometry_graph(&self, idx: usize) -> GeometryGraph<F>;
 
-    fn relate(&self, other: &impl Relate<F>) -> IntersectionMatrix {
-        RelateOperation::new(self.geometry_graph(0), other.geometry_graph(1))
-            .compute_intersection_matrix()
+    fn relate(&self, other: &impl Relate<F>) -> IntersectionMatrix
+    where
+        Self: Sized,
+    {
+        RelateOperation::new(self, other).compute_intersection_matrix()
     }
 }