diff --git a/geo/CHANGES.md b/geo/CHANGES.md index 32cd974d6..e1bd3f1ab 100644 --- a/geo/CHANGES.md +++ b/geo/CHANGES.md @@ -10,6 +10,8 @@ * * Fix `AffineTransform::compose` ordering to be conventional - such that the argument is applied *after* self. * +* Add `PreparedGeometry` to speed up repeated `Relate` operations. + * * Implement Frechet distance using linear algorithm to avoid `fatal runtime error: stack overflow` and improve overall performances. * * Bump `geo` MSRV to 1.74 and update CI diff --git a/geo/Cargo.toml b/geo/Cargo.toml index fee0dc04b..81a1b1bd8 100644 --- a/geo/Cargo.toml +++ b/geo/Cargo.toml @@ -83,6 +83,10 @@ harness = false name = "euclidean_distance" harness = false +[[bench]] +name = "prepared_geometry" +harness = false + [[bench]] name = "rotate" harness = false diff --git a/geo/benches/prepared_geometry.rs b/geo/benches/prepared_geometry.rs new file mode 100644 index 000000000..96253aa6b --- /dev/null +++ b/geo/benches/prepared_geometry.rs @@ -0,0 +1,65 @@ +use criterion::{criterion_group, criterion_main, Criterion}; +use geo::algorithm::Relate; +use geo::PreparedGeometry; +use geo_types::MultiPolygon; + +fn criterion_benchmark(c: &mut Criterion) { + c.bench_function("relate prepared polygons", |bencher| { + let plot_polygons: MultiPolygon = geo_test_fixtures::nl_plots(); + let zone_polygons = geo_test_fixtures::nl_zones(); + + bencher.iter(|| { + let mut intersects = 0; + let mut non_intersects = 0; + + let plot_polygons = plot_polygons + .iter() + .map(PreparedGeometry::from) + .collect::>(); + + let zone_polygons = zone_polygons + .iter() + .map(PreparedGeometry::from) + .collect::>(); + + for a in &plot_polygons { + for b in &zone_polygons { + if criterion::black_box(a.relate(b).is_intersects()) { + intersects += 1; + } else { + non_intersects += 1; + } + } + } + + assert_eq!(intersects, 974); + assert_eq!(non_intersects, 27782); + }); + }); + + c.bench_function("relate unprepared polygons", |bencher| { + let plot_polygons: MultiPolygon = geo_test_fixtures::nl_plots(); + let zone_polygons = geo_test_fixtures::nl_zones(); + + bencher.iter(|| { + let mut intersects = 0; + let mut non_intersects = 0; + + for a in &plot_polygons { + for b in &zone_polygons { + if criterion::black_box(a.relate(b).is_intersects()) { + intersects += 1; + } else { + non_intersects += 1; + } + } + } + + assert_eq!(intersects, 974); + assert_eq!(non_intersects, 27782); + }); + }); +} + +criterion_group!(benches, criterion_benchmark); +criterion_main!(benches); diff --git a/geo/src/algorithm/relate/geomgraph/edge.rs b/geo/src/algorithm/relate/geomgraph/edge.rs index f21561d4c..c4d372ac8 100644 --- a/geo/src/algorithm/relate/geomgraph/edge.rs +++ b/geo/src/algorithm/relate/geomgraph/edge.rs @@ -7,7 +7,7 @@ use std::collections::BTreeSet; /// An `Edge` represents a one dimensional line in a geometry. /// /// This is based on [JTS's `Edge` as of 1.18.1](https://github.com/locationtech/jts/blob/jts-1.18.1/modules/core/src/main/java/org/locationtech/jts/geomgraph/Edge.java) -#[derive(Debug)] +#[derive(Debug, PartialEq, Clone)] pub(crate) struct Edge { /// `coordinates` of the line geometry coords: Vec>, @@ -48,6 +48,13 @@ impl Edge { &mut self.label } + /// When comparing two prepared geometries, we cache each geometry's topology graph. Depending + /// on the order of the operation - `a.relate(b)` vs `b.relate(a)` - we may need to swap the + /// label. + pub fn swap_label_args(&mut self) { + self.label.swap_args() + } + pub fn coords(&self) -> &[Coord] { &self.coords } @@ -55,6 +62,7 @@ impl Edge { pub fn is_isolated(&self) -> bool { self.is_isolated } + pub fn mark_as_unisolated(&mut self) { self.is_isolated = false; } diff --git a/geo/src/algorithm/relate/geomgraph/edge_intersection.rs b/geo/src/algorithm/relate/geomgraph/edge_intersection.rs index b88908711..b8e4c8c04 100644 --- a/geo/src/algorithm/relate/geomgraph/edge_intersection.rs +++ b/geo/src/algorithm/relate/geomgraph/edge_intersection.rs @@ -6,7 +6,7 @@ use crate::{Coord, GeoFloat}; /// the start of the line segment) The intersection point must be precise. /// /// This is based on [JTS's EdgeIntersection as of 1.18.1](https://github.com/locationtech/jts/blob/jts-1.18.1/modules/core/src/main/java/org/locationtech/jts/geomgraph/EdgeIntersection.java) -#[derive(Debug)] +#[derive(Debug, Clone)] pub(crate) struct EdgeIntersection { coord: Coord, segment_index: usize, diff --git a/geo/src/algorithm/relate/geomgraph/geometry_graph.rs b/geo/src/algorithm/relate/geomgraph/geometry_graph.rs index 714b38e66..6746dc8ac 100644 --- a/geo/src/algorithm/relate/geomgraph/geometry_graph.rs +++ b/geo/src/algorithm/relate/geomgraph/geometry_graph.rs @@ -1,6 +1,7 @@ use super::{ index::{ - EdgeSetIntersector, RstarEdgeSetIntersector, SegmentIntersector, SimpleEdgeSetIntersector, + EdgeSetIntersector, RStarEdgeSetIntersector, Segment, SegmentIntersector, + SimpleEdgeSetIntersector, }, CoordNode, CoordPos, Direction, Edge, Label, LineIntersector, PlanarGraph, TopologyPosition, }; @@ -8,17 +9,18 @@ use super::{ use crate::HasDimensions; use crate::{Coord, GeoFloat, GeometryCow, Line, LineString, Point, Polygon}; +use rstar::{RTree, RTreeNum}; use std::cell::RefCell; use std::rc::Rc; -/// The computation of the [`IntersectionMatrix`] relies on the use of a -/// structure called a "topology graph". The topology graph contains [nodes](CoordNode) and -/// [edges](Edge) corresponding to the nodes and line segments of a [`Geometry`]. Each +/// The computation of the [`IntersectionMatrix`](crate::algorithm::relate::IntersectionMatrix) relies on the use of a +/// structure called a "topology graph". The topology graph contains nodes (CoordNode) and +/// edges (Edge) corresponding to the nodes and line segments of a [`Geometry`](crate::Geometry). Each /// node and edge in the graph is labeled with its topological location /// relative to the source geometry. /// /// Note that there is no requirement that points of self-intersection be a -/// vertex. Thus to obtain a correct topology graph, [`Geometry`] must be +/// vertex. Thus, to obtain a correct topology graph, [`Geometry`](crate::Geometry) must be /// self-noded before constructing their graphs. /// /// Two fundamental operations are supported by topology graphs: @@ -26,13 +28,16 @@ use std::rc::Rc; /// - Computing the intersections between the edges and nodes of two different graphs /// /// GeometryGraph is based on [JTS's `GeomGraph` as of 1.18.1](https://github.com/locationtech/jts/blob/jts-1.18.1/modules/core/src/main/java/org/locationtech/jts/geomgraph/GeometryGraph.java) -pub(crate) struct GeometryGraph<'a, F> +#[derive(Clone)] +pub struct GeometryGraph<'a, F> where F: GeoFloat, { arg_index: usize, - parent_geometry: &'a GeometryCow<'a, F>, + parent_geometry: GeometryCow<'a, F>, + tree: Option>>>, use_boundary_determination_rule: bool, + has_computed_self_nodes: bool, planar_graph: PlanarGraph, } @@ -44,44 +49,102 @@ impl GeometryGraph<'_, F> where F: GeoFloat, { - pub fn edges(&self) -> &[Rc>>] { + pub(crate) fn set_tree(&mut self, tree: Rc>>) { + self.tree = Some(tree); + } + + pub(crate) fn get_or_build_tree(&self) -> Rc>> { + self.tree + .clone() + .unwrap_or_else(|| Rc::new(self.build_tree())) + } + + pub(crate) fn build_tree(&self) -> RTree> { + let segments: Vec> = self + .edges() + .iter() + .enumerate() + .flat_map(|(edge_idx, edge)| { + let edge = RefCell::borrow(edge); + let start_of_final_segment: usize = edge.coords().len() - 1; + (0..start_of_final_segment).map(move |segment_idx| { + let p1 = edge.coords()[segment_idx]; + let p2 = edge.coords()[segment_idx + 1]; + Segment::new(edge_idx, segment_idx, p1, p2) + }) + }) + .collect(); + RTree::bulk_load(segments) + } + + pub(crate) fn assert_eq_graph(&self, other: &Self) { + assert_eq!(self.arg_index, other.arg_index); + assert_eq!( + self.use_boundary_determination_rule, + other.use_boundary_determination_rule + ); + assert_eq!(self.parent_geometry, other.parent_geometry); + self.planar_graph.assert_eq_graph(&other.planar_graph); + } + + pub(crate) fn clone_for_arg_index(&self, arg_index: usize) -> Self { + debug_assert!( + self.has_computed_self_nodes, + "should only be called after computing self nodes" + ); + let planar_graph = self + .planar_graph + .clone_for_arg_index(self.arg_index, arg_index); + Self { + arg_index, + parent_geometry: self.parent_geometry.clone(), + tree: self.tree.clone(), + use_boundary_determination_rule: self.use_boundary_determination_rule, + has_computed_self_nodes: true, + planar_graph, + } + } + + pub(crate) fn edges(&self) -> &[Rc>>] { self.planar_graph.edges() } - pub fn insert_edge(&mut self, edge: Edge) { + pub(crate) fn insert_edge(&mut self, edge: Edge) { self.planar_graph.insert_edge(edge) } - pub fn is_boundary_node(&self, coord: Coord) -> bool { + pub(crate) fn is_boundary_node(&self, coord: Coord) -> bool { self.planar_graph.is_boundary_node(self.arg_index, coord) } - pub fn add_node_with_coordinate(&mut self, coord: Coord) -> &mut CoordNode { + pub(crate) fn add_node_with_coordinate(&mut self, coord: Coord) -> &mut CoordNode { self.planar_graph.add_node_with_coordinate(coord) } - pub fn nodes_iter(&self) -> impl Iterator> { + pub(crate) fn nodes_iter(&self) -> impl Iterator> { self.planar_graph.nodes.iter() } } impl<'a, F> GeometryGraph<'a, F> where - F: GeoFloat, + F: GeoFloat + RTreeNum, { - pub fn new(arg_index: usize, parent_geometry: &'a GeometryCow) -> Self { + pub(crate) fn new(arg_index: usize, parent_geometry: GeometryCow<'a, F>) -> Self { let mut graph = GeometryGraph { arg_index, parent_geometry, use_boundary_determination_rule: true, + tree: None, + has_computed_self_nodes: false, planar_graph: PlanarGraph::new(), }; - graph.add_geometry(parent_geometry); + graph.add_geometry(&graph.parent_geometry.clone()); graph } - pub fn geometry(&self) -> &GeometryCow { - self.parent_geometry + pub(crate) fn geometry(&self) -> &GeometryCow { + &self.parent_geometry } /// Determine whether a component (node or edge) that appears multiple times in elements @@ -96,22 +159,11 @@ where } } - fn create_edge_set_intersector() -> Box> { - // PERF: faster algorithms exist. This one was chosen for simplicity of implementation and - // debugging - // Slow, but simple and good for debugging - // Box::new(SimpleEdgeSetIntersector::new()) - - // Should be much faster for sparse intersections, while not much slower than - // SimpleEdgeSetIntersector in the dense case - Box::new(RstarEdgeSetIntersector::new()) - } - fn boundary_nodes(&self) -> impl Iterator> { self.planar_graph.boundary_nodes(self.arg_index) } - pub fn add_geometry(&mut self, geometry: &GeometryCow) { + pub(crate) fn add_geometry(&mut self, geometry: &GeometryCow) { if geometry.is_empty() { return; } @@ -273,13 +325,13 @@ where /// assumed to be valid). /// /// `line_intersector` the [`LineIntersector`] to use to determine intersection - pub fn compute_self_nodes( - &mut self, - line_intersector: Box>, - ) -> SegmentIntersector { - let mut segment_intersector = SegmentIntersector::new(line_intersector, true); + pub(crate) fn compute_self_nodes(&mut self, line_intersector: Box>) { + if self.has_computed_self_nodes { + return; + } + self.has_computed_self_nodes = true; - let mut edge_set_intersector = Self::create_edge_set_intersector(); + let mut segment_intersector = SegmentIntersector::new(line_intersector, true); // optimize intersection search for valid Polygons and LinearRings let is_rings = match self.geometry() { @@ -290,18 +342,16 @@ where }; let check_for_self_intersecting_edges = !is_rings; + let edge_set_intersector = RStarEdgeSetIntersector; edge_set_intersector.compute_intersections_within_set( - self.edges(), + self, check_for_self_intersecting_edges, &mut segment_intersector, ); - self.add_self_intersection_nodes(); - - segment_intersector } - pub fn compute_edge_intersections( + pub(crate) fn compute_edge_intersections( &self, other: &GeometryGraph, line_intersector: Box>, @@ -312,10 +362,10 @@ where other.boundary_nodes().cloned().collect(), ); - let mut edge_set_intersector = Self::create_edge_set_intersector(); + let edge_set_intersector = RStarEdgeSetIntersector; edge_set_intersector.compute_intersections_between_sets( - self.edges(), - other.edges(), + self, + other, &mut segment_intersector, ); diff --git a/geo/src/algorithm/relate/geomgraph/index/edge_set_intersector.rs b/geo/src/algorithm/relate/geomgraph/index/edge_set_intersector.rs index 709eb2bf1..a1f81af66 100644 --- a/geo/src/algorithm/relate/geomgraph/index/edge_set_intersector.rs +++ b/geo/src/algorithm/relate/geomgraph/index/edge_set_intersector.rs @@ -1,4 +1,4 @@ -use super::super::Edge; +use super::super::{Edge, GeometryGraph}; use super::SegmentIntersector; use crate::{Coord, GeoFloat}; @@ -13,18 +13,18 @@ pub(crate) trait EdgeSetIntersector { /// `check_for_self_intersecting_edges`: if false, an edge is not checked for intersections with itself. /// `segment_intersector`: the SegmentIntersector to use fn compute_intersections_within_set( - &mut self, - edges: &[Rc>>], + &self, + graph: &GeometryGraph, check_for_self_intersecting_edges: bool, segment_intersector: &mut SegmentIntersector, ); /// Compute all intersections between two sets of edges, recording those intersections on /// the intersecting edges. - fn compute_intersections_between_sets( - &mut self, - edges0: &[Rc>>], - edges1: &[Rc>>], + fn compute_intersections_between_sets<'a>( + &self, + graph_0: &GeometryGraph<'a, F>, + graph_1: &GeometryGraph<'a, F>, segment_intersector: &mut SegmentIntersector, ); } diff --git a/geo/src/algorithm/relate/geomgraph/index/mod.rs b/geo/src/algorithm/relate/geomgraph/index/mod.rs index ae48b3db9..5695878e2 100644 --- a/geo/src/algorithm/relate/geomgraph/index/mod.rs +++ b/geo/src/algorithm/relate/geomgraph/index/mod.rs @@ -1,9 +1,13 @@ mod edge_set_intersector; +mod prepared_geometry; mod rstar_edge_set_intersector; +mod segment; mod segment_intersector; mod simple_edge_set_intersector; pub(crate) use edge_set_intersector::EdgeSetIntersector; -pub(crate) use rstar_edge_set_intersector::RstarEdgeSetIntersector; +pub use prepared_geometry::PreparedGeometry; +pub(crate) use rstar_edge_set_intersector::RStarEdgeSetIntersector; +pub(crate) use segment::Segment; pub(crate) use segment_intersector::SegmentIntersector; pub(crate) use simple_edge_set_intersector::SimpleEdgeSetIntersector; diff --git a/geo/src/algorithm/relate/geomgraph/index/prepared_geometry.rs b/geo/src/algorithm/relate/geomgraph/index/prepared_geometry.rs new file mode 100644 index 000000000..2fbc199e2 --- /dev/null +++ b/geo/src/algorithm/relate/geomgraph/index/prepared_geometry.rs @@ -0,0 +1,227 @@ +use super::Segment; +use crate::geometry::*; +use crate::relate::geomgraph::{GeometryGraph, RobustLineIntersector}; +use crate::GeometryCow; +use crate::{GeoFloat, Relate}; + +use std::cell::RefCell; +use std::rc::Rc; + +use rstar::{RTree, RTreeNum}; + +/// A `PreparedGeometry` can be more efficient than a plain Geometry when performing +/// multiple topological comparisons against the `PreparedGeometry`. +/// +/// ``` +/// use geo::{Relate, PreparedGeometry, wkt}; +/// +/// let polygon = wkt! { POLYGON((2.0 2.0,2.0 6.0,4.0 6.0)) }; +/// let touching_line = wkt! { LINESTRING(0.0 0.0,2.0 2.0) }; +/// let intersecting_line = wkt! { LINESTRING(0.0 0.0,3.0 3.0) }; +/// let contained_line = wkt! { LINESTRING(2.0 2.0,3.0 5.0) }; +/// +/// let prepared_polygon = PreparedGeometry::from(polygon); +/// assert!(prepared_polygon.relate(&touching_line).is_touches()); +/// assert!(prepared_polygon.relate(&intersecting_line).is_intersects()); +/// assert!(prepared_polygon.relate(&contained_line).is_contains()); +/// +/// ``` +pub struct PreparedGeometry<'a, F: GeoFloat + RTreeNum = f64> { + geometry_graph: GeometryGraph<'a, F>, +} + +mod conversions { + use crate::geometry_cow::GeometryCow; + use crate::relate::geomgraph::{GeometryGraph, RobustLineIntersector}; + use crate::{GeoFloat, PreparedGeometry}; + use geo_types::{ + Geometry, GeometryCollection, Line, LineString, MultiLineString, MultiPoint, MultiPolygon, + Point, Polygon, Rect, Triangle, + }; + use std::rc::Rc; + + impl<'a, F: GeoFloat> From> for PreparedGeometry<'a, F> { + fn from(point: Point) -> Self { + PreparedGeometry::from(GeometryCow::from(point)) + } + } + impl<'a, F: GeoFloat> From> for PreparedGeometry<'a, F> { + fn from(line: Line) -> Self { + PreparedGeometry::from(GeometryCow::from(line)) + } + } + impl<'a, F: GeoFloat> From> for PreparedGeometry<'a, F> { + fn from(line_string: LineString) -> Self { + PreparedGeometry::from(GeometryCow::from(line_string)) + } + } + impl<'a, F: GeoFloat> From> for PreparedGeometry<'a, F> { + fn from(polygon: Polygon) -> Self { + PreparedGeometry::from(GeometryCow::from(polygon)) + } + } + impl<'a, F: GeoFloat> From> for PreparedGeometry<'a, F> { + fn from(multi_point: MultiPoint) -> Self { + PreparedGeometry::from(GeometryCow::from(multi_point)) + } + } + impl<'a, F: GeoFloat> From> for PreparedGeometry<'a, F> { + fn from(multi_line_string: MultiLineString) -> Self { + PreparedGeometry::from(GeometryCow::from(multi_line_string)) + } + } + impl<'a, F: GeoFloat> From> for PreparedGeometry<'a, F> { + fn from(multi_polygon: MultiPolygon) -> Self { + PreparedGeometry::from(GeometryCow::from(multi_polygon)) + } + } + impl<'a, F: GeoFloat> From> for PreparedGeometry<'a, F> { + fn from(rect: Rect) -> Self { + PreparedGeometry::from(GeometryCow::from(rect)) + } + } + impl<'a, F: GeoFloat> From> for PreparedGeometry<'a, F> { + fn from(triangle: Triangle) -> Self { + PreparedGeometry::from(GeometryCow::from(triangle)) + } + } + impl<'a, F: GeoFloat> From> for PreparedGeometry<'a, F> { + fn from(geometry_collection: GeometryCollection) -> Self { + PreparedGeometry::from(GeometryCow::from(geometry_collection)) + } + } + impl<'a, F: GeoFloat> From> for PreparedGeometry<'a, F> { + fn from(geometry: Geometry) -> Self { + PreparedGeometry::from(GeometryCow::from(geometry)) + } + } + + impl<'a, F: GeoFloat> From<&'a Point> for PreparedGeometry<'a, F> { + fn from(point: &'a Point) -> Self { + PreparedGeometry::from(GeometryCow::from(point)) + } + } + impl<'a, F: GeoFloat> From<&'a Line> for PreparedGeometry<'a, F> { + fn from(line: &'a Line) -> Self { + PreparedGeometry::from(GeometryCow::from(line)) + } + } + impl<'a, F: GeoFloat> From<&'a LineString> for PreparedGeometry<'a, F> { + fn from(line_string: &'a LineString) -> Self { + PreparedGeometry::from(GeometryCow::from(line_string)) + } + } + impl<'a, F: GeoFloat> From<&'a Polygon> for PreparedGeometry<'a, F> { + fn from(polygon: &'a Polygon) -> Self { + PreparedGeometry::from(GeometryCow::from(polygon)) + } + } + impl<'a, F: GeoFloat> From<&'a MultiPoint> for PreparedGeometry<'a, F> { + fn from(multi_point: &'a MultiPoint) -> Self { + PreparedGeometry::from(GeometryCow::from(multi_point)) + } + } + impl<'a, F: GeoFloat> From<&'a MultiLineString> for PreparedGeometry<'a, F> { + fn from(multi_line_string: &'a MultiLineString) -> Self { + PreparedGeometry::from(GeometryCow::from(multi_line_string)) + } + } + impl<'a, F: GeoFloat> From<&'a MultiPolygon> for PreparedGeometry<'a, F> { + fn from(multi_polygon: &'a MultiPolygon) -> Self { + PreparedGeometry::from(GeometryCow::from(multi_polygon)) + } + } + impl<'a, F: GeoFloat> From<&'a GeometryCollection> for PreparedGeometry<'a, F> { + fn from(geometry_collection: &'a GeometryCollection) -> Self { + PreparedGeometry::from(GeometryCow::from(geometry_collection)) + } + } + impl<'a, F: GeoFloat> From<&'a Rect> for PreparedGeometry<'a, F> { + fn from(rect: &'a Rect) -> Self { + PreparedGeometry::from(GeometryCow::from(rect)) + } + } + impl<'a, F: GeoFloat> From<&'a Triangle> for PreparedGeometry<'a, F> { + fn from(triangle: &'a Triangle) -> Self { + PreparedGeometry::from(GeometryCow::from(triangle)) + } + } + + impl<'a, F: GeoFloat> From<&'a Geometry> for PreparedGeometry<'a, F> { + fn from(geometry: &'a Geometry) -> Self { + PreparedGeometry::from(GeometryCow::from(geometry)) + } + } + + impl<'a, F: GeoFloat> From> 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())); + + // 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 } + } + } +} + +impl<'a, F> PreparedGeometry<'a, F> +where + F: GeoFloat + RTreeNum, +{ + pub(crate) fn geometry(&self) -> &GeometryCow { + self.geometry_graph.geometry() + } +} + +impl Relate for PreparedGeometry<'_, F> { + /// Efficiently builds a [`GeometryGraph`] which can then be used for topological + /// computations. + fn geometry_graph(&self, arg_index: usize) -> GeometryGraph { + self.geometry_graph.clone_for_arg_index(arg_index) + } +} + +#[cfg(test)] +mod tests { + use super::*; + use crate::algorithm::Relate; + use crate::polygon; + + #[test] + fn relate() { + let p1 = polygon![(x: 0.0, y: 0.0), (x: 2.0, y: 0.0), (x: 1.0, y: 1.0)]; + let p2 = polygon![(x: 0.5, y: 0.0), (x: 2.0, y: 0.0), (x: 1.0, y: 1.0)]; + let prepared_1 = PreparedGeometry::from(&p1); + let prepared_2 = PreparedGeometry::from(&p2); + assert!(prepared_1.relate(&prepared_2).is_contains()); + assert!(prepared_2.relate(&prepared_1).is_within()); + } + + #[test] + fn prepared_with_unprepared() { + let p1 = polygon![(x: 0.0, y: 0.0), (x: 2.0, y: 0.0), (x: 1.0, y: 1.0)]; + let p2 = polygon![(x: 0.5, y: 0.0), (x: 2.0, y: 0.0), (x: 1.0, y: 1.0)]; + let prepared_1 = PreparedGeometry::from(&p1); + assert!(prepared_1.relate(&p2).is_contains()); + assert!(p2.relate(&prepared_1).is_within()); + } + + #[test] + fn swap_arg_index() { + let poly = polygon![(x: 0.0, y: 0.0), (x: 2.0, y: 0.0), (x: 1.0, y: 1.0)]; + let prepared_geom = PreparedGeometry::from(&poly); + + let poly_cow = GeometryCow::from(&poly); + + let cached_graph = prepared_geom.geometry_graph(0); + let fresh_graph = GeometryGraph::new(0, poly_cow.clone()); + cached_graph.assert_eq_graph(&fresh_graph); + + let cached_graph = prepared_geom.geometry_graph(1); + let fresh_graph = GeometryGraph::new(1, poly_cow); + cached_graph.assert_eq_graph(&fresh_graph); + } +} diff --git a/geo/src/algorithm/relate/geomgraph/index/rstar_edge_set_intersector.rs b/geo/src/algorithm/relate/geomgraph/index/rstar_edge_set_intersector.rs index 63c7b2b69..5c0190d33 100644 --- a/geo/src/algorithm/relate/geomgraph/index/rstar_edge_set_intersector.rs +++ b/geo/src/algorithm/relate/geomgraph/index/rstar_edge_set_intersector.rs @@ -1,105 +1,62 @@ -use super::super::Edge; -use super::{EdgeSetIntersector, SegmentIntersector}; +use super::super::{Edge, GeometryGraph}; +use super::{EdgeSetIntersector, Segment, SegmentIntersector}; use crate::GeoFloat; use std::cell::RefCell; use std::rc::Rc; -use rstar::RTree; +use rstar::{RTree, RTreeNum}; -pub(crate) struct RstarEdgeSetIntersector; +pub(crate) struct RStarEdgeSetIntersector; -impl RstarEdgeSetIntersector { - pub fn new() -> Self { - RstarEdgeSetIntersector - } -} - -struct Segment<'a, F: GeoFloat + rstar::RTreeNum> { - i: usize, - edge: &'a RefCell>, - envelope: rstar::AABB>, -} - -impl<'a, F> Segment<'a, F> -where - F: GeoFloat + rstar::RTreeNum, -{ - fn new(i: usize, edge: &'a RefCell>) -> Self { - use rstar::RTreeObject; - let p1 = edge.borrow().coords()[i]; - let p2 = edge.borrow().coords()[i + 1]; - Self { - i, - edge, - envelope: rstar::AABB::from_corners(p1, p2), - } - } -} - -impl<'a, F> rstar::RTreeObject for Segment<'a, F> -where - F: GeoFloat + rstar::RTreeNum, -{ - type Envelope = rstar::AABB>; - - fn envelope(&self) -> Self::Envelope { - self.envelope - } -} - -impl EdgeSetIntersector for RstarEdgeSetIntersector +impl EdgeSetIntersector for RStarEdgeSetIntersector where - F: GeoFloat + rstar::RTreeNum, + F: GeoFloat + RTreeNum, { fn compute_intersections_within_set( - &mut self, - edges: &[Rc>>], + &self, + graph: &GeometryGraph, check_for_self_intersecting_edges: bool, segment_intersector: &mut SegmentIntersector, ) { - let segments: Vec> = edges - .iter() - .flat_map(|edge| { - let start_of_final_segment: usize = RefCell::borrow(edge).coords().len() - 1; - (0..start_of_final_segment).map(|segment_i| Segment::new(segment_i, edge)) - }) - .collect(); - let tree = RTree::bulk_load(segments); - - for (edge0, edge1) in tree.intersection_candidates_with_other_tree(&tree) { - if check_for_self_intersecting_edges || edge0.edge.as_ptr() != edge1.edge.as_ptr() { - segment_intersector.add_intersections(edge0.edge, edge0.i, edge1.edge, edge1.i); + let edges = graph.edges(); + + let tree = graph.get_or_build_tree(); + for (segment_0, segment_1) in tree.intersection_candidates_with_other_tree(&tree) { + if check_for_self_intersecting_edges || segment_0.edge_idx != segment_1.edge_idx { + let edge_0 = &edges[segment_0.edge_idx]; + let edge_1 = &edges[segment_1.edge_idx]; + segment_intersector.add_intersections( + edge_0, + segment_0.segment_idx, + edge_1, + segment_1.segment_idx, + ); } } } - fn compute_intersections_between_sets( - &mut self, - edges0: &[Rc>>], - edges1: &[Rc>>], + fn compute_intersections_between_sets<'a>( + &self, + graph_0: &GeometryGraph<'a, F>, + graph_1: &GeometryGraph<'a, F>, segment_intersector: &mut SegmentIntersector, ) { - let segments0: Vec> = edges0 - .iter() - .flat_map(|edge| { - let start_of_final_segment: usize = RefCell::borrow(edge).coords().len() - 1; - (0..start_of_final_segment).map(|segment_i| Segment::new(segment_i, edge)) - }) - .collect(); - let tree_0 = RTree::bulk_load(segments0); - - let segments1: Vec> = edges1 - .iter() - .flat_map(|edge| { - let start_of_final_segment: usize = RefCell::borrow(edge).coords().len() - 1; - (0..start_of_final_segment).map(|segment_i| Segment::new(segment_i, edge)) - }) - .collect(); - let tree_1 = RTree::bulk_load(segments1); - - for (edge0, edge1) in tree_0.intersection_candidates_with_other_tree(&tree_1) { - segment_intersector.add_intersections(edge0.edge, edge0.i, edge1.edge, edge1.i); + let edges_0 = graph_0.edges(); + let edges_1 = graph_1.edges(); + + let tree_0 = graph_0.get_or_build_tree(); + let tree_1 = graph_1.get_or_build_tree(); + + for (segment_0, segment_1) in tree_0.intersection_candidates_with_other_tree(&tree_1) { + let edge_0 = &edges_0[segment_0.edge_idx]; + let edge_1 = &edges_1[segment_1.edge_idx]; + segment_intersector.add_intersections( + edge_0, + segment_0.segment_idx, + edge_1, + segment_1.segment_idx, + ); } } } diff --git a/geo/src/algorithm/relate/geomgraph/index/segment.rs b/geo/src/algorithm/relate/geomgraph/index/segment.rs new file mode 100644 index 000000000..11c5c2543 --- /dev/null +++ b/geo/src/algorithm/relate/geomgraph/index/segment.rs @@ -0,0 +1,33 @@ +use crate::Coord; +use crate::GeoFloat; + +#[derive(Debug, Clone)] +pub(crate) struct Segment { + pub edge_idx: usize, + pub segment_idx: usize, + pub envelope: rstar::AABB>, +} + +impl Segment +where + F: GeoFloat + rstar::RTreeNum, +{ + pub fn new(edge_idx: usize, segment_idx: usize, p1: Coord, p2: Coord) -> Self { + Self { + edge_idx, + segment_idx, + envelope: rstar::AABB::from_corners(p1, p2), + } + } +} + +impl rstar::RTreeObject for Segment +where + F: GeoFloat + rstar::RTreeNum, +{ + type Envelope = rstar::AABB>; + + fn envelope(&self) -> Self::Envelope { + self.envelope + } +} diff --git a/geo/src/algorithm/relate/geomgraph/index/simple_edge_set_intersector.rs b/geo/src/algorithm/relate/geomgraph/index/simple_edge_set_intersector.rs index db9cc9845..a73593446 100644 --- a/geo/src/algorithm/relate/geomgraph/index/simple_edge_set_intersector.rs +++ b/geo/src/algorithm/relate/geomgraph/index/simple_edge_set_intersector.rs @@ -1,4 +1,4 @@ -use super::super::Edge; +use super::super::{Edge, GeometryGraph}; use super::{EdgeSetIntersector, SegmentIntersector}; use crate::GeoFloat; @@ -13,7 +13,7 @@ impl SimpleEdgeSetIntersector { } fn compute_intersects( - &mut self, + &self, edge0: &Rc>>, edge1: &Rc>>, segment_intersector: &mut SegmentIntersector, @@ -30,11 +30,12 @@ impl SimpleEdgeSetIntersector { impl EdgeSetIntersector for SimpleEdgeSetIntersector { fn compute_intersections_within_set( - &mut self, - edges: &[Rc>>], + &self, + graph: &GeometryGraph, check_for_self_intersecting_edges: bool, segment_intersector: &mut SegmentIntersector, ) { + let edges = graph.edges(); for edge0 in edges.iter() { for edge1 in edges.iter() { if check_for_self_intersecting_edges || edge0.as_ptr() != edge1.as_ptr() { @@ -44,14 +45,17 @@ impl EdgeSetIntersector for SimpleEdgeSetIntersector { } } - fn compute_intersections_between_sets( - &mut self, - edges0: &[Rc>>], - edges1: &[Rc>>], + fn compute_intersections_between_sets<'a>( + &self, + graph_0: &GeometryGraph<'a, F>, + graph_1: &GeometryGraph<'a, F>, segment_intersector: &mut SegmentIntersector, ) { - for edge0 in edges0 { - for edge1 in edges1 { + let edges_0 = graph_0.edges(); + let edges_1 = graph_1.edges(); + + for edge0 in edges_0 { + for edge1 in edges_1 { self.compute_intersects(edge0, edge1, segment_intersector); } } diff --git a/geo/src/algorithm/relate/geomgraph/intersection_matrix.rs b/geo/src/algorithm/relate/geomgraph/intersection_matrix.rs index 1cd9f30ac..c365edf80 100644 --- a/geo/src/algorithm/relate/geomgraph/intersection_matrix.rs +++ b/geo/src/algorithm/relate/geomgraph/intersection_matrix.rs @@ -1,6 +1,7 @@ -use crate::{coordinate_position::CoordPos, dimensions::Dimensions}; +use crate::{coordinate_position::CoordPos, dimensions::Dimensions, GeoNum, GeometryCow}; use crate::geometry_cow::GeometryCow::Point; +use crate::relate::geomgraph::intersection_matrix::dimension_matcher::DimensionMatcher; use std::str::FromStr; /// Models a *Dimensionally Extended Nine-Intersection Model (DE-9IM)* matrix. @@ -104,10 +105,57 @@ impl std::fmt::Debug for IntersectionMatrix { } impl IntersectionMatrix { - pub fn empty() -> Self { + pub const fn empty() -> Self { IntersectionMatrix(LocationArray([LocationArray([Dimensions::Empty; 3]); 3])) } + pub(crate) const fn empty_disjoint() -> Self { + IntersectionMatrix(LocationArray([ + LocationArray([Dimensions::Empty, Dimensions::Empty, Dimensions::Empty]), + LocationArray([Dimensions::Empty, Dimensions::Empty, Dimensions::Empty]), + // since Geometries are finite and embedded in a 2-D space, + // the `(Outside, Outside)` element must always be 2-D + LocationArray([ + Dimensions::Empty, + Dimensions::Empty, + Dimensions::TwoDimensional, + ]), + ])) + } + + /// 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( + &mut self, + geometry_a: &GeometryCow, + geometry_b: &GeometryCow, + ) { + use crate::algorithm::dimensions::HasDimensions; + { + let dimensions = geometry_a.dimensions(); + if dimensions != Dimensions::Empty { + self.set(CoordPos::Inside, CoordPos::Outside, dimensions); + + let boundary_dimensions = geometry_a.boundary_dimensions(); + if boundary_dimensions != Dimensions::Empty { + self.set(CoordPos::OnBoundary, CoordPos::Outside, boundary_dimensions); + } + } + } + + { + let dimensions = geometry_b.dimensions(); + if dimensions != Dimensions::Empty { + self.set(CoordPos::Outside, CoordPos::Inside, dimensions); + + let boundary_dimensions = geometry_b.boundary_dimensions(); + if boundary_dimensions != Dimensions::Empty { + self.set(CoordPos::Outside, CoordPos::OnBoundary, boundary_dimensions); + } + } + } + } + /// Set `dimensions` of the cell specified by the positions. /// /// `position_a`: which position `dimensions` applies to within the first geometry diff --git a/geo/src/algorithm/relate/geomgraph/label.rs b/geo/src/algorithm/relate/geomgraph/label.rs index ee7ff64b7..85da82fed 100644 --- a/geo/src/algorithm/relate/geomgraph/label.rs +++ b/geo/src/algorithm/relate/geomgraph/label.rs @@ -13,7 +13,7 @@ use std::fmt; /// /// If the component has *no* incidence with one of the geometries, than the `Label`'s /// `TopologyPosition` for that geometry is called `empty`. -#[derive(Clone)] +#[derive(Clone, PartialEq)] pub(crate) struct Label { geometry_topologies: [TopologyPosition; 2], } @@ -29,6 +29,10 @@ impl fmt::Debug for Label { } impl Label { + pub fn swap_args(&mut self) { + self.geometry_topologies.swap(0, 1) + } + /// Construct an empty `Label` for relating a 1-D line or 0-D point to both geometries. pub fn empty_line_or_point() -> Label { Label { diff --git a/geo/src/algorithm/relate/geomgraph/mod.rs b/geo/src/algorithm/relate/geomgraph/mod.rs index 6405d6dc0..62af0bee0 100644 --- a/geo/src/algorithm/relate/geomgraph/mod.rs +++ b/geo/src/algorithm/relate/geomgraph/mod.rs @@ -8,7 +8,7 @@ pub(crate) use edge_end::{EdgeEnd, EdgeEndKey}; pub(crate) use edge_end_bundle::{EdgeEndBundle, LabeledEdgeEndBundle}; pub(crate) use edge_end_bundle_star::{EdgeEndBundleStar, LabeledEdgeEndBundleStar}; pub(crate) use edge_intersection::EdgeIntersection; -pub(crate) use geometry_graph::GeometryGraph; +pub use geometry_graph::GeometryGraph; pub(crate) use intersection_matrix::IntersectionMatrix; pub(crate) use label::Label; pub(crate) use line_intersector::{LineIntersection, LineIntersector}; diff --git a/geo/src/algorithm/relate/geomgraph/node.rs b/geo/src/algorithm/relate/geomgraph/node.rs index 550fca210..ed77d524a 100644 --- a/geo/src/algorithm/relate/geomgraph/node.rs +++ b/geo/src/algorithm/relate/geomgraph/node.rs @@ -1,7 +1,7 @@ use super::{CoordPos, Dimensions, EdgeEnd, EdgeEndBundleStar, IntersectionMatrix, Label}; use crate::{Coord, GeoFloat}; -#[derive(Debug, Clone)] +#[derive(Debug, Clone, PartialEq)] pub(crate) struct CoordNode where F: GeoFloat, @@ -11,6 +11,10 @@ where } impl CoordNode { + pub fn swap_label_args(&mut self) { + self.label.swap_args() + } + pub(crate) fn label(&self) -> &Label { &self.label } diff --git a/geo/src/algorithm/relate/geomgraph/node_map.rs b/geo/src/algorithm/relate/geomgraph/node_map.rs index 3ae78cbf9..b77e67565 100644 --- a/geo/src/algorithm/relate/geomgraph/node_map.rs +++ b/geo/src/algorithm/relate/geomgraph/node_map.rs @@ -6,6 +6,7 @@ use std::fmt; use std::marker::PhantomData; /// A map of nodes, indexed by the coordinate of the node +#[derive(Clone, PartialEq)] pub(crate) struct NodeMap where F: GeoFloat, @@ -16,7 +17,7 @@ where } /// Creates the node stored in `NodeMap` -pub(crate) trait NodeFactory { +pub(crate) trait NodeFactory: PartialEq { type Node; fn create_node(coordinate: Coord) -> Self::Node; } diff --git a/geo/src/algorithm/relate/geomgraph/planar_graph.rs b/geo/src/algorithm/relate/geomgraph/planar_graph.rs index 0acf8e405..454e410f1 100644 --- a/geo/src/algorithm/relate/geomgraph/planar_graph.rs +++ b/geo/src/algorithm/relate/geomgraph/planar_graph.rs @@ -7,6 +7,7 @@ use crate::{Coord, GeoFloat}; use std::cell::RefCell; use std::rc::Rc; +#[derive(Clone, PartialEq)] pub(crate) struct PlanarGraphNode; /// The basic node constructor does not allow for incident edges @@ -20,12 +21,44 @@ where } } +#[derive(Clone, PartialEq)] pub(crate) struct PlanarGraph { pub(crate) nodes: NodeMap, edges: Vec>>>, } impl PlanarGraph { + pub fn clone_for_arg_index(&self, from_arg_index: usize, to_arg_index: usize) -> Self { + let mut graph = Self { + nodes: self.nodes.clone(), + // deep copy edges + edges: self + .edges + .iter() + .map(|e| Rc::new(RefCell::new(e.borrow().clone()))) + .collect(), + }; + assert_eq!(from_arg_index, 0); + if from_arg_index != to_arg_index { + graph.swap_labels(); + } + graph + } + + fn swap_labels(&mut self) { + for node in self.nodes.iter_mut() { + node.swap_label_args(); + } + for edge in &mut self.edges { + edge.borrow_mut().swap_label_args(); + } + } + + pub fn assert_eq_graph(&self, other: &Self) { + assert_eq!(self.nodes, other.nodes); + assert_eq!(self.edges, other.edges); + } + pub fn edges(&self) -> &[Rc>>] { &self.edges } diff --git a/geo/src/algorithm/relate/geomgraph/topology_position.rs b/geo/src/algorithm/relate/geomgraph/topology_position.rs index 11cda7f57..17435ff56 100644 --- a/geo/src/algorithm/relate/geomgraph/topology_position.rs +++ b/geo/src/algorithm/relate/geomgraph/topology_position.rs @@ -14,7 +14,7 @@ use std::fmt; /// topological relationship attribute for the [`On`](Direction::On) position. /// /// See [`CoordPos`] for the possible values. -#[derive(Copy, Clone)] +#[derive(Copy, Clone, PartialEq)] pub(crate) enum TopologyPosition { Area { on: Option, diff --git a/geo/src/algorithm/relate/mod.rs b/geo/src/algorithm/relate/mod.rs index f6ba3f851..1959d254d 100644 --- a/geo/src/algorithm/relate/mod.rs +++ b/geo/src/algorithm/relate/mod.rs @@ -1,7 +1,10 @@ pub(crate) use edge_end_builder::EdgeEndBuilder; pub use geomgraph::intersection_matrix::IntersectionMatrix; +use relate_operation::RelateOperation; use crate::geometry::*; +pub use crate::relate::geomgraph::index::PreparedGeometry; +pub use crate::relate::geomgraph::GeometryGraph; use crate::{GeoFloat, GeometryCow}; mod edge_end_builder; @@ -51,66 +54,38 @@ mod relate_operation; /// ``` /// /// Note: `Relate` must not be called on geometries containing `NaN` coordinates. -pub trait Relate { - fn relate(&self, other: &T) -> IntersectionMatrix; -} +pub trait Relate { + /// Construct a [`GeometryGraph`] + fn geometry_graph(&self, arg_index: usize) -> GeometryGraph; -impl Relate> for GeometryCow<'_, F> { - fn relate(&self, other: &GeometryCow) -> IntersectionMatrix { - let mut relate_computer = relate_operation::RelateOperation::new(self, other); - relate_computer.compute_intersection_matrix() + fn relate(&self, other: &impl Relate) -> IntersectionMatrix { + RelateOperation::new(self.geometry_graph(0), other.geometry_graph(1)) + .compute_intersection_matrix() } } macro_rules! relate_impl { - ($k:ty, $t:ty) => { - relate_impl![($k, $t),]; - }; - ($(($k:ty, $t:ty),)*) => { + ($($t:ty ,)*) => { $( - impl Relate for $k { - fn relate(&self, other: &$t) -> IntersectionMatrix { - GeometryCow::from(self).relate(&GeometryCow::from(other)) + impl Relate for $t { + fn geometry_graph(&self, arg_index: usize) -> GeometryGraph { + GeometryGraph::new(arg_index, GeometryCow::from(self)) } } )* }; } -/// Call the given macro with every pair of inputs -/// -/// # Examples -/// -/// ```ignore -/// cartesian_pairs!(foo, [Bar, Baz, Qux]); -/// ``` -/// Is akin to calling: -/// ```ignore -/// foo![(Bar, Bar), (Bar, Baz), (Bar, Qux), (Baz, Bar), (Baz, Baz), (Baz, Qux), (Qux, Bar), (Qux, Baz), (Qux, Qux)]; -/// ``` -macro_rules! cartesian_pairs { - ($macro_name:ident, [$($a:ty),*]) => { - cartesian_pairs_helper! { [] [$($a,)*] [$($a,)*] [$($a,)*] $macro_name} - }; -} - -macro_rules! cartesian_pairs_helper { - // popped all a's - we're done. Use the accumulated output as the input to relate macro. - ([$($out_pairs:tt)*] [] [$($b:ty,)*] $init_b:tt $macro_name:ident) => { - $macro_name!{$($out_pairs)*} - }; - // finished one loop of b, pop next a and reset b - ($out_pairs:tt [$a_car:ty, $($a_cdr:ty,)*] [] $init_b:tt $macro_name:ident) => { - cartesian_pairs_helper!{$out_pairs [$($a_cdr,)*] $init_b $init_b $macro_name} - }; - // pop b through all of b with head of a - ([$($out_pairs:tt)*] [$a_car:ty, $($a_cdr:ty,)*] [$b_car:ty, $($b_cdr:ty,)*] $init_b:tt $macro_name:ident) => { - cartesian_pairs_helper!{[$($out_pairs)* ($a_car, $b_car),] [$a_car, $($a_cdr,)*] [$($b_cdr,)*] $init_b $macro_name} - }; -} - -// Implement Relate for every combination of Geometry. Alternatively we could do something like -// `impl Relate> for Into { }` -// but I don't know that we want to make GeometryCow public (yet?). -cartesian_pairs!(relate_impl, [Point, Line, LineString, Polygon, MultiPoint, MultiLineString, MultiPolygon, Rect, Triangle, GeometryCollection]); -relate_impl!(Geometry, Geometry); +relate_impl![ + Point, + Line, + LineString, + Polygon, + MultiPoint, + MultiLineString, + MultiPolygon, + Rect, + Triangle, + GeometryCollection, + Geometry, +]; diff --git a/geo/src/algorithm/relate/relate_operation.rs b/geo/src/algorithm/relate/relate_operation.rs index 4d8676e41..061bb2919 100644 --- a/geo/src/algorithm/relate/relate_operation.rs +++ b/geo/src/algorithm/relate/relate_operation.rs @@ -32,6 +32,7 @@ where isolated_edges: Vec>>>, } +#[derive(PartialEq)] pub(crate) struct RelateNodeFactory; impl NodeFactory for RelateNodeFactory where @@ -47,13 +48,10 @@ impl<'a, F> RelateOperation<'a, F> where F: GeoFloat, { - pub(crate) fn new( - geom_a: &'a GeometryCow<'a, F>, - geom_b: &'a GeometryCow<'a, F>, - ) -> RelateOperation<'a, F> { + pub(crate) fn new(graph_a: GeometryGraph<'a, F>, graph_b: GeometryGraph<'a, F>) -> Self { Self { - graph_a: GeometryGraph::new(0, geom_a), - graph_b: GeometryGraph::new(1, geom_b), + graph_a, + graph_b, nodes: NodeMap::new(), isolated_edges: vec![], line_intersector: RobustLineIntersector::new(), @@ -61,14 +59,7 @@ where } pub(crate) fn compute_intersection_matrix(&mut self) -> IntersectionMatrix { - let mut intersection_matrix = IntersectionMatrix::empty(); - // since Geometries are finite and embedded in a 2-D space, - // the `(Outside, Outside)` element must always be 2-D - intersection_matrix.set( - CoordPos::Outside, - CoordPos::Outside, - Dimensions::TwoDimensional, - ); + let mut intersection_matrix = IntersectionMatrix::empty_disjoint(); use crate::BoundingRect; use crate::Intersects; @@ -80,7 +71,8 @@ where if bounding_rect_a.intersects(&bounding_rect_b) => {} _ => { // since Geometries don't overlap, we can skip most of the work - self.compute_disjoint_intersection_matrix(&mut intersection_matrix); + intersection_matrix + .compute_disjoint(self.graph_a.geometry(), self.graph_b.geometry()); return intersection_matrix; } } @@ -293,44 +285,6 @@ where } } - /// If the Geometries are disjoint, we need to enter their dimension and boundary dimension in - /// the `Outside` rows in the IM - fn compute_disjoint_intersection_matrix(&self, intersection_matrix: &mut IntersectionMatrix) { - { - let geometry_a = self.graph_a.geometry(); - let dimensions = geometry_a.dimensions(); - if dimensions != Dimensions::Empty { - intersection_matrix.set(CoordPos::Inside, CoordPos::Outside, dimensions); - - let boundary_dimensions = geometry_a.boundary_dimensions(); - if boundary_dimensions != Dimensions::Empty { - intersection_matrix.set( - CoordPos::OnBoundary, - CoordPos::Outside, - boundary_dimensions, - ); - } - } - } - - { - let geometry_b = self.graph_b.geometry(); - let dimensions = geometry_b.dimensions(); - if dimensions != Dimensions::Empty { - intersection_matrix.set(CoordPos::Outside, CoordPos::Inside, dimensions); - - let boundary_dimensions = geometry_b.boundary_dimensions(); - if boundary_dimensions != Dimensions::Empty { - intersection_matrix.set( - CoordPos::Outside, - CoordPos::OnBoundary, - boundary_dimensions, - ); - } - } - } - } - fn update_intersection_matrix( &self, labeled_node_edges: Vec<(CoordNode, LabeledEdgeEndBundleStar)>, @@ -457,12 +411,8 @@ mod test { ] .into(); - let gc1 = GeometryCow::from(&square_a); - let gc2 = GeometryCow::from(&square_b); - let mut relate_computer = RelateOperation::new(&gc1, &gc2); - let intersection_matrix = relate_computer.compute_intersection_matrix(); assert_eq!( - intersection_matrix, + square_a.relate(&square_b), IntersectionMatrix::from_str("FF2FF1212").unwrap() ); } @@ -487,12 +437,8 @@ mod test { ] .into(); - let gca = GeometryCow::from(&square_a); - let gcb = GeometryCow::from(&square_b); - let mut relate_computer = RelateOperation::new(&gca, &gcb); - let intersection_matrix = relate_computer.compute_intersection_matrix(); assert_eq!( - intersection_matrix, + square_a.relate(&square_b), IntersectionMatrix::from_str("212FF1FF2").unwrap() ); } @@ -517,12 +463,8 @@ mod test { ] .into(); - let gca = &GeometryCow::from(&square_a); - let gcb = &GeometryCow::from(&square_b); - let mut relate_computer = RelateOperation::new(gca, gcb); - let intersection_matrix = relate_computer.compute_intersection_matrix(); assert_eq!( - intersection_matrix, + square_a.relate(&square_b), IntersectionMatrix::from_str("212101212").unwrap() ); } diff --git a/geo/src/geometry_cow.rs b/geo/src/geometry_cow.rs index 5cd1b2b51..108c24ab4 100644 --- a/geo/src/geometry_cow.rs +++ b/geo/src/geometry_cow.rs @@ -10,7 +10,7 @@ use std::borrow::Cow; /// This is a way to "upgrade" an inner type to something like a `Geometry` without `moving` it. /// /// As an example, see the [`Relate`] trait which uses `GeometryCow`. -#[derive(PartialEq, Debug, Hash)] +#[derive(PartialEq, Debug, Hash, Clone)] pub(crate) enum GeometryCow<'a, T> where T: CoordNum, @@ -103,3 +103,82 @@ impl<'a, T: CoordNum> From<&'a Triangle> for GeometryCow<'a, T> { GeometryCow::Triangle(Cow::Borrowed(triangle)) } } + +impl<'a, T: CoordNum> From> for GeometryCow<'a, T> { + fn from(point: Point) -> Self { + GeometryCow::Point(Cow::Owned(point)) + } +} + +impl<'a, T: CoordNum> From> for GeometryCow<'a, T> { + fn from(line_string: LineString) -> Self { + GeometryCow::LineString(Cow::Owned(line_string)) + } +} + +impl<'a, T: CoordNum> From> for GeometryCow<'a, T> { + fn from(line: Line) -> Self { + GeometryCow::Line(Cow::Owned(line)) + } +} + +impl<'a, T: CoordNum> From> for GeometryCow<'a, T> { + fn from(polygon: Polygon) -> Self { + GeometryCow::Polygon(Cow::Owned(polygon)) + } +} + +impl<'a, T: CoordNum> From> for GeometryCow<'a, T> { + fn from(multi_point: MultiPoint) -> GeometryCow<'a, T> { + GeometryCow::MultiPoint(Cow::Owned(multi_point)) + } +} + +impl<'a, T: CoordNum> From> for GeometryCow<'a, T> { + fn from(multi_line_string: MultiLineString) -> Self { + GeometryCow::MultiLineString(Cow::Owned(multi_line_string)) + } +} + +impl<'a, T: CoordNum> From> for GeometryCow<'a, T> { + fn from(multi_polygon: MultiPolygon) -> Self { + GeometryCow::MultiPolygon(Cow::Owned(multi_polygon)) + } +} + +impl<'a, T: CoordNum> From> for GeometryCow<'a, T> { + fn from(geometry_collection: GeometryCollection) -> Self { + GeometryCow::GeometryCollection(Cow::Owned(geometry_collection)) + } +} + +impl<'a, T: CoordNum> From> for GeometryCow<'a, T> { + fn from(rect: Rect) -> Self { + GeometryCow::Rect(Cow::Owned(rect)) + } +} + +impl<'a, T: CoordNum> From> for GeometryCow<'a, T> { + fn from(triangle: Triangle) -> Self { + GeometryCow::Triangle(Cow::Owned(triangle)) + } +} + +impl<'a, T: CoordNum> From> for GeometryCow<'a, T> { + fn from(geometry: Geometry) -> Self { + match geometry { + Geometry::Point(point) => GeometryCow::from(point), + Geometry::Line(line) => GeometryCow::from(line), + Geometry::LineString(line_string) => GeometryCow::from(line_string), + Geometry::Polygon(polygon) => GeometryCow::from(polygon), + Geometry::MultiPoint(multi_point) => GeometryCow::from(multi_point), + Geometry::MultiLineString(multi_line_string) => GeometryCow::from(multi_line_string), + Geometry::MultiPolygon(multi_polygon) => GeometryCow::from(multi_polygon), + Geometry::GeometryCollection(geometry_collection) => { + GeometryCow::from(geometry_collection) + } + Geometry::Rect(rect) => GeometryCow::from(rect), + Geometry::Triangle(triangle) => GeometryCow::from(triangle), + } + } +} diff --git a/geo/src/lib.rs b/geo/src/lib.rs index 66968ad15..f6ecfae66 100644 --- a/geo/src/lib.rs +++ b/geo/src/lib.rs @@ -222,6 +222,7 @@ pub use crate::algorithm::*; pub use crate::types::Closest; use std::cmp::Ordering; +pub use crate::relate::PreparedGeometry; pub use geo_types::{coord, line_string, point, polygon, wkt, CoordFloat, CoordNum}; pub mod geometry;