diff --git a/geo-benches/Cargo.toml b/geo-benches/Cargo.toml index ded3c54957..2421f49067 100644 --- a/geo-benches/Cargo.toml +++ b/geo-benches/Cargo.toml @@ -169,3 +169,8 @@ harness = false name = "voronoi" path = "src/voronoi.rs" harness = false + +[[bench]] +name = "validation" +path = "src/validation.rs" +harness = false diff --git a/geo-benches/src/validation.rs b/geo-benches/src/validation.rs new file mode 100644 index 0000000000..87ad6b881b --- /dev/null +++ b/geo-benches/src/validation.rs @@ -0,0 +1,66 @@ +//! Benchmarks for polygon validation, focusing on: +//! - Checkerboard patterns (stress test simply-connected interior detection) +//! - Simple and real-world polygon geometries + +use criterion::{BenchmarkId, Criterion, criterion_group, criterion_main}; +use geo::algorithm::Validation; +use geo::coord; +use geo::geometry::{LineString, Polygon}; +use geo_test_fixtures::checkerboard::create_checkerboard_polygon; + +fn create_simple_polygon(size: usize) -> Polygon { + let n = size; + let mut coords = Vec::with_capacity(n + 1); + for i in 0..n { + let angle = 2.0 * std::f64::consts::PI * (i as f64) / (n as f64); + coords.push(coord! { x: angle.cos() * 100.0, y: angle.sin() * 100.0 }); + } + coords.push(coords[0]); + Polygon::new(LineString::new(coords), vec![]) +} + +fn validation_benchmark(c: &mut Criterion) { + let mut group = c.benchmark_group("polygon_validation"); + + for size in [10, 100, 1000] { + let simple = create_simple_polygon(size); + group.bench_with_input( + BenchmarkId::new("simple_polygon", size), + &simple, + |b, poly| b.iter(|| criterion::black_box(poly.is_valid())), + ); + } + + for level in 0..=3 { + let checkerboard = create_checkerboard_polygon(level); + let num_holes = checkerboard.interiors().len(); + let label = format!("level_{level}_{num_holes}holes"); + group.bench_with_input( + BenchmarkId::new("checkerboard", &label), + &checkerboard, + |b, poly| b.iter(|| criterion::black_box(poly.is_valid())), + ); + } + + let east_baton_rouge: Polygon = geo_test_fixtures::east_baton_rouge(); + group.bench_with_input( + BenchmarkId::new("fixture", "east_baton_rouge"), + &east_baton_rouge, + |b, poly| b.iter(|| criterion::black_box(poly.is_valid())), + ); + + let nl_zones = geo_test_fixtures::nl_zones::(); + group.bench_with_input(BenchmarkId::new("fixture", "nl_zones"), &nl_zones, |b, mp| { + b.iter(|| criterion::black_box(mp.is_valid())) + }); + + let nl_plots = geo_test_fixtures::nl_plots_wgs84::(); + group.bench_with_input(BenchmarkId::new("fixture", "nl_plots"), &nl_plots, |b, mp| { + b.iter(|| criterion::black_box(mp.is_valid())) + }); + + group.finish(); +} + +criterion_group!(benches, validation_benchmark); +criterion_main!(benches); diff --git a/geo-test-fixtures/Cargo.toml b/geo-test-fixtures/Cargo.toml index fcbd33280c..3ccaab46fd 100644 --- a/geo-test-fixtures/Cargo.toml +++ b/geo-test-fixtures/Cargo.toml @@ -7,3 +7,4 @@ publish = false [dependencies] wkt = { version = "0.14.0" } geo-types = { path = "../geo-types" } +geo = { path = "../geo" } diff --git a/geo-test-fixtures/src/checkerboard.rs b/geo-test-fixtures/src/checkerboard.rs new file mode 100644 index 0000000000..6d8ca7a4ad --- /dev/null +++ b/geo-test-fixtures/src/checkerboard.rs @@ -0,0 +1,177 @@ +//! Checkerboard polygon generator for testing simply-connected interior validation. +//! +//! The checkerboard pattern is useful for testing cycle detection: each pair of +//! adjacent holes shares exactly one vertex, and these single-touch connections +//! form cycles in the ring adjacency graph, disconnecting the interior. +//! +//! # Pattern Structure +//! +//! ```text +//! 0 1 2 3 4 5 6 7 +//! 7 +-------------------------------+ +//! | | +//! 6 | +---+ +---+ +---+ | +//! | | # | | # | | # | | +//! 5 | +---*---*---*---*---+ | +//! | | # | | # | | +//! 4 | +---*---*---*---*---+ | # = holes +//! | | # | | # | | # | | * = shared vertices +//! 3 | +---*---*---*---*---+ | +//! | | # | | # | | +//! 2 | +---*---*---*---*---+ | +//! | | # | | # | | # | | +//! 1 | +---+ +---+ +---+ | +//! | | +//! 0 +-------------------------------+ +//! ``` + +use geo::Translate; +use geo_types::{coord, LineString, Polygon}; + +/// Create a box as a closed [`LineString`] (for use as exterior or hole). +pub fn box_ring(min_x: f64, min_y: f64, max_x: f64, max_y: f64) -> LineString { + LineString::new(vec![ + coord! { x: min_x, y: min_y }, + coord! { x: max_x, y: min_y }, + coord! { x: max_x, y: max_y }, + coord! { x: min_x, y: max_y }, + coord! { x: min_x, y: min_y }, + ]) +} + +/// Generate checkerboard hole positions for a given level. +/// +/// The pattern creates holes that share vertices at their corners. +/// Level 0 uses unit size 1, level 1 uses unit size 7, etc. +pub fn checkerboard_holes_at_level(level: usize) -> Vec> { + let base_sz: f64 = 7.0; + let sz = base_sz.powi(level as i32); + + let mut holes = Vec::new(); + + // Diagonal holes: (i, i) for i in 1..6 + for i in 1..6 { + let fi = i as f64; + holes.push(box_ring(fi * sz, fi * sz, (fi + 1.0) * sz, (fi + 1.0) * sz)); + } + + // Off-diagonal holes above the diagonal + for i in 1..4 { + let fi = i as f64; + holes.push(box_ring( + fi * sz, + (fi + 2.0) * sz, + (fi + 1.0) * sz, + (fi + 3.0) * sz, + )); + } + for i in 1..2 { + let fi = i as f64; + holes.push(box_ring( + fi * sz, + (fi + 4.0) * sz, + (fi + 1.0) * sz, + (fi + 5.0) * sz, + )); + } + + // Off-diagonal holes below the diagonal + for i in 1..4 { + let fi = i as f64; + holes.push(box_ring( + (fi + 2.0) * sz, + fi * sz, + (fi + 3.0) * sz, + (fi + 1.0) * sz, + )); + } + for i in 1..2 { + let fi = i as f64; + holes.push(box_ring( + (fi + 4.0) * sz, + fi * sz, + (fi + 5.0) * sz, + (fi + 1.0) * sz, + )); + } + + holes +} + +/// Create a checkerboard polygon with the given nesting level. +/// +/// - Level 0: Simple 7x7 checkerboard with 13 holes +/// - Level 1: 49x49 with level-0 checkerboard nested inside one of the "solid" squares +/// - Level 2: 343x343 with level-1 nested inside, etc. +/// +/// The nested checkerboards are placed at offset `(2*sz, 3*sz)` where +/// `sz` is the size of the outer checkerboard's unit cell. This places +/// the nested pattern inside the solid square at grid position (2,3). +/// +/// Returns the polygon and its expected area. +pub fn create_checkerboard(level: usize) -> (Polygon, f64) { + let base_sz: f64 = 7.0; + let sz = base_sz.powi((level + 1) as i32); + + let exterior = box_ring(0.0, 0.0, sz, sz); + let exterior_area = sz * sz; + + let mut all_holes: Vec> = Vec::new(); + let mut total_hole_area = 0.0; + + fn add_holes_recursive( + all_holes: &mut Vec>, + total_hole_area: &mut f64, + current_level: usize, + max_level: usize, + offset_x: f64, + offset_y: f64, + base_sz: f64, + ) { + let unit_sz = base_sz.powi((max_level - current_level) as i32); + + let holes = checkerboard_holes_at_level(max_level - current_level); + let hole_area = unit_sz * unit_sz; + + for hole in holes { + all_holes.push(hole.translate(offset_x, offset_y)); + *total_hole_area += hole_area; + } + + if current_level < max_level { + let next_offset_x = offset_x + 2.0 * unit_sz; + let next_offset_y = offset_y + 3.0 * unit_sz; + add_holes_recursive( + all_holes, + total_hole_area, + current_level + 1, + max_level, + next_offset_x, + next_offset_y, + base_sz, + ); + } + } + + add_holes_recursive( + &mut all_holes, + &mut total_hole_area, + 0, + level, + 0.0, + 0.0, + base_sz, + ); + + let polygon = Polygon::new(exterior, all_holes); + let expected_area = exterior_area - total_hole_area; + + (polygon, expected_area) +} + +/// Create a checkerboard polygon without computing expected area. +/// +/// Convenience wrapper for benchmarking where only the polygon is needed. +pub fn create_checkerboard_polygon(level: usize) -> Polygon { + create_checkerboard(level).0 +} diff --git a/geo-test-fixtures/src/lib.rs b/geo-test-fixtures/src/lib.rs index 558f589c34..8b9deea29e 100644 --- a/geo-test-fixtures/src/lib.rs +++ b/geo-test-fixtures/src/lib.rs @@ -3,6 +3,8 @@ use std::{path::PathBuf, str::FromStr}; use geo_types::{LineString, MultiPoint, MultiPolygon, Point, Polygon}; use wkt::{TryFromWkt, WktFloat}; +pub mod checkerboard; + pub fn louisiana() -> LineString where T: WktFloat + Default + FromStr, diff --git a/geo/CHANGES.md b/geo/CHANGES.md index cbc6e577f9..0c42dfc80b 100644 --- a/geo/CHANGES.md +++ b/geo/CHANGES.md @@ -2,6 +2,12 @@ ## Unreleased +- Add simply connected interior validation for polygons. Polygons with holes that touch at vertices in ways that disconnect the interior (e.g., two holes sharing 2+ vertices, or cycles of holes each sharing a vertex) are now detected as invalid via `Validation::is_valid()`. This aligns with OGC Simple Features and matches PostGIS behavior. + - +- Polygon validation now uses `PreparedGeometry` to cache R-tree structures for interior/exterior containment checks, improving validation speed for polygons with many holes. + - +- Update `i_overlay` to 4.4 and enable OGC-compliant polygon extraction for all boolean operations, fixing cases where holes sharing vertices produced invalid geometry. + - - Fix `CoordinatePosition` for `LineString` to handle dimensionally collapsed input e.g. `LINESTRING(0 0)` is treated like `POINT(0 0)`. - - Fix `CoordinatePosition` for `Triangle` to correctly return `CoordPos::OnBoundary` for coordinate within vertical segment. diff --git a/geo/Cargo.toml b/geo/Cargo.toml index abf932cee9..c6352753c8 100644 --- a/geo/Cargo.toml +++ b/geo/Cargo.toml @@ -37,7 +37,7 @@ proj = { version = "0.31.0", optional = true } robust = "1.1.0" rstar = "0.12.0" serde = { version = "1.0", optional = true, features = ["derive"] } -i_overlay = { version = "4.0.0, < 4.1.0", default-features = false } +i_overlay = { version = "4.4.0, < 4.5.0", default-features = false } sif-itree = "0.4.0" [target.'cfg(not(all(target_family = "wasm", target_os = "unknown")))'.dependencies] diff --git a/geo/src/algorithm/bool_ops/mod.rs b/geo/src/algorithm/bool_ops/mod.rs index 783812d265..e75eb4aa48 100644 --- a/geo/src/algorithm/bool_ops/mod.rs +++ b/geo/src/algorithm/bool_ops/mod.rs @@ -13,7 +13,7 @@ pub use i_overlay::core::fill_rule::FillRule; use i_overlay::core::overlay_rule::OverlayRule; use i_overlay::float::clip::FloatClip; use i_overlay::float::overlay::FloatOverlay; -use i_overlay::float::single::SingleFloatOverlay; +use i_overlay::float::overlay::OverlayOptions; use i_overlay::string::clip::ClipRule; /// Boolean Operations on geometry. @@ -96,7 +96,13 @@ pub trait BooleanOps { ) -> MultiPolygon { let subject = self.rings().map(ring_to_shape_path).collect::>(); let clip = other.rings().map(ring_to_shape_path).collect::>(); - let shapes = subject.overlay(&clip, op.into(), fill_rule); + let shapes = FloatOverlay::with_subj_and_clip_custom( + &subject, + &clip, + OverlayOptions::ogc(), + Default::default(), + ) + .overlay(op.into(), fill_rule); multi_polygon_from_shapes(shapes) } @@ -275,7 +281,9 @@ pub fn unary_union<'a, B: BooleanOps + 'a>( FillRule::Negative }; - let shapes = FloatOverlay::with_subj(&subject).overlay(OverlayRule::Subject, fill_rule); + let shapes = + FloatOverlay::with_subj_custom(&subject, OverlayOptions::ogc(), Default::default()) + .overlay(OverlayRule::Subject, fill_rule); multi_polygon_from_shapes(shapes) } diff --git a/geo/src/algorithm/bool_ops/tests.rs b/geo/src/algorithm/bool_ops/tests.rs index 0a4cd20d0b..852a298c9b 100644 --- a/geo/src/algorithm/bool_ops/tests.rs +++ b/geo/src/algorithm/bool_ops/tests.rs @@ -410,3 +410,76 @@ mod gh_issues { // The goal is just to get here without panic } } + +mod ogc_compliance { + use super::*; + + #[test] + fn test_disconnected_interior_split() { + // iOverlay ocg_tests test_0: Two L-shaped holes share vertices at + // (2,2) and (3,3), disconnecting the interior. OGC mode splits + // this into 2 valid polygons. + let square: Polygon = wkt!(POLYGON((0.0 0.0, 5.0 0.0, 5.0 5.0, 0.0 5.0, 0.0 0.0))); + let l_shapes: MultiPolygon = wkt!(MULTIPOLYGON( + ((1.0 2.0, 1.0 4.0, 3.0 4.0, 3.0 3.0, 2.0 3.0, 2.0 2.0, 1.0 2.0)), + ((2.0 1.0, 2.0 2.0, 3.0 2.0, 3.0 3.0, 4.0 3.0, 4.0 1.0, 2.0 1.0)) + )); + + let result = square.difference(&l_shapes); + + assert_eq!( + result.0.len(), + 2, + "OGC: disconnected interior must be split into 2 polygons" + ); + + let expected: MultiPolygon = wkt!(MULTIPOLYGON( + ( + (0.0 0.0, 5.0 0.0, 5.0 5.0, 0.0 5.0, 0.0 0.0), + (1.0 2.0, 1.0 4.0, 3.0 4.0, 3.0 3.0, 4.0 3.0, 4.0 1.0, 2.0 1.0, 2.0 2.0, 1.0 2.0) + ), + ((2.0 2.0, 3.0 2.0, 3.0 3.0, 2.0 3.0, 2.0 2.0)) + )); + let im = result.relate(&expected); + assert!( + im.is_equal_topo(), + "result: {}, expected: {}", + result.to_wkt(), + expected.to_wkt(), + ); + } + + #[test] + fn test_proper_holes_touching_boundary() { + // iOverlay ocg_tests test_6: U-shaped polygon with two small + // square holes touching the inner boundary. OGC mode produces + // proper interior rings (1 polygon, 3 contours). + let u_shape: Polygon = wkt!(POLYGON(( + 0.0 0.0, 5.0 0.0, 5.0 3.0, 3.0 3.0, + 3.0 2.0, 2.0 2.0, 2.0 3.0, 0.0 3.0, 0.0 0.0 + ))); + let squares: MultiPolygon = wkt!(MULTIPOLYGON( + ((1.0 1.0, 2.0 1.0, 2.0 2.0, 1.0 2.0, 1.0 1.0)), + ((3.0 1.0, 4.0 1.0, 4.0 2.0, 3.0 2.0, 3.0 1.0)) + )); + + let result = u_shape.difference(&squares); + + assert_eq!(result.0.len(), 1, "result should be a single polygon"); + let poly = &result.0[0]; + assert_eq!(poly.interiors().len(), 2, "polygon should have 2 holes"); + + let expected: MultiPolygon = wkt!(MULTIPOLYGON(( + (0.0 0.0, 5.0 0.0, 5.0 3.0, 3.0 3.0, 3.0 2.0, 2.0 2.0, 2.0 3.0, 0.0 3.0, 0.0 0.0), + (1.0 1.0, 1.0 2.0, 2.0 2.0, 2.0 1.0, 1.0 1.0), + (3.0 1.0, 3.0 2.0, 4.0 2.0, 4.0 1.0, 3.0 1.0) + ))); + let im = result.relate(&expected); + assert!( + im.is_equal_topo(), + "result: {}, expected: {}", + result.to_wkt(), + expected.to_wkt(), + ); + } +} diff --git a/geo/src/algorithm/buffer.rs b/geo/src/algorithm/buffer.rs index 47fd17df9d..ca59414fdf 100644 --- a/geo/src/algorithm/buffer.rs +++ b/geo/src/algorithm/buffer.rs @@ -541,7 +541,7 @@ mod tests { coord!(x: 3.0, y: 0.0).into(), coord!(x: -1.0, y: 2.0).into(), ]; - let style = BufferStyle::new(2.0).line_cap(LineCap::Custom(arrow_cap)); + let style = BufferStyle::new(2.0).line_cap(LineCap::Custom(arrow_cap.into())); let actual = line_string.buffer_with_style(style); let expected_output = wkt! { MULTIPOLYGON(((-6.0 0.0,2.0 -4.0,0.0 -2.0,10.0 -2.0,10.3901806473732 -1.9615705609321594,10.765366852283478 -1.8477590680122375,11.111140459775925 -1.6629392206668854,11.414213567972183 -1.4142135679721832,11.662939220666885 -1.1111404597759247,11.847759068012238 -0.7653668522834778,11.96157056093216 -0.39018064737319946,12.0 0.0,12.0 10.0,14.0 8.0,10.0 16.0,6.0 8.0,8.0 10.0,8.0 2.0,0.0 2.0,2.0 4.0,-6.0 0.0))) @@ -708,7 +708,7 @@ mod tests { coord!(x: 3.0, y: 0.0).into(), coord!(x: -1.0, y: 2.0).into(), ]; - let style = BufferStyle::new(2.0).line_cap(LineCap::Custom(arrow_cap)); + let style = BufferStyle::new(2.0).line_cap(LineCap::Custom(arrow_cap.into())); let actual = point.buffer_with_style(style); let expected_output = wkt! { POLYGON EMPTY diff --git a/geo/src/algorithm/validation/geometry_collection.rs b/geo/src/algorithm/validation/geometry_collection.rs index 216a651df7..38c7da79a7 100644 --- a/geo/src/algorithm/validation/geometry_collection.rs +++ b/geo/src/algorithm/validation/geometry_collection.rs @@ -90,11 +90,6 @@ mod tests { assert_eq!( errors[1].to_string(), - "geometry at index 3 is invalid: interior ring at index 0 is not contained within the polygon's exterior" - ); - - assert_eq!( - errors[2].to_string(), "geometry at index 3 is invalid: exterior ring and interior ring at index 0 intersect on a line" ); } diff --git a/geo/src/algorithm/validation/polygon.rs b/geo/src/algorithm/validation/polygon.rs index 026f3ce8fa..88d213baa3 100644 --- a/geo/src/algorithm/validation/polygon.rs +++ b/geo/src/algorithm/validation/polygon.rs @@ -1,8 +1,11 @@ use super::{CoordIndex, RingRole, Validation, utils}; use crate::coordinate_position::CoordPos; use crate::dimensions::Dimensions; -use crate::{GeoFloat, HasDimensions, Polygon, Relate}; +use crate::relate::geomgraph::GeometryGraph; +use crate::{Coord, GeoFloat, HasDimensions, Polygon, PreparedGeometry, Relate}; +use std::cell::RefCell; +use std::collections::HashSet; use std::fmt; /// A [`Polygon`] must follow these rules to be valid: @@ -10,9 +13,7 @@ use std::fmt; /// - [x] boundary rings do not cross /// - [x] boundary rings may touch at points but only as a tangent (i.e. not in a line) /// - [x] interior rings are contained in the exterior ring -/// - [ ] the polygon interior is simply connected (i.e. the rings must not touch in a way that splits the polygon into more than one part) -/// -/// Note: the simple connectivity of the interior is not checked by this implementation. +/// - [x] the polygon interior is simply connected (i.e. the rings must not touch in a way that splits the polygon into more than one part) #[derive(Debug, Clone, PartialEq)] pub enum InvalidPolygon { /// A ring must have at least 4 points to be valid. Note that, in order to close the ring, the first and final points will be identical. @@ -27,6 +28,16 @@ pub enum InvalidPolygon { IntersectingRingsOnALine(RingRole, RingRole), /// A valid polygon's rings must not intersect one another. In this case, the intersection is 2-dimensional. IntersectingRingsOnAnArea(RingRole, RingRole), + /// The polygon interior is not simply connected because rings touch in a way that + /// disconnects the interior into multiple regions. + /// + /// Per OGC Simple Feature Specification (ISO 19125-1), section 6.1.11.1: + /// "The interior of every Surface is a connected point set." + /// + /// This can occur when: + /// - Two rings share two or more vertices (creating a "corridor" that cuts through the interior) + /// - Rings form a cycle of single-vertex touches that encloses part of the interior + InteriorNotSimplyConnected(RingRole, RingRole), } impl fmt::Display for InvalidPolygon { @@ -50,6 +61,12 @@ impl fmt::Display for InvalidPolygon { InvalidPolygon::IntersectingRingsOnAnArea(ring_1, ring_2) => { write!(f, "{ring_1} and {ring_2} intersect on an area") } + InvalidPolygon::InteriorNotSimplyConnected(ring_1, ring_2) => { + write!( + f, + "{ring_1} and {ring_2} touch in a way that disconnects the polygon interior" + ) + } } } } @@ -99,14 +116,32 @@ impl Validation for Polygon { } } + // Skip interior checks if there are no non-empty interiors + let has_interiors = self.interiors().iter().any(|ring| !ring.is_empty()); + if !has_interiors { + return Ok(()); + } + + // Use PreparedGeometry for the exterior to cache its R-tree, avoiding + // graph reconstruction for each interior containment check. let polygon_exterior = Polygon::new(self.exterior().clone(), vec![]); + let prepared_exterior = PreparedGeometry::from(&polygon_exterior); + + // Track ring pairs that already have intersection errors so the + // simply-connected check can skip them (avoids duplicate reporting). + // Keys use GeometryGraph edge indices: 0 = exterior, N = interior N-1. + let mut errored_edge_pairs: HashSet<(usize, usize)> = HashSet::new(); for (interior_1_idx, interior_1) in self.interiors().iter().enumerate() { let ring_role_1 = RingRole::Interior(interior_1_idx); + let edge_idx_1 = interior_1_idx + 1; if interior_1.is_empty() { continue; } - let exterior_vs_interior = polygon_exterior.relate(interior_1); + + let interior_1_as_poly = Polygon::new(interior_1.clone(), vec![]); + let prepared_interior_1 = PreparedGeometry::from(&interior_1_as_poly); + let exterior_vs_interior = prepared_exterior.relate(&prepared_interior_1); if !exterior_vs_interior.is_contains() { handle_validation_error(InvalidPolygon::InteriorRingNotContainedInExteriorRing( @@ -115,29 +150,32 @@ impl Validation for Polygon { } // Interior ring and exterior ring may only touch at point (not as a line) - // and not cross - if exterior_vs_interior.get(CoordPos::OnBoundary, CoordPos::Inside) + if exterior_vs_interior.get(CoordPos::OnBoundary, CoordPos::OnBoundary) == Dimensions::OneDimensional { + errored_edge_pairs.insert((0, edge_idx_1)); handle_validation_error(InvalidPolygon::IntersectingRingsOnALine( RingRole::Exterior, ring_role_1, ))?; } - // PERF: consider using PreparedGeometry - let interior_1_as_poly = Polygon::new(interior_1.clone(), vec![]); - for (interior_2_idx, interior_2) in self.interiors().iter().enumerate().skip(interior_1_idx + 1) { let ring_role_2 = RingRole::Interior(interior_2_idx); + let edge_idx_2 = interior_2_idx + 1; + if interior_2.is_empty() { + continue; + } + let interior_2_as_poly = Polygon::new(interior_2.clone(), vec![]); - let intersection_matrix = interior_1_as_poly.relate(&interior_2_as_poly); + let intersection_matrix = prepared_interior_1.relate(&interior_2_as_poly); if intersection_matrix.get(CoordPos::Inside, CoordPos::Inside) == Dimensions::TwoDimensional { + errored_edge_pairs.insert((edge_idx_1, edge_idx_2)); handle_validation_error(InvalidPolygon::IntersectingRingsOnAnArea( ring_role_1, ring_role_2, @@ -146,6 +184,7 @@ impl Validation for Polygon { if intersection_matrix.get(CoordPos::OnBoundary, CoordPos::OnBoundary) == Dimensions::OneDimensional { + errored_edge_pairs.insert((edge_idx_1, edge_idx_2)); handle_validation_error(InvalidPolygon::IntersectingRingsOnALine( ring_role_1, ring_role_2, @@ -153,10 +192,246 @@ impl Validation for Polygon { } } } + + // Check that the interior is simply connected. + let prepared_polygon = PreparedGeometry::from(self); + if let Some((edge_a, edge_b)) = check_interior_simply_connected_from_graph( + &prepared_polygon.geometry_graph, + &errored_edge_pairs, + ) { + let role_a = edge_index_to_ring_role(edge_a); + let role_b = edge_index_to_ring_role(edge_b); + handle_validation_error(InvalidPolygon::InteriorNotSimplyConnected(role_a, role_b))?; + } + Ok(()) } } +/// Convert a GeometryGraph edge index to a [`RingRole`]. +/// +/// Edge 0 is the exterior ring; edge N (N >= 1) is interior ring N-1. +fn edge_index_to_ring_role(edge_idx: usize) -> RingRole { + if edge_idx == 0 { + RingRole::Exterior + } else { + RingRole::Interior(edge_idx - 1) + } +} + +/// Check that the polygon interior is simply connected using the GeometryGraph. +/// +/// Extracts touch-point information from the pre-computed GeometryGraph. The +/// interior is disconnected when rings touch in a way that creates separate +/// regions. This occurs when: +/// - Two rings share 2+ touch points at different coordinates +/// - Rings form a cycle through distinct single touch points +/// (graph-theoretically: the ring adjacency graph has a cycle whose edges +/// correspond to distinct coordinates, meaning the cycle encloses area) +/// +/// Multiple rings meeting at a *single* coordinate does NOT disconnect +/// the interior — the connected regions can still reach each other around +/// the shared point. +/// +/// Ring pairs in `skip_pairs` are excluded (they already have intersection +/// errors reported by the caller). +/// +/// Returns `None` if the interior is simply connected, or `Some((edge_a, edge_b))` +/// identifying the ring pair that causes disconnection. +fn check_interior_simply_connected_from_graph( + graph: &GeometryGraph, + skip_pairs: &HashSet<(usize, usize)>, +) -> Option<(usize, usize)> { + let edges = graph.edges(); + if edges.len() < 2 { + return None; + } + + // Collect all intersection points with their edge index and vertex status. + struct IntersectionInfo { + coord: Coord, + edge_idx: usize, + is_vertex: bool, + } + + let mut all_intersections: Vec> = Vec::new(); + for (edge_idx, edge) in edges.iter().enumerate() { + let edge = RefCell::borrow(edge); + let coords = edge.coords(); + for ei in edge.edge_intersections() { + let coord = ei.coordinate(); + let is_vertex = coords.contains(&coord); + all_intersections.push(IntersectionInfo { + coord, + edge_idx, + is_vertex, + }); + } + } + + // Sort by (coordinate, edge_idx) using total_cmp for consistent ordering. + all_intersections.sort_by(|a, b| { + a.coord + .x + .total_cmp(&b.coord.x) + .then_with(|| a.coord.y.total_cmp(&b.coord.y)) + .then_with(|| a.edge_idx.cmp(&b.edge_idx)) + }); + + // Group all intersections sharing the same coordinate and build: + // - coord_edges: (coord, edge_a, edge_b) for each touch, used by cycle detection + // - ring_pair_seen: if a pair is seen a second time, it shares 2+ touch + // points at different coordinates, which always disconnects the interior + let mut coord_edges: Vec<(Coord, usize, usize)> = Vec::new(); + let mut ring_pair_seen: HashSet<(usize, usize)> = HashSet::new(); + let mut unique_edges: Vec<(usize, bool)> = Vec::new(); + let mut i = 0; + while i < all_intersections.len() { + let current_coord = all_intersections[i].coord; + + // Find the range of intersections sharing this coordinate. + let mut j = i + 1; + while j < all_intersections.len() && all_intersections[j].coord == current_coord { + j += 1; + } + let group = &all_intersections[i..j]; + + // Deduplicate edges within the group. The same edge can appear multiple + // times with different vertex flags, so we merge them to ensure we only + // count distinct edges and know if any of them is a vertex. + unique_edges.clear(); + for info in group { + if let Some(last) = unique_edges.last_mut() { + if last.0 == info.edge_idx { + last.1 |= info.is_vertex; + continue; + } + } + unique_edges.push((info.edge_idx, info.is_vertex)); + } + + // A "touch" requires 2+ distinct edges, at least one having the point + // as a vertex. Skip pairs the caller already reported as intersecting. + if unique_edges.len() >= 2 && unique_edges.iter().any(|(_, is_v)| *is_v) { + for (idx_a, &(edge_a, _)) in unique_edges.iter().enumerate() { + for &(edge_b, _) in unique_edges.iter().skip(idx_a + 1) { + let key = (edge_a, edge_b); + if !skip_pairs.contains(&key) { + if !ring_pair_seen.insert(key) { + return Some(key); + } + coord_edges.push((current_coord, edge_a, edge_b)); + } + } + } + } + + i = j; + } + + check_ring_touches_disconnect_interior(&mut coord_edges, edges.len()) +} + +/// Detect cycles through distinct coordinates in the ring touch graph. +/// +/// Each entry in `coord_edges` is `(coord, ring_a, ring_b)` representing a +/// single-touch between two rings at the given coordinate. The caller has +/// already handled the case where any pair shares 2+ touch points. +/// +/// A cycle of ring touches only disconnects the interior if the cycle passes +/// through at least two distinct coordinates (enclosing area). A cycle where +/// every edge shares the same coordinate — e.g. three holes all meeting at +/// one point — does not enclose any area and is harmless. +/// +/// To distinguish these cases, edges are sorted and grouped by coordinate, +/// then processed with two Union-Find structures: +/// +/// - **Global UF**: accumulates connectivity across all coordinate groups +/// processed so far. "Are these two rings connected by ANY path of touches?" +/// - **Local UF**: reset for each coordinate group. "Are these two rings +/// connected through touches at THIS coordinate only?" +/// +/// For each edge `(u, v)` in the current coordinate group: +/// - `!global.connected(u, v)`: first connection between these components — no +/// cycle at all. **Safe.** +/// - `global.connected(u, v) && local.connected(u, v)`: the rings are already +/// connected within this coordinate group — a same-coordinate cycle. **Safe.** +/// - `global.connected(u, v) && !local.connected(u, v)`: the rings are +/// reachable through a previously processed (different) coordinate, but not +/// through this one — closing a multi-coordinate cycle. **Invalid.** +fn check_ring_touches_disconnect_interior( + coord_edges: &mut [(Coord, usize, usize)], + num_rings: usize, +) -> Option<(usize, usize)> { + coord_edges.sort_by(|a, b| { + a.0.x + .total_cmp(&b.0.x) + .then_with(|| a.0.y.total_cmp(&b.0.y)) + }); + + let mut global = UnionFind::new(num_rings); + let mut group_start = 0; + while group_start < coord_edges.len() { + let mut group_end = group_start + 1; + while group_end < coord_edges.len() + && coord_edges[group_end].0 == coord_edges[group_start].0 + { + group_end += 1; + } + + let mut local = UnionFind::new(num_rings); + for &(_, u, v) in &coord_edges[group_start..group_end] { + if global.find(u) == global.find(v) && local.find(u) != local.find(v) { + return Some((u.min(v), u.max(v))); + } + global.union(u, v); + local.union(u, v); + } + + group_start = group_end; + } + + None +} + +struct UnionFind { + parent: Vec, + rank: Vec, +} + +impl UnionFind { + fn new(n: usize) -> Self { + Self { + parent: (0..n).collect(), + rank: vec![0; n], + } + } + + fn find(&mut self, mut x: usize) -> usize { + while self.parent[x] != x { + self.parent[x] = self.parent[self.parent[x]]; + x = self.parent[x]; + } + x + } + + fn union(&mut self, x: usize, y: usize) { + let rx = self.find(x); + let ry = self.find(y); + if rx == ry { + return; + } + if self.rank[rx] < self.rank[ry] { + self.parent[rx] = ry; + } else if self.rank[rx] > self.rank[ry] { + self.parent[ry] = rx; + } else { + self.parent[ry] = rx; + self.rank[rx] += 1; + } + } +} + #[cfg(test)] mod tests { use super::*; diff --git a/geo/src/algorithm/validation/tests.rs b/geo/src/algorithm/validation/tests.rs index 626adcc27d..48928d111b 100644 --- a/geo/src/algorithm/validation/tests.rs +++ b/geo/src/algorithm/validation/tests.rs @@ -2,3 +2,300 @@ fn jts_validation_tests() { jts_test_runner::assert_jts_tests_succeed("*Valid*"); } + +mod simply_connected_interior { + //! Tests for simply connected interior validation. + //! + //! OGC Simple Feature Specification (ISO 19125-1), section 6.1.11.1 states: + //! "The interior of every Surface is a connected point set." + //! + //! These tests verify that we correctly detect polygons with disconnected + //! interiors, which can occur when: + //! - Two holes share 2+ vertices (creating a "corridor") + //! - Rings form a cycle of single-vertex touches that encloses part of the interior + + use crate::algorithm::validation::polygon::InvalidPolygon; + use crate::algorithm::validation::{RingRole, Validation}; + use crate::coord; + use crate::geometry::{LineString, Polygon}; + use geo_test_fixtures::checkerboard::{box_ring, create_checkerboard}; + + /// Two L-shaped holes sharing vertices at (2,2) and (3,3). + /// + /// Simplest case: two holes share exactly 2 vertices, creating a + /// "corridor" that cuts through the interior. + #[test] + fn two_holes_sharing_two_vertices() { + let exterior = box_ring(0.0, 0.0, 5.0, 5.0); + + let top_l = LineString::new(vec![ + coord! { x: 1.0, y: 2.0 }, + coord! { x: 2.0, y: 2.0 }, + coord! { x: 2.0, y: 3.0 }, + coord! { x: 3.0, y: 3.0 }, + coord! { x: 3.0, y: 4.0 }, + coord! { x: 1.0, y: 4.0 }, + coord! { x: 1.0, y: 2.0 }, + ]); + + let bottom_l = LineString::new(vec![ + coord! { x: 2.0, y: 1.0 }, + coord! { x: 4.0, y: 1.0 }, + coord! { x: 4.0, y: 3.0 }, + coord! { x: 3.0, y: 3.0 }, + coord! { x: 3.0, y: 2.0 }, + coord! { x: 2.0, y: 2.0 }, + coord! { x: 2.0, y: 1.0 }, + ]); + + let polygon = Polygon::new(exterior, vec![top_l, bottom_l]); + + let errors = polygon.validation_errors(); + assert_eq!(errors.len(), 1); + assert_eq!( + errors[0], + InvalidPolygon::InteriorNotSimplyConnected( + RingRole::Interior(0), + RingRole::Interior(1), + ), + ); + } + + /// Checkerboard level 0: 13 holes in a checkerboard pattern where + /// adjacent holes share single vertices at grid intersections. + #[test] + fn checkerboard_level_0() { + let (polygon, expected_area) = create_checkerboard(0); + + assert_eq!(polygon.interiors().len(), 13); + + use crate::algorithm::Area; + let actual_area = polygon.unsigned_area(); + assert!( + (actual_area - expected_area).abs() < 1e-10, + "Area mismatch: expected {expected_area}, got {actual_area}", + ); + + let errors = polygon.validation_errors(); + assert_eq!(errors.len(), 1); + assert!( + matches!(errors[0], InvalidPolygon::InteriorNotSimplyConnected(_, _)), + "Expected InteriorNotSimplyConnected, got: {:?}", + errors[0], + ); + } + + /// Checkerboard level 1: nested checkerboard with level-0 inside one + /// of the "filled" squares. + #[test] + fn checkerboard_level_1() { + let (polygon, expected_area) = create_checkerboard(1); + + assert_eq!(polygon.interiors().len(), 26); + + use crate::algorithm::Area; + let actual_area = polygon.unsigned_area(); + assert!( + (actual_area - expected_area).abs() < 1e-10, + "Area mismatch: expected {expected_area}, got {actual_area}", + ); + + let errors = polygon.validation_errors(); + assert_eq!(errors.len(), 1); + assert!( + matches!(errors[0], InvalidPolygon::InteriorNotSimplyConnected(_, _)), + "Expected InteriorNotSimplyConnected, got: {:?}", + errors[0], + ); + } + + /// Two holes sharing exactly ONE vertex is valid (interior remains connected). + #[test] + fn holes_sharing_one_vertex_is_valid() { + let exterior = box_ring(0.0, 0.0, 6.0, 6.0); + let hole1 = box_ring(1.0, 1.0, 3.0, 3.0); + let hole2 = box_ring(3.0, 3.0, 5.0, 5.0); + + let polygon = Polygon::new(exterior, vec![hole1, hole2]); + assert!( + polygon.is_valid(), + "Polygon with two holes sharing 1 vertex should be valid", + ); + } + + /// Non-touching holes: always valid. + #[test] + fn separate_holes_is_valid() { + let exterior = box_ring(0.0, 0.0, 10.0, 10.0); + let hole1 = box_ring(1.0, 1.0, 3.0, 3.0); + let hole2 = box_ring(5.0, 5.0, 7.0, 7.0); + let hole3 = box_ring(1.0, 6.0, 3.0, 8.0); + + let polygon = Polygon::new(exterior, vec![hole1, hole2, hole3]); + assert!( + polygon.is_valid(), + "Polygon with separate non-touching holes should be valid", + ); + } + + /// Three triangular holes meeting at a single vertex: valid. + /// + /// When three holes share a single vertex but don't share any edges, + /// the interior is still connected around the perimeter. + #[test] + fn three_holes_meeting_at_vertex_valid() { + let exterior = box_ring(-1.1, -1.1, 1.1, 1.1); + + let hole_a = LineString::new(vec![ + coord! { x: 0.0, y: 0.0 }, + coord! { x: 0.5, y: 0.8660254037844386 }, + coord! { x: 1.0, y: 0.0 }, + coord! { x: 0.0, y: 0.0 }, + ]); + let hole_b = LineString::new(vec![ + coord! { x: 0.0, y: 0.0 }, + coord! { x: -1.0, y: 0.0 }, + coord! { x: -0.5, y: 0.8660254037844386 }, + coord! { x: 0.0, y: 0.0 }, + ]); + let hole_c = LineString::new(vec![ + coord! { x: 0.0, y: 0.0 }, + coord! { x: 0.5, y: -0.8660254037844386 }, + coord! { x: -0.5, y: -0.8660254037844386 }, + coord! { x: 0.0, y: 0.0 }, + ]); + + let polygon = Polygon::new(exterior, vec![hole_a, hole_b, hole_c]); + assert!( + polygon.is_valid(), + "Polygon with 3 holes meeting at one vertex (wedge pattern) should be valid", + ); + } + + /// Four triangular holes forming a cycle through single touches. + /// + /// Each pair shares exactly one vertex at a different coordinate, + /// forming cycle A-B-C-D-A. This encloses the center region. + #[test] + fn four_holes_forming_cycle() { + let exterior = box_ring(-10.0, -10.0, 10.0, 10.0); + + // A: bottom, shares (-8,-8) with B, (8,-8) with D + let hole_a = LineString::new(vec![ + coord! { x: -8.0, y: -8.0 }, + coord! { x: 0.0, y: -4.0 }, + coord! { x: 8.0, y: -8.0 }, + coord! { x: -8.0, y: -8.0 }, + ]); + // B: left, shares (-8,8) with C, (-8,-8) with A + let hole_b = LineString::new(vec![ + coord! { x: -8.0, y: 8.0 }, + coord! { x: -4.0, y: 0.0 }, + coord! { x: -8.0, y: -8.0 }, + coord! { x: -8.0, y: 8.0 }, + ]); + // C: top, shares (8,8) with D, (-8,8) with B + let hole_c = LineString::new(vec![ + coord! { x: 8.0, y: 8.0 }, + coord! { x: 0.0, y: 4.0 }, + coord! { x: -8.0, y: 8.0 }, + coord! { x: 8.0, y: 8.0 }, + ]); + // D: right, shares (8,-8) with A, (8,8) with C + let hole_d = LineString::new(vec![ + coord! { x: 8.0, y: -8.0 }, + coord! { x: 4.0, y: 0.0 }, + coord! { x: 8.0, y: 8.0 }, + coord! { x: 8.0, y: -8.0 }, + ]); + + let polygon = Polygon::new(exterior, vec![hole_a, hole_b, hole_c, hole_d]); + + let errors = polygon.validation_errors(); + assert_eq!(errors.len(), 1); + assert!( + matches!(errors[0], InvalidPolygon::InteriorNotSimplyConnected(_, _)), + "Expected InteriorNotSimplyConnected, got: {:?}", + errors[0], + ); + } + + /// Three triangular holes forming a 3-node cycle through distinct touch points. + /// + /// Each pair shares exactly one vertex at a different coordinate: + /// - H0-H1 share (10, 2) + /// - H1-H2 share (16, 10) + /// - H2-H0 share (4, 10) + /// + /// This forms a cycle H0-H1-H2-H0 through three distinct coordinates, + /// enclosing the central triangle and disconnecting the interior. + #[test] + fn three_holes_forming_triangle_cycle() { + let exterior = box_ring(0.0, 0.0, 20.0, 20.0); + + // H0: bottom-left, shares (4,10) with H2 and (10,2) with H1 + let hole_0 = LineString::new(vec![ + coord! { x: 4.0, y: 10.0 }, + coord! { x: 2.0, y: 2.0 }, + coord! { x: 10.0, y: 2.0 }, + coord! { x: 4.0, y: 10.0 }, + ]); + // H1: bottom-right, shares (10,2) with H0 and (16,10) with H2 + let hole_1 = LineString::new(vec![ + coord! { x: 10.0, y: 2.0 }, + coord! { x: 18.0, y: 2.0 }, + coord! { x: 16.0, y: 10.0 }, + coord! { x: 10.0, y: 2.0 }, + ]); + // H2: top, shares (16,10) with H1 and (4,10) with H0 + let hole_2 = LineString::new(vec![ + coord! { x: 16.0, y: 10.0 }, + coord! { x: 10.0, y: 18.0 }, + coord! { x: 4.0, y: 10.0 }, + coord! { x: 16.0, y: 10.0 }, + ]); + + let polygon = Polygon::new(exterior, vec![hole_0, hole_1, hole_2]); + + let errors = polygon.validation_errors(); + assert_eq!(errors.len(), 1); + assert!( + matches!(errors[0], InvalidPolygon::InteriorNotSimplyConnected(_, _)), + "Expected InteriorNotSimplyConnected, got: {:?}", + errors[0], + ); + } + + /// Two holes forming a chain with vertex-on-edge touches to exterior. + /// + /// H0 touches exterior at (0,5), H1 touches exterior at (20,5), + /// both share (10,5). This creates a chain through the exterior with + /// distinct coordinates, disconnecting the interior. + #[test] + fn hole_chain_with_vertex_on_edge_touch() { + let exterior = box_ring(0.0, 0.0, 20.0, 15.0); + + let hole_0 = LineString::new(vec![ + coord! { x: 0.0, y: 5.0 }, + coord! { x: 10.0, y: 5.0 }, + coord! { x: 5.0, y: 10.0 }, + coord! { x: 0.0, y: 5.0 }, + ]); + let hole_1 = LineString::new(vec![ + coord! { x: 10.0, y: 5.0 }, + coord! { x: 20.0, y: 5.0 }, + coord! { x: 15.0, y: 10.0 }, + coord! { x: 10.0, y: 5.0 }, + ]); + + let polygon = Polygon::new(exterior, vec![hole_0, hole_1]); + + let errors = polygon.validation_errors(); + assert_eq!(errors.len(), 1); + assert!( + matches!(errors[0], InvalidPolygon::InteriorNotSimplyConnected(_, _)), + "Expected InteriorNotSimplyConnected, got: {:?}", + errors[0], + ); + } +}