diff --git a/geo/src/algorithm/mod.rs b/geo/src/algorithm/mod.rs index 791ad98fe..4aca03233 100644 --- a/geo/src/algorithm/mod.rs +++ b/geo/src/algorithm/mod.rs @@ -258,6 +258,10 @@ pub mod triangulate_spade; #[cfg(feature = "spade")] pub use triangulate_spade::TriangulateSpade; +/// Operations on invalid polygons to make them valid +mod validity; +pub use validity::Validify; + /// Vector Operations for 2D coordinates mod vector_ops; pub use vector_ops::Vector2DOps; diff --git a/geo/src/algorithm/validity/impls.rs b/geo/src/algorithm/validity/impls.rs new file mode 100644 index 000000000..8277c43cc --- /dev/null +++ b/geo/src/algorithm/validity/impls.rs @@ -0,0 +1,32 @@ +use super::*; + +use geo_types::{MultiPolygon, Polygon}; + +use crate::GeoFloat; + +/// Splits self intersecting polygons into individual valid polygons +/// +/// Self intersecting here includes banana polygons +/// +/// There are no real guarantees in terms of minimality we can give since the algorithms is a +/// heuristic. In practice it turned out to work quiet well though. +/// +/// # Validity +/// +/// The algo will return non-sense if the input is non-sense. This mainly means: +/// +/// - a polygon including a weird, disconnected linestring +impl Validify for Polygon { + type ValidResult = MultiPolygon; + fn split_into_valid(&self) -> Self::ValidResult { + let mp = MultiPolygon::new(vec![self.clone()]); + split::split_invalid_multipolygon(&mp) + } +} + +impl Validify for MultiPolygon { + type ValidResult = MultiPolygon; + fn split_into_valid(&self) -> Self::ValidResult { + split::split_invalid_multipolygon(self) + } +} diff --git a/geo/src/algorithm/validity/mod.rs b/geo/src/algorithm/validity/mod.rs new file mode 100644 index 000000000..29a647e07 --- /dev/null +++ b/geo/src/algorithm/validity/mod.rs @@ -0,0 +1,20 @@ +mod impls; +mod split; + +/// Trait which includes a bunch of methods for dealing with invalid polygons. +/// +/// There are a few implicit assumptions that are made when dealing with geo's types. You can read +/// more about the specific assumptions of a kind of geometry in it's own documentation. +/// +/// Some of geo's algorithms act in unexpected ways if those validity assumptions are not +/// respected. An example of this is `.make_ccw_winding()` on [`LineStrings`] which might fail if +/// the linestring doesn't include unique points. +pub trait Validify { + type ValidResult; + + /// splits invalid geometries into valid ones + fn split_into_valid(&self) -> Self::ValidResult; +} + +#[cfg(test)] +mod split_tests; diff --git a/geo/src/algorithm/validity/split/banana_utils.rs b/geo/src/algorithm/validity/split/banana_utils.rs new file mode 100644 index 000000000..1f7e3e782 --- /dev/null +++ b/geo/src/algorithm/validity/split/banana_utils.rs @@ -0,0 +1,582 @@ +use std::ops::RangeInclusive; + +use geo_types::*; + +use crate::{ClosestPoint, Contains, CoordsIter, GeoFloat, LinesIter}; + +use super::types::{fold_closest, ClosestPointInfo, ClosestPointPreciseInfo, ConnectionKind}; +use super::utils::{const_true, filter_points_not_creating_intersections, is_point_traversed_once}; + +/// Banana ranges contain index ranges of the parts of a banana. This handles all cases: +/// +/// - banana through connection point +/// ┌─────────────┐ +/// │ / \ │ +/// │ / \ │ +/// │ / \ │ +/// │ ─────── │ +/// │ │ +/// │ │ +/// │ │ +/// └─────────────┘ +/// - banana through connection line +/// ┌──────┬──────┐ +/// │ │ │ +/// │ │ │ +/// │ / \ │ +/// │ / \ │ +/// │ / \ │ +/// │ ─────── │ +/// │ │ +/// └─────────────┘ +/// - recursive bananas (bananas in bananas) +#[derive(Debug)] +pub struct BananaRanges { + /// range of indices, indexing into an exterior linestring + /// + /// the range points to one side of the banana + pub min_range: RangeInclusive, + /// range of indices, indexing into an exterior linestring + /// + /// the range points to the opposite side of the banana than the `min_range` field + /// + /// if this field here is none, it's assumed that every index not present in the first range is + /// included in the second range + pub max_range: Option>, +} + +/// given +/// +/// - a predicate for the index of the iterated point +/// - a predicate for the point itself +/// +/// returns all coords in the polygons exterior that fulfill these predicates +fn filter_out_coords_by( + index_predicate: impl Fn(&usize) -> bool, + point_pred: impl Fn(Point) -> bool, +) -> impl Fn(&LineString) -> Vec> { + move |ls| { + ls.points() + .enumerate() + .filter_map(|(i, p)| (index_predicate(&i) && point_pred(p)).then_some(p.0)) + .collect::>() + } +} + +/// this function checks whether the polygon is a banana polygon and if so, it'll try to find a +/// minimal connecting line which either: +/// +/// - connects a hole to the banana +/// - splits the banana (🍨🍌 yum) into two polygons, reducing the banananess* +/// +/// * we can't prove that a split banana isn't banana anymore (it could be a whole bunch of bananas +/// in the same poly), but continuous splitting should eventually split the polygon into a polygon +/// without bananas +pub(crate) fn find_closest_lines_for_banana( + p: &Polygon, +) -> Option> { + // find banana ranges + let banana_range = find_banana_ranges(p.exterior())?; + + let BananaRanges { + min_range, + max_range, + } = banana_range; + + let lines = p.lines_iter().collect::>(); + let traversed_once_pred = |p| is_point_traversed_once(&lines, p); + let ext = p.exterior(); + + // get valid coords based on banana ranges + let inner_candidates = + filter_out_coords_by(|i| min_range.contains(i), traversed_once_pred)(ext); + let outer_exclusion_range = max_range.unwrap_or(min_range); + let outer_candidates = + filter_out_coords_by(|i| !outer_exclusion_range.contains(i), traversed_once_pred)(ext); + + // get best connection possible + find_closest_lines_for_banana_inner(p, inner_candidates, outer_candidates) +} + +pub(crate) fn find_and_split_outer_banana(p: &Polygon) -> Option> { + // find banana ranges + let banana_range = find_banana_ranges(p.exterior())?; + + let BananaRanges { + min_range, + max_range, + } = banana_range; + + let ext = p.exterior(); + + // get valid coords based on banana ranges + let inner_candidates = filter_out_coords_by(|i| min_range.contains(i), const_true)(ext); + let outer_exclusion_range = max_range.unwrap_or(min_range); + let mut outer_candidates = filter_out_coords_by( + |i| !outer_exclusion_range.contains(i) || i == outer_exclusion_range.start(), + const_true, + )(ext); + outer_candidates.dedup(); + + // check if neither polygon includes the other and if so, just return the split apart polys + let poly1 = Polygon::new(LineString::new(inner_candidates), vec![]); + let poly2 = Polygon::new(LineString::new(outer_candidates), vec![]); + + let is_outer_banana = !poly1.contains(&poly2) && !poly2.contains(&poly1); + + is_outer_banana.then(|| MultiPolygon::new(vec![poly1, poly2])) +} + +/// for given source (start) and target (end) coords, find the smallest connection line between the +/// banana polygon and either: +/// +/// - itself +/// - a interior hole +fn find_closest_lines_for_banana_inner( + polygon: &Polygon, + source_coords: Vec>, + target_coords: Vec>, +) -> Option> { + // Define some helping closures + let lines = polygon + .lines_iter() + .filter(|l| l.start != l.end) + .collect::>(); + let target_linestring = LineString::new(target_coords); + let iter_targets = || { + std::iter::once(&target_linestring) + .chain(polygon.interiors()) + .enumerate() + }; + + let find_closest_for = |point_in_self: Point| { + iter_targets() + .map(|(id, ls)| { + let ps = filter_points_not_creating_intersections(&lines, point_in_self)(ls); + (id, ps) + }) + .map(|(id, ps)| (id, ps.closest_point(&point_in_self))) + .map(|(id, point_in_other)| ClosestPointInfo { + point_in_other, + point_in_self, + from_linestring: ConnectionKind::from_banana_index(id), + }) + .fold(None, fold_closest) + }; + + source_coords + .into_iter() + .map(Point::from) + .filter_map(find_closest_for) + .fold(None, fold_closest) + .and_then(ClosestPointPreciseInfo::from_unprecise) +} + +/// for one given banana polygon, traverses through the polygon and returns index ranges which mark +/// one side of the banana +fn find_banana_ranges(exterior: &LineString) -> Option { + let (_, ranges) = exterior.coords_iter().enumerate().fold( + (vec![], vec![]), + |(mut state_coords, mut ranges), (end, coord)| { + // add new range if the current coord is closing some opened range + ranges.extend( + state_coords + .iter() + .position(|&other| other == coord) + .filter(|start| end - start + 1 < exterior.0.len()) + .map(|start| start..=end), + ); + // keep track of existing coords + state_coords.push(coord); + (state_coords, ranges) + }, + ); + + let range_len = |range: &&RangeInclusive| *range.end() - *range.start(); + let min_range = ranges.iter().min_by_key(range_len).cloned()?; + + // based on the min range, try to find the closest range that's bigger than the min range + let max_range = get_max_range_from_min_range(exterior, &min_range); + + Some(BananaRanges { + min_range, + max_range, + }) +} + +/// given some minimum banana range, this function calculates the banana range of the other side of +/// a possible banana connection linestring. +/// +/// In case there is no banana linestring, this function returns None +fn get_max_range_from_min_range( + exterior: &LineString, + min_range: &RangeInclusive, +) -> Option> { + let start = *min_range.start(); + let end = *min_range.end(); + + let ext_iter = || exterior.0.iter(); + let ext_iter_indexed = || ext_iter().enumerate(); + + // get all points the are occuring more than once + let duplicates = ext_iter() + .filter(|c| ext_iter().skip(1).filter(|o| o == c).count() > 1) + .collect::>(); + + // in case of a banana connection linestring, go along that linestring and find the point + // that is furthest away from the min_range + fn try_take_last_duplicate_in<'a, F: GeoFloat + 'a>( + iter: impl Iterator)>, + duplicates: &[&Coord], + ) -> Option { + iter.take_while(|(_, c)| duplicates.contains(c)) + .last() + .map(|(i, _)| i) + } + + let rev_iter_at_start = ext_iter_indexed().rev().skip(exterior.0.len() - start); + let normal_iter_at_end = ext_iter_indexed().skip(end); + + // note that the early returns here make sure that we always have matching points on both + // sides. This prevents us creating a range of duplicate points on the one side of the + // connection line string while the other side actually features no duplicates + // + // in this case we just use the inverse min_range later on + let max_start = try_take_last_duplicate_in(rev_iter_at_start, duplicates.as_slice())?; + let max_end = try_take_last_duplicate_in(normal_iter_at_end, duplicates.as_slice())?; + + Some(max_start..=max_end) +} + +/// simple predicate to check whether a polygon is a banana polygon +pub(crate) fn is_banana_polygon(poly: &Polygon) -> bool { + let ext = poly.exterior(); + let lines = ext.lines().collect::>(); + ext.points().any(|p| !is_point_traversed_once(&lines, p)) +} + +#[cfg(test)] +mod banana_tests { + use pretty_env_logger::env_logger::*; + + use super::*; + + fn gen_basic_inner_banana() -> Polygon { + polygon!( + Coord { x: 0.0, y: 0.0 }, + Coord { x: 2.0, y: 0.0 }, + Coord { x: 2.0, y: 1.0 }, // banana start + Coord { x: 1.0, y: 1.0 }, + Coord { x: 1.0, y: 3.0 }, + Coord { x: 3.0, y: 3.0 }, + Coord { x: 3.0, y: 1.0 }, + Coord { x: 2.0, y: 1.0 }, // banana end + Coord { x: 2.0, y: 0.0 }, + Coord { x: 4.0, y: 0.0 }, + Coord { x: 4.0, y: 4.0 }, + Coord { x: 0.0, y: 4.0 }, + ) + } + + fn gen_long_connection_inner_banana() -> Polygon { + polygon!( + Coord { x: 0.0, y: 0.0 }, + Coord { x: 2.0, y: 0.0 }, + Coord { x: 2.25, y: 0.25 }, + Coord { x: 1.75, y: 0.75 }, + Coord { x: 2.0, y: 1.0 }, // banana start + Coord { x: 1.0, y: 1.0 }, + Coord { x: 1.0, y: 3.0 }, + Coord { x: 3.0, y: 3.0 }, + Coord { x: 3.0, y: 1.0 }, + Coord { x: 2.0, y: 1.0 }, // banana end + Coord { x: 1.75, y: 0.75 }, + Coord { x: 2.25, y: 0.25 }, + Coord { x: 2.0, y: 0.0 }, + Coord { x: 4.0, y: 0.0 }, + Coord { x: 4.0, y: 4.0 }, + Coord { x: 0.0, y: 4.0 }, + ) + } + + fn gen_basic_outer_banana() -> Polygon { + polygon!( + Coord { x: 0.0, y: 0.0 }, + Coord { x: 2.0, y: 0.0 }, + Coord { x: 2.0, y: -1.0 }, // banana start + Coord { x: 1.0, y: -1.0 }, + Coord { x: 1.0, y: -3.0 }, + Coord { x: 3.0, y: -3.0 }, + Coord { x: 3.0, y: -1.0 }, + Coord { x: 2.0, y: -1.0 }, // banana end + Coord { x: 2.0, y: 0.0 }, + Coord { x: 4.0, y: 0.0 }, + Coord { x: 4.0, y: 4.0 }, + Coord { x: 0.0, y: 4.0 }, + ) + } + + fn gen_long_connection_outer_banana() -> Polygon { + polygon!( + Coord { x: 0.0, y: 0.0 }, + Coord { x: 2.0, y: 0.0 }, + Coord { x: 2.25, y: -0.25 }, + Coord { x: 1.75, y: -0.75 }, + Coord { x: 2.0, y: -1.0 }, // banana start + Coord { x: 1.0, y: -1.0 }, + Coord { x: 1.0, y: -3.0 }, + Coord { x: 3.0, y: -3.0 }, + Coord { x: 3.0, y: -1.0 }, + Coord { x: 2.0, y: -1.0 }, // banana end + Coord { x: 1.75, y: -0.75 }, + Coord { x: 2.25, y: -0.25 }, + Coord { x: 2.0, y: 0.0 }, + Coord { x: 4.0, y: 0.0 }, + Coord { x: 4.0, y: 4.0 }, + Coord { x: 0.0, y: 4.0 }, + ) + } + + fn gen_recursive_banana() -> Polygon { + polygon!( + Coord { x: 0.0, y: 0.0 }, + Coord { x: 2.0, y: 0.0 }, + Coord { x: 2.0, y: 1.0 }, // level 1 banana start + Coord { x: 1.0, y: 1.0 }, + Coord { x: 1.0, y: 3.0 }, + Coord { x: 3.0, y: 3.0 }, + Coord { x: 3.0, y: 2.0 }, + Coord { x: 4.0, y: 2.0 }, // level 2 banana start + Coord { x: 4.0, y: 3.0 }, + Coord { x: 6.0, y: 3.0 }, + Coord { x: 6.0, y: 1.0 }, + Coord { x: 4.0, y: 1.0 }, + Coord { x: 4.0, y: 2.0 }, // level 2 banana end + Coord { x: 3.0, y: 2.0 }, + Coord { x: 3.0, y: 1.0 }, + Coord { x: 2.0, y: 1.0 }, // level 1 banana end + Coord { x: 2.0, y: 0.0 }, + Coord { x: 10.0, y: 0.0 }, + Coord { x: 10.0, y: 10.0 }, + Coord { x: 0.0, y: 10.0 }, + ) + } + + fn run_finder_test_for( + gen: impl Fn() -> Polygon, + expected_inner_range: RangeInclusive, + expected_outer_range: Option>, + ) { + _ = try_init(); + + let banana = gen(); + let banana_connection = find_banana_ranges(banana.exterior()); + + assert!(banana_connection.is_some()); + let banana_connection = banana_connection.unwrap(); + + assert_eq!(banana_connection.min_range, expected_inner_range); + assert_eq!( + banana.exterior().0[*banana_connection.min_range.start()], + banana.exterior().0[*banana_connection.min_range.end()] + ); + assert_eq!(banana_connection.max_range, expected_outer_range); + assert_eq!( + banana_connection + .max_range + .as_ref() + .map(|range| banana.exterior().0[*range.start()]), + banana_connection + .max_range + .as_ref() + .map(|range| banana.exterior().0[*range.end()]), + ); + } + + mod found_banana { + use super::*; + + #[test] + fn basic_inner_banana_found() { + run_finder_test_for(gen_basic_inner_banana, 2..=7, Some(1..=8)); + } + + #[test] + fn long_inner_banana_found() { + run_finder_test_for(gen_long_connection_inner_banana, 4..=9, Some(1..=12)); + } + + #[test] + fn basic_outer_banana_found() { + run_finder_test_for(gen_basic_outer_banana, 2..=7, Some(1..=8)); + } + + #[test] + fn long_outer_banana_found() { + run_finder_test_for(gen_long_connection_outer_banana, 4..=9, Some(1..=12)); + } + + #[test] + fn recursive_banana_found() { + run_finder_test_for(gen_recursive_banana, 7..=12, Some(6..=13)); + } + } + + fn run_fixing_test_inner( + gen: impl Fn() -> Polygon, + expected_closest: Option>, + ) { + _ = try_init(); + let banana = gen(); + let fixed = find_closest_lines_for_banana(&banana); + + assert_eq!(fixed, expected_closest); + } + + fn run_fixing_test_outer( + gen: impl Fn() -> Polygon, + expected_closest: Option>, + ) { + _ = try_init(); + let banana = gen(); + let fixed = find_and_split_outer_banana(&banana); + + assert_eq!(fixed, expected_closest); + } + + #[test] + fn basic_inner_banana_fixed() { + run_fixing_test_inner( + gen_basic_inner_banana, + Some(ClosestPointPreciseInfo { + from_linestring: ConnectionKind::Exterior, + point_in_self: Point::new(1.0, 1.0), + point_in_other: Point::new(0.0, 0.0), + }), + ); + } + + #[test] + fn long_inner_banana_fixed() { + run_fixing_test_inner( + gen_long_connection_inner_banana, + Some(ClosestPointPreciseInfo { + from_linestring: ConnectionKind::Exterior, + point_in_self: Point::new(1.0, 1.0), + point_in_other: Point::new(0.0, 0.0), + }), + ); + } + + #[test] + fn basic_outer_banana_fixed() { + let expected = [ + vec![ + Coord { x: 2.0, y: -1.0 }, + Coord { x: 1.0, y: -1.0 }, + Coord { x: 1.0, y: -3.0 }, + Coord { x: 3.0, y: -3.0 }, + Coord { x: 3.0, y: -1.0 }, + ], + vec![ + Coord { x: 0.0, y: 0.0 }, + Coord { x: 2.0, y: 0.0 }, + Coord { x: 4.0, y: 0.0 }, + Coord { x: 4.0, y: 4.0 }, + Coord { x: 0.0, y: 4.0 }, + ], + ] + .map(|coords| Polygon::new(LineString::new(coords), vec![])) + .to_vec(); + let expected = MultiPolygon::new(expected); + run_fixing_test_outer(gen_basic_outer_banana, Some(expected)); + } + + #[test] + fn long_outer_banana_fixed() { + let expected = [ + vec![ + Coord { x: 2.0, y: -1.0 }, + Coord { x: 1.0, y: -1.0 }, + Coord { x: 1.0, y: -3.0 }, + Coord { x: 3.0, y: -3.0 }, + Coord { x: 3.0, y: -1.0 }, + ], + vec![ + Coord { x: 0.0, y: 0.0 }, + Coord { x: 2.0, y: 0.0 }, + Coord { x: 4.0, y: 0.0 }, + Coord { x: 4.0, y: 4.0 }, + Coord { x: 0.0, y: 4.0 }, + ], + ] + .map(|coords| Polygon::new(LineString::new(coords), vec![])) + .to_vec(); + let expected = MultiPolygon::new(expected); + run_fixing_test_outer(gen_long_connection_outer_banana, Some(expected)); + } + + #[test] + fn real1() { + let real = || { + polygon! { + Coord { x: 0.0, y: 0.0 }, + Coord { x: 10.0, y: 0.0 }, + Coord { x: 10.0, y: 3.0 }, + Coord { x: 0.0, y: 3.0 }, + Coord { x: 0.0, y: 0.0 }, + Coord { x: 4.0, y: 1.0 }, + Coord { x: 4.0, y: 2.0 }, + Coord { x: 4.25, y: 2.0 }, + Coord { x: 4.75, y: 2.0 }, + Coord { x: 5.0, y: 2.0 }, + Coord { x: 5.0, y: 1.0 }, + Coord { x: 4.75, y: 1.0 }, + Coord { x: 4.75, y: 2.0 }, + Coord { x: 4.25, y: 2.0 }, + Coord { x: 4.25, y: 1.0 }, + Coord { x: 4.0, y: 1.0 }, + Coord { x: 0.0, y: 0.0 }, + } + }; + run_fixing_test_inner( + real, + Some(ClosestPointPreciseInfo { + from_linestring: ConnectionKind::Exterior, + point_in_self: Point::new(0.0, 3.0), + point_in_other: Point::new(4.0, 2.0), + }), + ); + } + + #[test] + fn real2() { + let real = || { + polygon! { + Coord { x: 4.0, y: 2.0 }, + Coord { x: 4.25, y: 2.0 }, + Coord { x: 4.75, y: 2.0 }, + Coord { x: 5.0, y: 2.0 }, + Coord { x: 5.0, y: 1.0 }, + Coord { x: 4.75, y: 1.0 }, + Coord { x: 4.75, y: 2.0 }, + Coord { x: 4.25, y: 2.0 }, + Coord { x: 4.25, y: 1.0 }, + Coord { x: 4.0, y: 1.0 }, + Coord { x: 0.0, y: 0.0 }, + Coord { x: 10.0, y: 0.0 }, + Coord { x: 10.0, y: 3.0 }, + Coord { x: 0.0, y: 3.0 }, + Coord { x: 4.0, y: 2.0 }, + } + }; + run_fixing_test_inner( + real, + Some(ClosestPointPreciseInfo { + from_linestring: ConnectionKind::Exterior, + point_in_self: Point::new(4.75, 1.0), + point_in_other: Point::new(4.25, 1.0), + }), + ); + } +} diff --git a/geo/src/algorithm/validity/split/mod.rs b/geo/src/algorithm/validity/split/mod.rs new file mode 100644 index 000000000..009b4425b --- /dev/null +++ b/geo/src/algorithm/validity/split/mod.rs @@ -0,0 +1,277 @@ +mod banana_utils; +mod types; +mod utils; + +use geo_types::*; + +use banana_utils::{find_and_split_outer_banana, find_closest_lines_for_banana, is_banana_polygon}; +use types::{ + fold_closest, fold_closest_precise, ClosestPointInfo, ClosestPointPreciseInfo, ConnectionKind, + PolyAndClosest, Retry, +}; +use utils::{ + create_linestring_between_points, filter_points_not_creating_intersections, + find_and_remove_polygon, is_contact_point_of_exterior_interior, is_point_traversed_once, + iter_between, prepare_winding, reassociate_holes, remove_nth_hole, remove_two_holes, +}; + +use crate::{ClosestPoint, GeoFloat, LinesIter}; + +/// find's two connections between the first hole of a polygon and the rest of the polygon in the +/// following ordering: +/// +/// - Let H: be the hole-linestring of the first hole of the poly +/// - Let R: be the rest of the polygon (exterior + every but the first hole) +/// +/// 1. find a point in `H` with the minimum distance to any point `R` and determine the connection to +/// this other point +/// 2. find a point in `H` with the maximum distance to the first point in the hole-linestring found +/// in 1. +/// 3. find the minimum distance connection between the point from 2. and any other point in `R` +fn find_closest_lines_for_first_hole( + polygon: &Polygon, +) -> Option> { + // Define some helping closures + let lines = polygon + .lines_iter() + .filter(|l| l.start != l.end) + .collect::>(); + let first_hole = polygon.interiors().first().cloned()?; + let iter_targets = || { + std::iter::once(polygon.exterior()) + .chain(polygon.interiors().iter().skip(1)) + .enumerate() + }; + // exclusion set is to prevent the algo from choosing the same point twice since that would + // potentially to banana polygons which are not considered valid + let find_closest_for = |point_in_self: Point| { + iter_targets() + .map(|(id, linestring)| { + let ps = + filter_points_not_creating_intersections(&lines, point_in_self)(linestring); + (id, ps) + }) + .map(|(id, target_points)| (id, target_points.closest_point(&point_in_self))) + .map(|(id, point_in_other)| ClosestPointInfo { + point_in_other, + point_in_self, + from_linestring: ConnectionKind::from_normal_index(id), + }) + .fold(None, fold_closest) + }; + + let valid_point_filter = |p: &Point| { + // this is to give the algo a chance to just cut the poly at the contact point + let is_contact = is_contact_point_of_exterior_interior(polygon, p); + // these are all other points that guarantee that the connection isn't made between another + // connection part + let is_traversed_once = is_point_traversed_once(&lines, *p); + is_contact || is_traversed_once + }; + + let closest_connection = first_hole + .points() + .filter(valid_point_filter) + .filter_map(find_closest_for) + .fold(None, fold_closest) + .and_then(ClosestPointPreciseInfo::from_unprecise)?; + + Some(closest_connection) +} + +/// based on the location of the two closest points found, dispatch to functions that use these +/// points to take one step into splitting up the polygon +fn handle_closest_points_found(args: PolyAndClosest) -> MultiPolygon { + handle_closest_points_found_normal_cases(args) + .unwrap_or_else(handle_closest_points_found_banana_cases) +} + +fn handle_closest_points_found_normal_cases(args: PolyAndClosest) -> Retry { + // dispatch to specialized functions + match args.closest.from_linestring { + // we're splitting up the polygon into two new polygons + ConnectionKind::Exterior => connect_hole_to_ext(args), + // in any other case it's enough to make progress by merging one hole + ConnectionKind::Interior(_) => merge_holes(args), + } +} + +fn handle_closest_points_found_banana_cases( + args: PolyAndClosest, +) -> MultiPolygon { + // dispatch to specialized functions + match args.closest.from_linestring { + // we're splitting up the polygon into two new polygons + ConnectionKind::Exterior => make_banana_split(args), + // in any other case it's enough to make progress by merging one hole + ConnectionKind::Interior(_) => connect_hole_to_banana(args), + } +} + +/// hmm, yummy 🍨🍌🍨 +fn make_banana_split( + PolyAndClosest { poly, closest }: PolyAndClosest, +) -> MultiPolygon { + let ext = poly.exterior().clone(); + + let ext_iter = iter_between(ext); + + let create_poly_with_idxs = |[a, b]: [Point; 2]| { + let mut ext_coords = ext_iter([a, b]); + ext_coords.dedup(); + Polygon::new(LineString::new(ext_coords), vec![]) + }; + + let p1 = create_poly_with_idxs([closest.point_in_self, closest.point_in_other]); + let p2 = create_poly_with_idxs([closest.point_in_other, closest.point_in_self]); + + let new_polys = reassociate_holes([p1, p2], poly.interiors().iter().cloned()); + + MultiPolygon::new(new_polys.to_vec()) +} + +fn connect_hole_to_banana( + PolyAndClosest { poly, closest }: PolyAndClosest, +) -> MultiPolygon { + // get the index of the second hole (which isn't the first hole) + let (c, index) = match closest.from_linestring { + ConnectionKind::Interior(index) => (closest, index), + _ => unreachable!("got ext connection when hole was expected"), + }; + + let ext = poly.exterior().clone(); + + let (poly_without_first_hole, hole) = + remove_nth_hole(&poly, index).expect("this hole exists since it produced the index"); + + let ext_iter = iter_between(ext); + let hole_iter = iter_between(hole); + + let create_poly_with_idxs = |[a, b]: [Point; 2]| { + let ext_part_1 = ext_iter([a, a]); + let ext_part_2 = hole_iter([b, b]); + let mut ext_coords = [ext_part_1.clone(), ext_part_2.clone()].concat(); + ext_coords.dedup(); + + Polygon::new(LineString::new(ext_coords), vec![]) + }; + + // ┌──────────x───────────┐ + // │ │ │ + // │ ┌───x───┐ │ + // │ │ │ │ + // │ │ │ │ + // │ └───────┘ │ + // │ │ + // └──────────────────────┘ + + let mut result_poly = create_poly_with_idxs([c.point_in_self, c.point_in_other]); + + for hole in poly_without_first_hole.interiors().iter().cloned() { + result_poly.interiors_push(hole) + } + + MultiPolygon::new(vec![result_poly]) +} + +/// connects a hole with two lines to the exterior. In this process, it's guaranteed that we split +/// up the polygon into two new polygons. The rest of the holes will be re-associated with the two +/// new polygons by containment tests +fn connect_hole_to_ext( + PolyAndClosest { poly, closest }: PolyAndClosest, +) -> Retry { + // winding is important here. It doesn't have to be exactly like this but the exterior and the + // hole have to have opposite winding orders for this to work + // + // If you get why take a look at the sketch further below and play the scenario through. + + let ext = poly.exterior().clone(); + + let (poly_without_first_hole, hole) = + remove_nth_hole(&poly, 0).ok_or(PolyAndClosest { poly, closest })?; + + let ext_iter = iter_between(ext); + let hole_iter = iter_between(hole); + + // ┌──────────x───────────┐ + // │ │ │ + // │ ┌───x───┐ │ + // │ │ │ │ + // │ │ │ │ + // │ └───────┘ │ + // │ │ + // └──────────────────────┘ + + let result_poly_exterior = create_linestring_between_points(ext_iter, hole_iter)([ + closest.point_in_other, + closest.point_in_self, + ]); + let mut result_poly = Polygon::new(result_poly_exterior, vec![]); + + for hole in poly_without_first_hole.interiors().iter().cloned() { + result_poly.interiors_push(hole) + } + + Ok(MultiPolygon::new(vec![result_poly])) +} + +/// This merges two holes by simply connecting them via a line +fn merge_holes(PolyAndClosest { poly, closest }: PolyAndClosest) -> Retry { + // get the index of the other hole (which isn't the first hole) + let (c, index) = match closest.from_linestring { + ConnectionKind::Interior(index) => (closest, index), + _ => return Err(PolyAndClosest { poly, closest }), + }; + + let (mut p, other_hole, first_hole) = + remove_two_holes(&poly, [index, 0]).ok_or(PolyAndClosest { poly, closest })?; + + let first_iter = iter_between(first_hole); + let other_iter = iter_between(other_hole); + + let hole = create_linestring_between_points(first_iter, other_iter)([ + c.point_in_self, + c.point_in_other, + ]); + + p.interiors_push(hole); + + Ok(MultiPolygon::new(vec![p])) +} + +/// splits a MultiPolygon with holes into MultiPolygon without holes. +pub fn split_invalid_multipolygon(mp: &MultiPolygon) -> MultiPolygon { + let mut mp = prepare_winding(mp.clone()); + // take one step, this ends in finite time since it'll eliminate exactly one hole per loop + // iteration. This means we'll run this exactly `num_holes` times + while let Some(p) = find_and_remove_polygon(&mut mp, |p| { + !p.interiors().is_empty() || is_banana_polygon(p) + }) { + if let Some(new_mp) = find_and_split_outer_banana(&p) { + mp.0.extend(new_mp); + continue; + } + + let closest_connection = find_closest_lines_for_banana(&p) + .into_iter() + .chain(find_closest_lines_for_first_hole(&p)) + .fold(None, fold_closest_precise); + let new_mp = if let Some(closest) = closest_connection { + handle_closest_points_found(PolyAndClosest { + poly: p.clone(), + closest, + }) + } else { + match remove_nth_hole(&p, 0) { + Some((p, _)) => MultiPolygon::new(vec![p]), + None => { + log::error!("Didn't find a hole although one should exist. Returning mid algo"); + return mp; + } + } + }; + + mp.0.extend(new_mp); + } + mp +} diff --git a/geo/src/algorithm/validity/split/types.rs b/geo/src/algorithm/validity/split/types.rs new file mode 100644 index 000000000..5f37e4eb3 --- /dev/null +++ b/geo/src/algorithm/validity/split/types.rs @@ -0,0 +1,176 @@ +use geo_types::*; + +use crate::{Closest, EuclideanDistance, GeoFloat}; + +/// Mostly the same as `ClosestPointPreciseInfo` with some more info on the closest point useful +/// for comparison of results but generally too unergonomic to deal with +/// +/// (since we would need to get the point out of the `Closest` context a lot of times, causing +/// unwraps all over the place) +#[derive(Debug, Clone, Copy)] +pub(crate) struct ClosestPointInfo { + pub(crate) from_linestring: ConnectionKind, + pub(crate) point_in_self: Point, + /// closest point to the `point_in_hole` point (with some extra information, `Closest` context) + pub(crate) point_in_other: Closest, +} + +/// This struct holds all the data that is needed to describe the closest connection between a +/// point in the first hole of the polygon and another point either in the exterior or another hole +/// of the polygon +#[derive(Debug, Clone, Copy, PartialEq, Eq)] +pub(crate) struct ClosestPointPreciseInfo { + /// this field classifies the target (exterior or other hole). If it's another hole it also + /// provides closer information which index it has in the `.interiors()` return vector + pub(crate) from_linestring: ConnectionKind, + /// the point in the first hole of the polygon. This data is needed to draw the connecting line + pub(crate) point_in_self: Point, + /// closest point to the `point_in_hole` point + pub(crate) point_in_other: Point, +} + +impl ClosestPointPreciseInfo { + // conversion function, fails if closest value was Indeterminate + pub(crate) fn from_unprecise(value: ClosestPointInfo) -> Option { + Some(ClosestPointPreciseInfo { + from_linestring: value.from_linestring, + point_in_self: value.point_in_self, + point_in_other: unpack_closest(value)?, + }) + } +} + +pub(crate) struct PolyAndClosest { + pub(crate) poly: Polygon, + pub(crate) closest: ClosestPointPreciseInfo, +} + +/// type that allows a retry with the same arguments in the case of an error +pub(crate) type Retry = Result, PolyAndClosest>; + +/// This struct holds all the data that is needed to describe the closest connection between a +/// point in the first hole of the polygon and another point either in the exterior or another hole +/// of the polygon +#[derive(Debug, Clone, PartialEq, Eq)] +pub(crate) struct ClosestPointPreciseInfoOuterBanana { + // indices of points of first poly + pub(crate) half_1: Vec, + // indices of points of second poly + pub(crate) half_2: Vec, +} + +/// Human readable enum to: +/// +/// - distinguish between exterior and interior cases +/// - hold information which interior we're talking about in the interior case +#[derive(Debug, Clone, Copy, PartialEq, Eq)] +pub(crate) enum ConnectionKind { + Exterior, + Interior(usize), +} + +impl ConnectionKind { + /// This function handles incoming indices as follows + /// + /// - index 0 => exterior + /// - index 1..n => hole 1..n + /// + /// Note that this excludes the first hole (index 0) + pub(crate) fn from_normal_index(n: usize) -> Self { + if n == 0 { + Self::Exterior + } else { + // in this case we enumerate all but the first + Self::Interior(n) + } + } + + /// This function handles incoming indices as follows + /// + /// - index 0 => exterior + /// - index 1..(n+1) => hole 0..n + /// + /// Note that this includes the first hole (index 0) + pub(crate) fn from_banana_index(n: usize) -> Self { + if n == 0 { + Self::Exterior + } else { + Self::Interior(n.saturating_sub(1)) + } + } +} + +/// returns: +/// +/// - None for `Indeterminate` cases +/// - `Some(input)` for valid cases +pub(crate) fn filter_out_indeterminate( + c: ClosestPointInfo, +) -> Option> { + match c.point_in_other { + Closest::Intersection(_) => Some(c), + Closest::SinglePoint(_) => Some(c), + Closest::Indeterminate => None, + } +} + +/// returns: +/// +/// - None for `Indeterminate` cases +/// - `Some(point)` for valid cases +pub(crate) fn unpack_closest(c: ClosestPointInfo) -> Option> { + match c.point_in_other { + Closest::Intersection(p) => Some(p), + Closest::SinglePoint(p) => Some(p), + Closest::Indeterminate => None, + } +} + +// copied from the impl of Closest::best_of_two and slightly adjusted here +/// compares two closest points and returns the best of the two, i.e. the one which is closer to +/// it's target +pub(crate) fn best_of_two( + a: ClosestPointInfo, + b: ClosestPointInfo, +) -> ClosestPointInfo { + let inner_pa = match a.point_in_other { + Closest::Indeterminate => return b, + Closest::Intersection(_) => return a, + Closest::SinglePoint(l) => l, + }; + let inner_pb = match b.point_in_other { + Closest::Indeterminate => return a, + Closest::Intersection(_) => return b, + Closest::SinglePoint(r) => r, + }; + + if inner_pa.euclidean_distance(&a.point_in_self) + <= inner_pb.euclidean_distance(&b.point_in_self) + { + a + } else { + b + } +} + +/// helper function to reduce a list of closest points to the minimum +pub(crate) fn fold_closest( + acc: Option>, + new: ClosestPointInfo, +) -> Option> { + let new_best = acc.map_or(new, |old| best_of_two(old, new)); + filter_out_indeterminate(new_best) +} + +/// helper function to reduce a list of closest points to the minimum +pub(crate) fn fold_closest_precise( + acc: Option>, + new: ClosestPointPreciseInfo, +) -> Option> { + acc.map_or(Some(new), |old| { + let new_dist = new.point_in_other.euclidean_distance(&new.point_in_self); + let old_dist = old.point_in_other.euclidean_distance(&old.point_in_self); + let res = if new_dist < old_dist { new } else { old }; + Some(res) + }) +} diff --git a/geo/src/algorithm/validity/split/utils.rs b/geo/src/algorithm/validity/split/utils.rs new file mode 100644 index 000000000..9f66fc8bb --- /dev/null +++ b/geo/src/algorithm/validity/split/utils.rs @@ -0,0 +1,210 @@ +use geo_types::*; + +use crate::line_intersection::line_intersection; +use crate::{Area, Contains, GeoFloat, LineIntersection, Scale, Winding}; + +/// given a linestring, returns a function which takes two points as an argument +/// +/// The returned function will return all the points between the two argument points also including +/// these two argument points at the start and end +pub(crate) fn iter_between( + ls: LineString, +) -> impl Fn([Point; 2]) -> Vec> { + move |[a, b]| { + // start + std::iter::once(a.0) + // between + .chain( + ls.0.iter() + .cycle() + .skip_while(move |&&c| c != a.0) + .skip(1) + .take_while(move |&&c| c != b.0) + .cloned(), + ) + // end + .chain(std::iter::once(b.0)) + .collect::>() + } +} + +/// Creates a `LineString` between two points with custom strategies for each half. +/// +/// A strategy is a function that collects all the points between it's two argument points. +/// +/// # Parameters +/// - `first_half_maker`: Closure for the first half strategy. +/// - `second_half_maker`: Closure for the second half strategy. +/// +/// # Returns +/// A closure generating a `LineString` from two points using the specified strategies. +pub(crate) fn create_linestring_between_points( + first_half_maker: FN, + second_half_maker: FN, +) -> impl Fn([Point; 2]) -> LineString +where + F: GeoFloat, + FN: Fn([Point; 2]) -> Vec>, +{ + move |[first_p, second_p]: [Point; 2]| { + let ext_part_1 = first_half_maker([first_p; 2]); + let ext_part_2 = second_half_maker([second_p; 2]); + let mut ext_coords = [ext_part_1.clone(), ext_part_2.clone()].concat(); + ext_coords.dedup(); + + LineString::new(ext_coords) + } +} + +/// given N polygons and a collection of potential holes, looks at which holes fit into which +/// polygon and assignes the holes correctly +pub(crate) fn reassociate_holes( + mut polys: [Polygon; N], + holes: impl IntoIterator>, +) -> [Polygon; N] { + holes.into_iter().for_each(|hole| { + if let Some(poly) = polys + .iter_mut() + .filter(|poly| poly.unsigned_area() > F::zero()) + .find(|poly| { + hole.0 + .iter() + .all(|p| poly.scale(F::from(0.99).unwrap()).contains(p)) + }) + { + poly.interiors_push(hole); + } + }); + polys +} + +/// remove the specified hole from the polygon and return the +/// +/// - the polygon without that hole +/// - the hole (optionally, if it exists) +pub(crate) fn remove_nth_hole( + p: &Polygon, + i: usize, +) -> Option<(Polygon, LineString)> { + let hole = p.interiors().iter().nth(i).cloned()?; + let p = Polygon::new( + p.exterior().clone(), + p.interiors() + .iter() + .take(i) + .chain(p.interiors().iter().skip(i + 1)) + .cloned() + .collect::>(), + ); + Some((p, hole)) +} + +/// finds a polygon fulfilling the predicate given to this function and removes it from the +/// MultiPolygon +pub(crate) fn find_and_remove_polygon( + mp: &mut MultiPolygon, + predicate: impl Fn(&Polygon) -> bool, +) -> Option> { + let idx = mp.iter().position(predicate)?; + Some(mp.0.remove(idx)) +} + +/// remove the two specified holes from the polygon and return the +/// +/// - the polygon without that hole +/// - the holes (optionally, if two exists) +/// +/// this function additionally takes care of removing the holes in the right order, to prevent +/// index invalidation +pub(crate) fn remove_two_holes( + poly: &Polygon, + idxs: [usize; 2], +) -> Option<(Polygon, LineString, LineString)> { + let [a, b] = idxs; + + let max = a.max(b); + let min = a.min(b); + + let (poly, other_hole) = remove_nth_hole(poly, max)?; + let (poly, first_hole) = remove_nth_hole(&poly, min)?; + + Some((poly, other_hole, first_hole)) +} + +/// predicate to check if a point is both present in the exterior and interior +pub(crate) fn is_contact_point_of_exterior_interior( + polygon: &Polygon, + point: &Point, +) -> bool { + polygon.exterior().contains(point) && polygon.interiors().iter().any(|int| int.contains(point)) +} + +/// predicate that checks whether the given point is only traveresd once +pub(crate) fn is_point_traversed_once(lines: &[Line], point: Point) -> bool { + lines + .iter() + .filter(|l| l.start == point.0 || l.end == point.0) + .count() + == 2 +} + +/// given a set of lines and a start point, returns a predicate that checks if for a given endpoint +/// a new line would intersect one of the given lines +pub(crate) fn exclude_points_forming_invalid_lines( + lines: &[Line], + point_in_hole: Point, +) -> impl Fn(&Point) -> bool + '_ { + move |p| { + // filter out len zero lines as they are invalid + if *p == point_in_hole { + return true; + } + + let new_line = Line::new(p.0, point_in_hole.0).scale(F::from(0.99).unwrap()); + !lines.iter().any(|old_line| { + let maybe_intersection = line_intersection(*old_line, new_line); + maybe_intersection + .filter(|kind| match kind { + // is_proper ensures they don't intersect in endpoints only + LineIntersection::SinglePoint { is_proper, .. } => *is_proper, + LineIntersection::Collinear { .. } => true, + }) + .is_some() + }) + } +} + +/// creates a closure that maps [`exclude_lines_intersecting`] over a range of points and only +/// return the valid ones +pub(crate) fn filter_points_not_creating_intersections( + lines: &[Line], + point_in_hole: Point, +) -> impl Fn(&LineString) -> MultiPoint + '_ { + move |ls| { + let points = ls + .points() + .filter(exclude_points_forming_invalid_lines(lines, point_in_hole)) + .collect::>(); + MultiPoint::new(points) + } +} + +/// The polygon splitting algorithm expects it's input polygon to be in a certain shape +/// +/// - the exterior needs to be ccw +/// - all interiors need to be cw +/// +/// Otherwise connecting two linestrings or splitting one up would fail in unexpected ways. That's +/// why this function is preparing these invariants. +pub(crate) fn prepare_winding(mut mp: MultiPolygon) -> MultiPolygon { + mp.iter_mut().for_each(|p| { + p.exterior_mut(|ls| ls.make_ccw_winding()); + p.interiors_mut(|lss| lss.iter_mut().for_each(|ls| ls.make_cw_winding())); + }); + mp +} + +/// function that always returns true +pub(crate) fn const_true(_: T) -> bool { + true +} diff --git a/geo/src/algorithm/validity/split_tests.rs b/geo/src/algorithm/validity/split_tests.rs new file mode 100644 index 000000000..b3f43b735 --- /dev/null +++ b/geo/src/algorithm/validity/split_tests.rs @@ -0,0 +1,481 @@ +use super::*; +use geo_types::*; + +use crate::{polygon, CoordsIter, GeoFloat, MapCoords, Winding}; +use pretty_env_logger::env_logger::*; + +fn one() -> Coord { + Coord { x: 1.0, y: 1.0 } +} + +fn test_consecutive_points_unique(mp: &MultiPolygon) { + mp.0.iter().for_each(|p| { + p.exterior() + .0 + .windows(2) + .map(|w| [w[0], w[1]]) + .for_each(|[a, b]| { + assert_ne!(a, b); + }); + }); +} + +#[test] +fn basic_stacked_rectangles_split() { + _ = try_init(); + let big_rect = Rect::new(Coord::zero(), one() * 3.0); + let small_rect = Rect::new(one(), one() * 2.0); + + let p = Polygon::new( + big_rect.to_polygon().exterior().clone(), + vec![small_rect.to_polygon().exterior().clone()], + ); + + let split = p.split_into_valid(); + + assert_eq!(split.0.len(), 2); + test_consecutive_points_unique(&split); +} + +#[test] +fn two_rectangles_inside_separate() { + _ = try_init(); + let big_rect = Rect::new(Coord::zero(), Coord { x: 10.0, y: 3.0 }); + let small_rect_1 = Rect::new(one(), one() * 2.0); + let small_rect_2 = Rect::new(Coord { x: 8.0, y: 1.0 }, Coord { x: 9.0, y: 2.0 }); + + let p = Polygon::new( + big_rect.to_polygon().exterior().clone(), + [small_rect_1, small_rect_2] + .map(|r| r.to_polygon().exterior().clone()) + .to_vec(), + ); + + let split = p.split_into_valid(); + + assert_eq!(split.0.len(), 3); + test_consecutive_points_unique(&split); +} + +#[test] +fn two_rectangles_inside_together() { + _ = try_init(); + let big_rect = Rect::new(Coord::zero(), Coord { x: 10.0, y: 3.0 }); + let small_rect_1 = Rect::new(Coord { x: 4.0, y: 1.0 }, Coord { x: 4.25, y: 2.0 }); + let small_rect_2 = Rect::new(Coord { x: 4.75, y: 1.0 }, Coord { x: 5.0, y: 2.0 }); + + let p = Polygon::new( + big_rect.to_polygon().exterior().clone(), + [small_rect_1, small_rect_2] + .map(|r| r.to_polygon().exterior().clone()) + .to_vec(), + ); + + let split = p.split_into_valid(); + + assert_eq!(split.0.len(), 3); + test_consecutive_points_unique(&split); +} + +#[test] +fn two_rectangles_inside_ext_and_hole_connecting() { + _ = try_init(); + let big_rect = Rect::new(Coord::zero(), Coord { x: 10.0, y: 3.0 }); + let small_rect_1 = Rect::new(Coord { x: 0.5, y: 1.0 }, Coord { x: 0.75, y: 2.0 }); + let small_rect_2 = Rect::new(Coord { x: 1.0, y: 1.0 }, Coord { x: 2.0, y: 2.0 }); + + let p = Polygon::new( + big_rect.to_polygon().exterior().clone(), + [small_rect_1, small_rect_2] + .map(|r| r.to_polygon().exterior().clone()) + .to_vec(), + ); + + let split = p.split_into_valid(); + + assert_eq!(split.0.len(), 3); + test_consecutive_points_unique(&split); +} + +#[test] +fn funny_case_two_holes() { + _ = try_init(); + + let p = polygon! { + exterior: [ + (x: 0.0, y: 0.0), + (x: 8.0, y: 0.0), + (x: 8.0, y: 9.0), + (x: -1.0, y: 9.0), + (x: 2.0, y: 4.0), + ], + interiors: [ + [ + (x: 4.0, y: 4.0), + (x: 6.0, y: 5.0), + (x: 3.0, y: 5.0), + ], + [ + (x: 3.0, y: 6.0), + (x: 7.0, y: 5.0), + (x: 7.0, y: 6.0), + ], + ] + } + .map_coords(|c| Coord { x: c.x, y: -c.y }); + + let split = p.split_into_valid(); + + assert_eq!(split.0.len(), 2); + test_consecutive_points_unique(&split); +} + +#[test] +fn poly_within_a_poly() { + _ = try_init(); + let big_rect = Rect::new(Coord::zero(), one() * 10.0); + let big_rect_hole = Rect::new(one(), one() * 9.0); + let small_rect = Rect::new(one() * 2.0, one() * 8.0); + let small_rect_hole = Rect::new(one() * 3.0, one() * 7.0); + + let outer = Polygon::new( + big_rect.to_polygon().exterior().clone(), + vec![big_rect_hole.to_polygon().exterior().clone()], + ); + let inner = Polygon::new( + small_rect.to_polygon().exterior().clone(), + vec![small_rect_hole.to_polygon().exterior().clone()], + ); + + let p = MultiPolygon::new(vec![outer, inner]); + + let split = p.split_into_valid(); + + assert_eq!(split.0.len(), 4); + test_consecutive_points_unique(&split); +} + +#[test] +fn poly_within_a_poly_within_a_poly() { + _ = try_init(); + let big_rect = Rect::new(Coord::zero(), one() * 10.0); + let big_rect_hole = Rect::new(one(), one() * 9.0); + let middle_rect = Rect::new(one() * 2.0, one() * 8.0); + let middle_rect_hole = Rect::new(one() * 3.0, one() * 7.0); + let small_rect = Rect::new(one() * 4.0, one() * 6.0); + let small_rect_hole = Rect::new(one() * 4.5, one() * 5.5); + + let outer = Polygon::new( + big_rect.to_polygon().exterior().clone(), + vec![big_rect_hole.to_polygon().exterior().clone()], + ); + let middle = Polygon::new( + middle_rect.to_polygon().exterior().clone(), + vec![middle_rect_hole.to_polygon().exterior().clone()], + ); + let inner = Polygon::new( + small_rect.to_polygon().exterior().clone(), + vec![small_rect_hole.to_polygon().exterior().clone()], + ); + + let p = MultiPolygon::new(vec![outer, middle, inner]); + + let split = p.split_into_valid(); + + assert_eq!(split.0.len(), 6); + test_consecutive_points_unique(&split); +} + +// if this were the case, this test would create a banana polygon which is basically not +// hole-less +#[test] +fn connection_to_exterior_isnt_same_point() { + _ = try_init(); + let big_rect = Rect::new(Coord::zero(), one() * 100.0); + let tri1 = Triangle::new( + Coord { x: 49.0, y: 1.0 }, + Coord { x: 49.5, y: 1.0 }, + Coord { x: 49.5, y: 25.0 }, + ); + let tri2 = Triangle::new( + Coord { y: 49.0, x: 1.0 }, + Coord { y: 49.5, x: 1.0 }, + Coord { y: 49.5, x: 25.0 }, + ); + + let p = Polygon::new( + big_rect.to_polygon().exterior().clone(), + [tri1, tri2] + .map(|t| t.to_polygon().exterior().clone()) + .to_vec(), + ); + + let split = p.split_into_valid(); + + assert_eq!(split.0.len(), 2); + test_consecutive_points_unique(&split); +} + +#[test] +fn connection_to_exterior_isnt_same_point_real_data() { + let p = polygon! { + exterior: [ + (x: -32.24332, y: -20.182356), + (x: -9.921326, y: -43.82522), + (x: -2.090634, y: -10.18567), + (x: -24.399603, y: -6.548216), + (x: -32.24332, y: -20.182356) + ], + interiors: [ + [ + (x: -15.338705, y: -25.178802), + (x: -16.906475, y: -22.304386), + (x: -14.203756, y: -21.644667), + (x: -15.338705, y: -25.178802), + ], + [ + (x: -14.332268, y: -16.889647), + (x: -25.32022, y: -16.11856), + (x: -18.89452, y: -13.0395775), + (x: -14.332268, y: -16.889647), + ], + ], + }; + + let split = p.split_into_valid(); + + assert_eq!(split.0.len(), 2); + test_consecutive_points_unique(&split); +} + +#[test] +fn connecting_line_doesnt_intersect_other_line() { + let p = polygon! { + exterior: [ + Coord { x: -16.599487, y: -15.243347, }, + Coord { x: 22.438654, y: -20.59876, }, + Coord { x: 13.0709915, y: 11.69747, }, + Coord { x: -18.707924, y: 4.3485994, }, + Coord { x: -18.472523, y: -10.658596, }, + Coord { x: -16.599487, y: -15.243347, }, + ], + interiors: [ + [ + Coord { x: -14.979336, y: -7.160144, }, + Coord { x: -6.7293396, y: -0.86014414, }, + Coord { x: 1.4206578, y: -5.460142, }, + Coord { x: -14.979336, y: -7.160144, }, + ], + [ + Coord { x: -12.000986, y: -11.615011, }, + Coord { x: 1.9206573, y: -8.860142, }, + Coord { x: -2.417613, y: -14.433651, }, + Coord { x: -12.000986, y: -11.615011, }, + ], + ], + }; + + let split = p.split_into_valid(); + + assert_eq!(split.0.len(), 2); + test_consecutive_points_unique(&split); +} + +#[test] +fn winding_sanity_check() { + // I think this caused a bug and if not we don't need to prove it again + + let mut t1 = LineString::new( + Triangle::from([ + Coord { + x: 8.237429, + y: 9.415966, + }, + Coord { + x: -4.9125605, + y: 8.0659685, + }, + Coord { + x: -0.6625643, + y: 11.465963, + }, + ]) + .to_array() + .to_vec(), + ); + let mut t2 = LineString::new( + Triangle::from([ + Coord { + x: -10.012558, + y: 1.6159716, + }, + Coord { + x: 4.0338764, + y: 1.6316788, + }, + Coord { + x: -1.9625626, + y: -1.4340262, + }, + ]) + .to_array() + .to_vec(), + ); + + t1.close(); + t1.make_cw_winding(); + t2.close(); + t2.make_cw_winding(); + + let ls = t1.coords_iter().chain(t2.coords_iter()).collect::>(); + + let p = Polygon::new(LineString::new(ls), vec![]); + + let mut ls = p.exterior().clone(); + assert!(!ls.is_cw()); + ls.make_cw_winding(); + // !!! 😮 !!! + assert!( + !ls.is_cw(), + "!!! IF THIS TEST FAILS THAT'S ACTUALLY A GOOD SIGN SINCE THE ASSERTION IS WEIRD !!!" + ); +} + +#[test] +fn holes_arent_used_twice() { + _ = try_init(); + let p = polygon! { + exterior: [ + Coord { x: -15.516099, y: 0.4316782, }, + Coord { x: -10.416106, y: -6.7683125, }, + Coord { x: 11.183868, y: -2.5683184, }, + Coord { x: 16.634357, y: 16.267822, }, + Coord { x: -13.116101, y: 14.781662, }, + Coord { x: -15.516099, y: 0.4316782, }, + ], + interiors: [ + [ + Coord { x: 8.237429, y: 9.415966, }, + Coord { x: -4.9125605, y: 8.0659685, }, + Coord { x: -0.6625643, y: 11.465963, }, + Coord { x: 8.237429, y: 9.415966, }, + ], + [ + Coord { x: -10.012558, y: 1.6159716, }, + Coord { x: 4.0338764, y: 1.6316788, }, + Coord { x: -1.9625626, y: -1.4340262, }, + Coord { x: -10.012558, y: 1.6159716, }, + ], + ], + }; + + let split = p.split_into_valid(); + + assert_eq!(split.0.len(), 2); + test_consecutive_points_unique(&split); +} + +#[test] +fn no_crossing_connecting_lines() { + _ = try_init(); + let p = polygon! { + exterior: [ + Coord { x: 30.180992, y: 131.80675, }, + Coord { x: -406.34613, y: -18.371037, }, + Coord { x: -30.24865, y: -56.683075, }, + Coord { x: 43.076523, y: -34.146957, }, + Coord { x: 43.471027, y: 12.403208, }, + Coord { x: 30.180992, y: 131.80675, }, + ], + interiors: [[ + Coord { x: -149.56546, y: -4.3648252, }, + Coord { x: -57.746937, y: 21.313248, }, + Coord { x: -59.631237, y: -8.495655, }, + Coord { x: -149.56546, y: -4.3648252, }, + ]], + }; + + let split = p.split_into_valid(); + + assert_eq!(split.0.len(), 2); + test_consecutive_points_unique(&split); +} + +#[test] +fn interior_one_vertex_on_exterior() { + _ = try_init(); + let p = polygon!( + exterior: [ + Coord { x: 0.0, y: 0.0 }, + Coord { x: 50.0, y: 0.0 }, + Coord { x: 100.0, y: 0.0 }, + Coord { x: 50.0, y: 100.0 } + ], + interiors: [[ + Coord { x: 50.0, y: 0.0 }, + Coord { x: 52.5, y: 2.5 }, + Coord { x: 47.5, y: 2.5 }, + ]], + ); + + let split = p.split_into_valid(); + + assert_eq!(split.0.len(), 2); + test_consecutive_points_unique(&split); +} + +#[test] +fn point_and_line_connection() { + _ = try_init(); + let p = polygon! { + Coord { x: -17.5161, y: 5.790348, }, + Coord { x: -3.7623396, y: -3.2200394, }, + Coord { x: -17.965897, y: -24.153437, }, + Coord { x: 25.049488, y: -16.02196, }, + Coord { x: 29.76677, y: 8.099451, }, + Coord { x: 5.395004, y: 8.099451, }, + Coord { x: 5.616444, y: 4.0554075, }, + Coord { x: 10.43572, y: -0.892385, }, + Coord { x: 3.3031864, y: -4.8763237, }, + Coord { x: 4.5291767, y: -8.794404, }, + Coord { x: 10.440827, y: -12.521313, }, + Coord { x: -5.6928024, y: -16.63536, }, + Coord { x: 4.5291767, y: -8.794404, }, + Coord { x: 3.3031864, y: -4.8763237, }, + Coord { x: -3.7623396, y: -3.2200394, }, + Coord { x: 5.616444, y: 4.0554075, }, + Coord { x: 5.395004, y: 8.099451, }, + Coord { x: -17.5161, y: 5.790348, }, + }; + + let split = p.split_into_valid(); + + assert_eq!(split.0.len(), 3); + test_consecutive_points_unique(&split); +} + +// this was found but we didn't really care about this case at the time. The issue here is that +// the algo splits the polygon into two polygons, which is good. But one of the return polygons +// is banana like which is kind of invalid and trips up algorithms +#[test] +fn interior_one_vertex_on_exterior_found_acidentally() { + _ = try_init(); + let p = polygon!( + exterior: [ + Coord { x: 0.0, y: 0.0 }, + Coord { x: 55.0, y: 0.0 }, + Coord { x: 100.0, y: 0.0 }, + Coord { x: 50.0, y: 100.0 } + ], + interiors: [[ + Coord { x: 50.0, y: 0.0 }, + Coord { x: 52.5, y: 2.5 }, + Coord { x: 47.5, y: 2.5 }, + ]], + ); + + let split = p.split_into_valid(); + + assert_eq!(split.0.len(), 2); + test_consecutive_points_unique(&split); +}