From 2691ef9386041571d7246431f2677bf1551bb122 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Stephan=20H=C3=BCgel?= Date: Fri, 17 Jun 2022 15:49:35 +0100 Subject: [PATCH 1/5] First pass at densification Remove stray JSON --- geo/src/algorithm/line_interpolate_point.rs | 160 +++++++++++++++++++- 1 file changed, 159 insertions(+), 1 deletion(-) diff --git a/geo/src/algorithm/line_interpolate_point.rs b/geo/src/algorithm/line_interpolate_point.rs index 8f9cebc1c..a7d3733b2 100644 --- a/geo/src/algorithm/line_interpolate_point.rs +++ b/geo/src/algorithm/line_interpolate_point.rs @@ -1,7 +1,137 @@ +use geo_types::{MultiLineString, MultiPolygon}; + use crate::coords_iter::CoordsIter; use std::ops::AddAssign; -use crate::{CoordFloat, EuclideanLength, Line, LineString, Point}; +use crate::{CoordFloat, EuclideanLength, Line, LineString, Point, Polygon}; + +/// Return a new linear geometry containing both existing and new interpolated coordinates with +/// a maximum distance of `max_distance` between them. +/// +/// # Examples +/// ``` +/// use geo::{coord, Line, LineString}; +/// use geo::line_interpolate_point::Densify; +/// +/// let line: Line = Line::new(coord! {x: 0.0, y: 6.0}, coord! {x: 1.0, y: 8.0}); +/// let correct: LineString = vec![[0.0, 6.0], [0.5, 7.0], [1.0, 8.0]].into(); +/// let max_dist = 2.0; +/// let densified = line.densify(max_dist); +/// assert_eq!(densified, correct); +///``` +pub trait Densify { + type Output; + + fn densify(&self, max_distance: F) -> Self::Output; +} + +// Helper for densification trait +fn densify_line(line: Line, container: &mut Vec>, max_distance: T) { + container.push(line.start_point()); + let num_segments = (line.euclidean_length() / max_distance).ceil() - T::one(); + if num_segments > T::zero() { + // distance "unit" for this line segment + let frac = T::one() / (T::from(num_segments).unwrap() + T::one()); + // conversion to int should be OK because we've already called ceil + (0..num_segments.to_u32().unwrap()) + .enumerate() + .for_each(|(seg, _)| { + // multiply distance unit by 1-indexed segment number to get correct offset + let np = line + .line_interpolate_point(frac * (T::from(seg).unwrap() + T::one())) + .unwrap(); + container.push(np); + }) + } +} + +impl Densify for MultiPolygon +where + T: CoordFloat, + Line: EuclideanLength, + LineString: EuclideanLength, +{ + type Output = MultiPolygon; + + fn densify(&self, max_distance: T) -> Self::Output { + MultiPolygon::new( + self.iter() + .map(|polygon| polygon.densify(max_distance)) + .collect(), + ) + } +} + +impl Densify for Polygon +where + T: CoordFloat, + Line: EuclideanLength, + LineString: EuclideanLength, +{ + type Output = Polygon; + + fn densify(&self, max_distance: T) -> Self::Output { + let densified_exterior = self.exterior().densify(max_distance); + let densified_interiors = self + .interiors() + .iter() + .map(|ring| ring.densify(max_distance)) + .collect(); + Polygon::new(densified_exterior, densified_interiors) + } +} + +impl Densify for MultiLineString +where + T: CoordFloat, + Line: EuclideanLength, + LineString: EuclideanLength, +{ + type Output = MultiLineString; + + fn densify(&self, max_distance: T) -> Self::Output { + MultiLineString::new( + self.iter() + .map(|linestring| linestring.densify(max_distance)) + .collect(), + ) + } +} + +impl Densify for LineString +where + T: CoordFloat, + Line: EuclideanLength, + LineString: EuclideanLength, +{ + type Output = LineString; + + fn densify(&self, max_distance: T) -> Self::Output { + let mut new_line = vec![]; + self.lines() + .for_each(|line| densify_line(line, &mut new_line, max_distance)); + // we're done, push the last coordinate on to finish + new_line.push(self.points().last().unwrap()); + LineString::from(new_line) + } +} + +impl Densify for Line +where + T: CoordFloat, + Line: EuclideanLength, + LineString: EuclideanLength, +{ + type Output = LineString; + + fn densify(&self, max_distance: T) -> Self::Output { + let mut new_line = vec![]; + densify_line(*self, &mut new_line, max_distance); + // we're done, push the last coordinate on to finish + new_line.push(self.end_point()); + LineString::from(new_line) + } +} /// Returns an option of the point that lies a given fraction along the line. /// @@ -109,6 +239,7 @@ where #[cfg(test)] mod test { + use super::*; use crate::{coord, point}; use crate::{ClosestPoint, LineLocatePoint}; @@ -295,4 +426,31 @@ mod test { _ => panic!("The closest point should be a SinglePoint"), // example chosen to not be an intersection }; } + #[test] + fn test_linestring_densify() { + let linestring: LineString = + vec![[-1.0, 0.0], [0.0, 0.0], [0.0, 6.0], [1.0, 8.0]].into(); + let correct: LineString = vec![ + [-1.0, 0.0], + [0.0, 0.0], + [0.0, 2.0], + [0.0, 4.0], + [0.0, 6.0], + [0.5, 7.0], + [1.0, 8.0], + ] + .into(); + let max_dist = 2.0; + let densified = linestring.densify(max_dist); + assert_eq!(densified, correct); + } + + #[test] + fn test_line_densify() { + let line: Line = Line::new(coord! {x: 0.0, y: 6.0}, coord! {x: 1.0, y: 8.0}); + let correct: LineString = vec![[0.0, 6.0], [0.5, 7.0], [1.0, 8.0]].into(); + let max_dist = 2.0; + let densified = line.densify(max_dist); + assert_eq!(densified, correct); + } } From b9cd493c91a692f1f137cbf4df9b60c5cca83860 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Stephan=20H=C3=BCgel?= Date: Fri, 17 Jun 2022 19:26:42 +0100 Subject: [PATCH 2/5] Move to own module, simplify logic, document Also add remaining geometries --- geo/src/algorithm/densify.rs | 193 ++++++++++++++++++++ geo/src/algorithm/line_interpolate_point.rs | 160 +--------------- geo/src/algorithm/mod.rs | 4 + geo/src/lib.rs | 1 + 4 files changed, 199 insertions(+), 159 deletions(-) create mode 100644 geo/src/algorithm/densify.rs diff --git a/geo/src/algorithm/densify.rs b/geo/src/algorithm/densify.rs new file mode 100644 index 000000000..dabd5eef0 --- /dev/null +++ b/geo/src/algorithm/densify.rs @@ -0,0 +1,193 @@ +use crate::{ + CoordFloat, EuclideanLength, Line, LineInterpolatePoint, LineString, MultiLineString, + MultiPolygon, Point, Polygon, Rect, Triangle, +}; + +/// Return a new linear geometry containing both existing and new interpolated coordinates with +/// a maximum distance of `max_distance` between them. +/// +/// Note: `max_distance` must be greater than 0. +/// +/// # Examples +/// ``` +/// use geo::{coord, Line, LineString}; +/// use geo::line_interpolate_point::Densify; +/// +/// let line: Line = Line::new(coord! {x: 0.0, y: 6.0}, coord! {x: 1.0, y: 8.0}); +/// let correct: LineString = vec![[0.0, 6.0], [0.5, 7.0], [1.0, 8.0]].into(); +/// let max_dist = 2.0; +/// let densified = line.densify(max_dist); +/// assert_eq!(densified, correct); +///``` +pub trait Densify { + type Output; + + fn densify(&self, max_distance: F) -> Self::Output; +} + +// Helper for densification trait +fn densify_line(line: Line, container: &mut Vec>, max_distance: T) { + assert!(max_distance > T::zero()); + container.push(line.start_point()); + let num_segments = (line.euclidean_length() / max_distance) + .ceil() + .to_u64() + .unwrap(); + // distance "unit" for this line segment + let frac = T::one() / T::from(num_segments).unwrap(); + for segment_idx in 1..num_segments { + let ratio = frac * T::from(segment_idx).unwrap(); + let interpolated_point = line + .line_interpolate_point(ratio) + .expect("ratio should be between 0..1"); + container.push(interpolated_point); + } +} + +impl Densify for MultiPolygon +where + T: CoordFloat, + Line: EuclideanLength, + LineString: EuclideanLength, +{ + type Output = MultiPolygon; + + fn densify(&self, max_distance: T) -> Self::Output { + MultiPolygon::new( + self.iter() + .map(|polygon| polygon.densify(max_distance)) + .collect(), + ) + } +} + +impl Densify for Polygon +where + T: CoordFloat, + Line: EuclideanLength, + LineString: EuclideanLength, +{ + type Output = Polygon; + + fn densify(&self, max_distance: T) -> Self::Output { + let densified_exterior = self.exterior().densify(max_distance); + let densified_interiors = self + .interiors() + .iter() + .map(|ring| ring.densify(max_distance)) + .collect(); + Polygon::new(densified_exterior, densified_interiors) + } +} + +impl Densify for MultiLineString +where + T: CoordFloat, + Line: EuclideanLength, + LineString: EuclideanLength, +{ + type Output = MultiLineString; + + fn densify(&self, max_distance: T) -> Self::Output { + MultiLineString::new( + self.iter() + .map(|linestring| linestring.densify(max_distance)) + .collect(), + ) + } +} + +impl Densify for LineString +where + T: CoordFloat, + Line: EuclideanLength, + LineString: EuclideanLength, +{ + type Output = LineString; + + fn densify(&self, max_distance: T) -> Self::Output { + let mut new_line = vec![]; + self.lines() + .for_each(|line| densify_line(line, &mut new_line, max_distance)); + // we're done, push the last coordinate on to finish + new_line.push(self.points().last().unwrap()); + LineString::from(new_line) + } +} + +impl Densify for Line +where + T: CoordFloat, + Line: EuclideanLength, + LineString: EuclideanLength, +{ + type Output = LineString; + + fn densify(&self, max_distance: T) -> Self::Output { + let mut new_line = vec![]; + densify_line(*self, &mut new_line, max_distance); + // we're done, push the last coordinate on to finish + new_line.push(self.end_point()); + LineString::from(new_line) + } +} + +impl Densify for Triangle +where + T: CoordFloat, + Line: EuclideanLength, + LineString: EuclideanLength, +{ + type Output = Polygon; + + fn densify(&self, max_distance: T) -> Self::Output { + self.to_polygon().densify(max_distance) + } +} + +impl Densify for Rect +where + T: CoordFloat, + Line: EuclideanLength, + LineString: EuclideanLength, +{ + type Output = Polygon; + + fn densify(&self, max_distance: T) -> Self::Output { + self.to_polygon().densify(max_distance) + } +} + +#[cfg(test)] +mod tests { + use super::*; + use crate::coord; + + #[test] + fn test_linestring_densify() { + let linestring: LineString = + vec![[-1.0, 0.0], [0.0, 0.0], [0.0, 6.0], [1.0, 8.0]].into(); + let correct: LineString = vec![ + [-1.0, 0.0], + [0.0, 0.0], + [0.0, 2.0], + [0.0, 4.0], + [0.0, 6.0], + [0.5, 7.0], + [1.0, 8.0], + ] + .into(); + let max_dist = 2.0; + let densified = linestring.densify(max_dist); + assert_eq!(densified, correct); + } + + #[test] + fn test_line_densify() { + let line: Line = Line::new(coord! {x: 0.0, y: 6.0}, coord! {x: 1.0, y: 8.0}); + let correct: LineString = vec![[0.0, 6.0], [0.5, 7.0], [1.0, 8.0]].into(); + let max_dist = 2.0; + let densified = line.densify(max_dist); + assert_eq!(densified, correct); + } +} diff --git a/geo/src/algorithm/line_interpolate_point.rs b/geo/src/algorithm/line_interpolate_point.rs index a7d3733b2..fdf2cdd2a 100644 --- a/geo/src/algorithm/line_interpolate_point.rs +++ b/geo/src/algorithm/line_interpolate_point.rs @@ -1,138 +1,7 @@ -use geo_types::{MultiLineString, MultiPolygon}; - use crate::coords_iter::CoordsIter; +use crate::{CoordFloat, EuclideanLength, Line, LineString, Point}; use std::ops::AddAssign; -use crate::{CoordFloat, EuclideanLength, Line, LineString, Point, Polygon}; - -/// Return a new linear geometry containing both existing and new interpolated coordinates with -/// a maximum distance of `max_distance` between them. -/// -/// # Examples -/// ``` -/// use geo::{coord, Line, LineString}; -/// use geo::line_interpolate_point::Densify; -/// -/// let line: Line = Line::new(coord! {x: 0.0, y: 6.0}, coord! {x: 1.0, y: 8.0}); -/// let correct: LineString = vec![[0.0, 6.0], [0.5, 7.0], [1.0, 8.0]].into(); -/// let max_dist = 2.0; -/// let densified = line.densify(max_dist); -/// assert_eq!(densified, correct); -///``` -pub trait Densify { - type Output; - - fn densify(&self, max_distance: F) -> Self::Output; -} - -// Helper for densification trait -fn densify_line(line: Line, container: &mut Vec>, max_distance: T) { - container.push(line.start_point()); - let num_segments = (line.euclidean_length() / max_distance).ceil() - T::one(); - if num_segments > T::zero() { - // distance "unit" for this line segment - let frac = T::one() / (T::from(num_segments).unwrap() + T::one()); - // conversion to int should be OK because we've already called ceil - (0..num_segments.to_u32().unwrap()) - .enumerate() - .for_each(|(seg, _)| { - // multiply distance unit by 1-indexed segment number to get correct offset - let np = line - .line_interpolate_point(frac * (T::from(seg).unwrap() + T::one())) - .unwrap(); - container.push(np); - }) - } -} - -impl Densify for MultiPolygon -where - T: CoordFloat, - Line: EuclideanLength, - LineString: EuclideanLength, -{ - type Output = MultiPolygon; - - fn densify(&self, max_distance: T) -> Self::Output { - MultiPolygon::new( - self.iter() - .map(|polygon| polygon.densify(max_distance)) - .collect(), - ) - } -} - -impl Densify for Polygon -where - T: CoordFloat, - Line: EuclideanLength, - LineString: EuclideanLength, -{ - type Output = Polygon; - - fn densify(&self, max_distance: T) -> Self::Output { - let densified_exterior = self.exterior().densify(max_distance); - let densified_interiors = self - .interiors() - .iter() - .map(|ring| ring.densify(max_distance)) - .collect(); - Polygon::new(densified_exterior, densified_interiors) - } -} - -impl Densify for MultiLineString -where - T: CoordFloat, - Line: EuclideanLength, - LineString: EuclideanLength, -{ - type Output = MultiLineString; - - fn densify(&self, max_distance: T) -> Self::Output { - MultiLineString::new( - self.iter() - .map(|linestring| linestring.densify(max_distance)) - .collect(), - ) - } -} - -impl Densify for LineString -where - T: CoordFloat, - Line: EuclideanLength, - LineString: EuclideanLength, -{ - type Output = LineString; - - fn densify(&self, max_distance: T) -> Self::Output { - let mut new_line = vec![]; - self.lines() - .for_each(|line| densify_line(line, &mut new_line, max_distance)); - // we're done, push the last coordinate on to finish - new_line.push(self.points().last().unwrap()); - LineString::from(new_line) - } -} - -impl Densify for Line -where - T: CoordFloat, - Line: EuclideanLength, - LineString: EuclideanLength, -{ - type Output = LineString; - - fn densify(&self, max_distance: T) -> Self::Output { - let mut new_line = vec![]; - densify_line(*self, &mut new_line, max_distance); - // we're done, push the last coordinate on to finish - new_line.push(self.end_point()); - LineString::from(new_line) - } -} - /// Returns an option of the point that lies a given fraction along the line. /// /// If the given fraction is @@ -426,31 +295,4 @@ mod test { _ => panic!("The closest point should be a SinglePoint"), // example chosen to not be an intersection }; } - #[test] - fn test_linestring_densify() { - let linestring: LineString = - vec![[-1.0, 0.0], [0.0, 0.0], [0.0, 6.0], [1.0, 8.0]].into(); - let correct: LineString = vec![ - [-1.0, 0.0], - [0.0, 0.0], - [0.0, 2.0], - [0.0, 4.0], - [0.0, 6.0], - [0.5, 7.0], - [1.0, 8.0], - ] - .into(); - let max_dist = 2.0; - let densified = linestring.densify(max_dist); - assert_eq!(densified, correct); - } - - #[test] - fn test_line_densify() { - let line: Line = Line::new(coord! {x: 0.0, y: 6.0}, coord! {x: 1.0, y: 8.0}); - let correct: LineString = vec![[0.0, 6.0], [0.5, 7.0], [1.0, 8.0]].into(); - let max_dist = 2.0; - let densified = line.densify(max_dist); - assert_eq!(densified, correct); - } } diff --git a/geo/src/algorithm/mod.rs b/geo/src/algorithm/mod.rs index 66c08f781..3bfd7ed55 100644 --- a/geo/src/algorithm/mod.rs +++ b/geo/src/algorithm/mod.rs @@ -191,3 +191,7 @@ pub mod sweep; /// Boolean Ops; pub mod bool_ops; pub use bool_ops::{BooleanOps, OpType}; + +/// Densify linear geometry components +pub mod densify; +pub use densify::Densify; diff --git a/geo/src/lib.rs b/geo/src/lib.rs index 3803d262d..c15702934 100644 --- a/geo/src/lib.rs +++ b/geo/src/lib.rs @@ -136,6 +136,7 @@ //! - **[`HaversineIntermediate`](HaversineIntermediate)**: //! - **[`Proj`](Proj)**: Project geometries with the `proj` crate (requires the `use-proj` feature) //! - **[`ChaikinSmoothing`](ChaikinSmoothing)**: Smoothen `LineString`, `Polygon`, `MultiLineString` and `MultiPolygon` using Chaikins algorithm. +//! - **[`Densify`](Densify)**: Densify linear geometry components by interpolating points //! //! # Features //! From 2c643c310fa7b5d56f4c92f6e590f4fea1d7a503 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Stephan=20H=C3=BCgel?= Date: Fri, 17 Jun 2022 19:33:51 +0100 Subject: [PATCH 3/5] Fix export in example --- geo/src/algorithm/densify.rs | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/geo/src/algorithm/densify.rs b/geo/src/algorithm/densify.rs index dabd5eef0..a057ae88c 100644 --- a/geo/src/algorithm/densify.rs +++ b/geo/src/algorithm/densify.rs @@ -11,7 +11,7 @@ use crate::{ /// # Examples /// ``` /// use geo::{coord, Line, LineString}; -/// use geo::line_interpolate_point::Densify; +/// use geo::Densify; /// /// let line: Line = Line::new(coord! {x: 0.0, y: 6.0}, coord! {x: 1.0, y: 8.0}); /// let correct: LineString = vec![[0.0, 6.0], [0.5, 7.0], [1.0, 8.0]].into(); From 9b35062d9699332f67b424f8c9e1d1abe85d7129 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Stephan=20H=C3=BCgel?= Date: Sat, 18 Jun 2022 00:49:31 +0100 Subject: [PATCH 4/5] Update changelog --- geo/CHANGES.md | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/geo/CHANGES.md b/geo/CHANGES.md index fe3543e2c..7f6a1e8c8 100644 --- a/geo/CHANGES.md +++ b/geo/CHANGES.md @@ -1,5 +1,9 @@ # Changes +## Unreleased +* Add densification algorithm for linear geometry components + * + ## 0.21.0 * Boolean operations for `Polygon`s and `MultiPolygon`s: intersect, union, xor, From 0336b607a0c253e52bde188c6b358568f648b0a7 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Stephan=20H=C3=BCgel?= Date: Sat, 18 Jun 2022 09:58:34 +0100 Subject: [PATCH 5/5] Add Polygon densify test --- geo/src/algorithm/densify.rs | 46 +++++++++++++++++++++++++++++++++++- 1 file changed, 45 insertions(+), 1 deletion(-) diff --git a/geo/src/algorithm/densify.rs b/geo/src/algorithm/densify.rs index a057ae88c..8bfd422b0 100644 --- a/geo/src/algorithm/densify.rs +++ b/geo/src/algorithm/densify.rs @@ -161,7 +161,51 @@ where #[cfg(test)] mod tests { use super::*; - use crate::coord; + use crate::{coord, Coordinate}; + + #[test] + fn test_polygon_densify() { + let linestring: LineString = + vec![[-5.0, 0.0], [0.0, 5.0], [5.0, 0.0], [-5.0, 0.0]].into(); + let interior: LineString = + vec![[-3.0, 0.0], [0.0, 3.0], [3.0, 0.0], [-3.0, 0.0]].into(); + let polygon = Polygon::new(linestring, vec![interior]); + let correct_ext: LineString = LineString(vec![ + Coordinate { x: -5.0, y: 0.0 }, + Coordinate { x: -3.75, y: 1.25 }, + Coordinate { x: -2.5, y: 2.5 }, + Coordinate { x: -1.25, y: 3.75 }, + Coordinate { x: 0.0, y: 5.0 }, + Coordinate { x: 1.25, y: 3.75 }, + Coordinate { x: 2.5, y: 2.5 }, + Coordinate { x: 3.75, y: 1.25 }, + Coordinate { x: 5.0, y: 0.0 }, + Coordinate { x: 3.0, y: 0.0 }, + Coordinate { x: 1.0, y: 0.0 }, + Coordinate { + x: -1.0000000000000009, + y: 0.0, + }, + Coordinate { x: -3.0, y: 0.0 }, + Coordinate { x: -5.0, y: 0.0 }, + ]); + let correct_int: LineString = LineString(vec![ + Coordinate { x: -3.0, y: 0.0 }, + Coordinate { x: -2.0, y: 1.0 }, + Coordinate { x: -1.0, y: 2.0 }, + Coordinate { x: 0.0, y: 3.0 }, + Coordinate { x: 1.0, y: 2.0 }, + Coordinate { x: 2.0, y: 1.0 }, + Coordinate { x: 3.0, y: 0.0 }, + Coordinate { x: 1.0, y: 0.0 }, + Coordinate { x: -1.0, y: 0.0 }, + Coordinate { x: -3.0, y: 0.0 }, + ]); + let correct_polygon = Polygon::new(correct_ext, vec![correct_int]); + let max_dist = 2.0; + let densified = polygon.densify(max_dist); + assert_eq!(densified, correct_polygon); + } #[test] fn test_linestring_densify() {