geo/algorithm/
closest_point.rs

1use crate::algorithm::{EuclideanLength, Intersects};
2use crate::geometry::*;
3use crate::Closest;
4use crate::GeoFloat;
5
6use std::iter;
7
8/// Find the closest `Point` between a given geometry and an input `Point`.
9/// The closest point may intersect the geometry, be a single
10/// point, or be indeterminate, as indicated by the value of the returned enum.
11///
12/// # Examples
13///
14/// We have a horizontal line which goes through `(-50, 0) -> (50, 0)`,
15/// and want to find the closest point to the point `(0, 100)`.
16/// Drawn on paper, the point on the line which is closest to `(0, 100)` is the origin (0, 0).
17///
18/// ```rust
19/// # use geo::ClosestPoint;
20/// # use geo::{Point, Line, Closest};
21/// let p: Point<f32> = Point::new(0.0, 100.0);
22/// let horizontal_line: Line<f32> = Line::new(Point::new(-50.0, 0.0), Point::new(50.0, 0.0));
23///
24/// let closest = horizontal_line.closest_point(&p);
25/// assert_eq!(closest, Closest::SinglePoint(Point::new(0.0, 0.0)));
26/// ```
27pub trait ClosestPoint<F: GeoFloat, Rhs = Point<F>> {
28    /// Find the closest point between `self` and `p`.
29    fn closest_point(&self, p: &Rhs) -> Closest<F>;
30}
31
32impl<'a, F, C> ClosestPoint<F> for &'a C
33where
34    C: ClosestPoint<F>,
35    F: GeoFloat,
36{
37    fn closest_point(&self, p: &Point<F>) -> Closest<F> {
38        (*self).closest_point(p)
39    }
40}
41
42impl<F: GeoFloat> ClosestPoint<F> for Point<F> {
43    fn closest_point(&self, p: &Self) -> Closest<F> {
44        if self == p {
45            Closest::Intersection(*self)
46        } else {
47            Closest::SinglePoint(*self)
48        }
49    }
50}
51
52#[allow(clippy::many_single_char_names)]
53impl<F: GeoFloat> ClosestPoint<F> for Line<F> {
54    fn closest_point(&self, p: &Point<F>) -> Closest<F> {
55        let line_length = self.euclidean_length();
56        if line_length == F::zero() {
57            // if we've got a zero length line, technically the entire line
58            // is the closest point...
59            return Closest::Indeterminate;
60        }
61
62        // For some line AB, there is some point, C, which will be closest to
63        // P. The line AB will be perpendicular to CP.
64        //
65        // Line equation: P = start + t * (end - start)
66
67        let direction_vector = Point::from(self.end - self.start);
68        let to_p = Point::from(p.0 - self.start);
69
70        let t = to_p.dot(direction_vector) / direction_vector.dot(direction_vector);
71
72        // check the cases where the closest point is "outside" the line
73        if t < F::zero() {
74            return Closest::SinglePoint(self.start.into());
75        } else if t > F::one() {
76            return Closest::SinglePoint(self.end.into());
77        }
78
79        let x = direction_vector.x();
80        let y = direction_vector.y();
81        let c = Point::from(self.start + (t * x, t * y).into());
82
83        if self.intersects(p) {
84            Closest::Intersection(c)
85        } else {
86            Closest::SinglePoint(c)
87        }
88    }
89}
90
91/// A generic function which takes some iterator of points and gives you the
92/// "best" `Closest` it can find. Where "best" is the first intersection or
93/// the `Closest::SinglePoint` which is closest to `p`.
94///
95/// If the iterator is empty, we get `Closest::Indeterminate`.
96fn closest_of<C, F, I>(iter: I, p: Point<F>) -> Closest<F>
97where
98    F: GeoFloat,
99    I: IntoIterator<Item = C>,
100    C: ClosestPoint<F>,
101{
102    let mut best = Closest::Indeterminate;
103
104    for element in iter {
105        let got = element.closest_point(&p);
106        best = got.best_of_two(&best, p);
107        if matches!(best, Closest::Intersection(_)) {
108            // short circuit - nothing can be closer than an intersection
109            return best;
110        }
111    }
112
113    best
114}
115
116impl<F: GeoFloat> ClosestPoint<F> for LineString<F> {
117    fn closest_point(&self, p: &Point<F>) -> Closest<F> {
118        closest_of(self.lines(), *p)
119    }
120}
121
122impl<F: GeoFloat> ClosestPoint<F> for Polygon<F> {
123    fn closest_point(&self, p: &Point<F>) -> Closest<F> {
124        if self.intersects(p) {
125            return Closest::Intersection(*p);
126        }
127        let prospectives = self.interiors().iter().chain(iter::once(self.exterior()));
128        closest_of(prospectives, *p)
129    }
130}
131
132impl<F: GeoFloat> ClosestPoint<F> for Coord<F> {
133    fn closest_point(&self, p: &Point<F>) -> Closest<F> {
134        Point::from(*self).closest_point(p)
135    }
136}
137
138impl<F: GeoFloat> ClosestPoint<F> for Triangle<F> {
139    fn closest_point(&self, p: &Point<F>) -> Closest<F> {
140        if self.intersects(p) {
141            return Closest::Intersection(*p);
142        }
143        closest_of(self.to_lines(), *p)
144    }
145}
146
147impl<F: GeoFloat> ClosestPoint<F> for Rect<F> {
148    fn closest_point(&self, p: &Point<F>) -> Closest<F> {
149        if self.intersects(p) {
150            return Closest::Intersection(*p);
151        }
152        closest_of(self.to_lines(), *p)
153    }
154}
155
156impl<F: GeoFloat> ClosestPoint<F> for MultiPolygon<F> {
157    fn closest_point(&self, p: &Point<F>) -> Closest<F> {
158        closest_of(self.iter(), *p)
159    }
160}
161
162impl<F: GeoFloat> ClosestPoint<F> for MultiPoint<F> {
163    fn closest_point(&self, p: &Point<F>) -> Closest<F> {
164        closest_of(self.iter(), *p)
165    }
166}
167
168impl<F: GeoFloat> ClosestPoint<F> for MultiLineString<F> {
169    fn closest_point(&self, p: &Point<F>) -> Closest<F> {
170        closest_of(self.iter(), *p)
171    }
172}
173
174impl<F: GeoFloat> ClosestPoint<F> for GeometryCollection<F> {
175    fn closest_point(&self, p: &Point<F>) -> Closest<F> {
176        closest_of(self.iter(), *p)
177    }
178}
179impl<F: GeoFloat> ClosestPoint<F> for Geometry<F> {
180    crate::geometry_delegate_impl! {
181        fn closest_point(&self, p: &Point<F>) -> Closest<F>;
182    }
183}
184
185#[cfg(test)]
186mod tests {
187    use super::*;
188    use crate::algorithm::{Contains, Translate};
189    use crate::{point, polygon};
190
191    /// Create a test which checks that we get `$should_be` when trying to find
192    /// the closest distance between `$p` and the line `(0, 0) -> (100, 100)`.
193    macro_rules! closest {
194        (intersects: $name:ident, $p:expr) => {
195            closest!($name, $p => Closest::Intersection($p.into()));
196        };
197        ($name:ident, $p:expr => $should_be:expr) => {
198            #[test]
199            fn $name() {
200                let line: Line<f32> = Line::from([(0., 0.), (100.0, 100.0)]);
201                let p: Point<f32> = $p.into();
202                let should_be: Closest<f32> = $should_be;
203
204                let got = line.closest_point(&p);
205                assert_eq!(got, should_be);
206            }
207        };
208    }
209
210    closest!(intersects: start_point, (0.0, 0.0));
211    closest!(intersects: end_point, (100.0, 100.0));
212    closest!(intersects: mid_point, (50.0, 50.0));
213    closest!(in_line_far_away, (1000.0, 1000.0) => Closest::SinglePoint(Point::new(100.0, 100.0)));
214    closest!(perpendicular_from_50_50, (0.0, 100.0) => Closest::SinglePoint(Point::new(50.0, 50.0)));
215
216    fn a_square(width: f32) -> LineString<f32> {
217        LineString::from(vec![
218            (0.0, 0.0),
219            (width, 0.0),
220            (width, width),
221            (0.0, width),
222            (0.0, 0.0),
223        ])
224    }
225
226    #[test]
227    fn zero_length_line_is_indeterminate() {
228        let line: Line<f32> = Line::from([(0.0, 0.0), (0.0, 0.0)]);
229        let p: Point<f32> = Point::new(100.0, 100.0);
230        let should_be: Closest<f32> = Closest::Indeterminate;
231
232        let got = line.closest_point(&p);
233        assert_eq!(got, should_be);
234    }
235
236    #[test]
237    fn line_string_with_single_element_behaves_like_line() {
238        let points = vec![(0.0, 0.0), (100.0, 100.0)];
239        let line_string = LineString::<f32>::from(points.clone());
240        let line = Line::new(points[0], points[1]);
241
242        let some_random_points = vec![
243            point!(x: 0.0, y: 0.0),
244            point!(x: 100.0, y: 100.0),
245            point!(x: 1000.0, y: 1000.0),
246            point!(x: 100.0, y: 0.0),
247            point!(x: 50.0, y: 50.0),
248            point!(x: 1234.567, y: -987.6543),
249        ];
250
251        for p in some_random_points {
252            assert_eq!(
253                line_string.closest_point(&p),
254                line.closest_point(&p),
255                "closest point to: {:?}",
256                p
257            );
258        }
259    }
260
261    #[test]
262    fn empty_line_string_is_indeterminate() {
263        let ls = LineString::new(Vec::new());
264        let p = Point::new(0.0, 0.0);
265
266        let got = ls.closest_point(&p);
267        assert_eq!(got, Closest::Indeterminate);
268    }
269
270    /// A polygon with 2 holes in it.
271    fn holy_polygon() -> Polygon<f32> {
272        let square: LineString<f32> = a_square(100.0);
273        let ring_1 = a_square(20.0).translate(10.0, 10.0);
274        let ring_2 = a_square(10.0).translate(70.0, 60.0);
275        Polygon::new(square, vec![ring_1, ring_2])
276    }
277
278    #[test]
279    fn polygon_without_rings_and_point_outside_is_same_as_linestring() {
280        let poly = holy_polygon();
281        let p = Point::new(1000.0, 12345.678);
282        assert!(
283            !poly.exterior().contains(&p),
284            "`p` should be outside the polygon!"
285        );
286
287        let poly_closest = poly.closest_point(&p);
288        let exterior_closest = poly.exterior().closest_point(&p);
289
290        assert_eq!(poly_closest, exterior_closest);
291    }
292
293    #[test]
294    fn polygon_with_point_on_interior_ring() {
295        let poly = holy_polygon();
296        let p = poly.interiors()[0].0[3];
297        let should_be = Closest::Intersection(p.into());
298
299        let got = poly.closest_point(&p.into());
300
301        assert_eq!(got, should_be);
302    }
303
304    #[test]
305    fn polygon_with_point_near_interior_ring() {
306        let poly = holy_polygon();
307        let p = point!(x: 17.0, y: 33.0);
308        assert!(poly.intersects(&p), "sanity check");
309
310        assert_eq!(Closest::Intersection(p), poly.closest_point(&p));
311    }
312
313    #[test]
314    fn polygon_with_interior_point() {
315        let square = polygon![
316            (x: 0.0, y: 0.0),
317            (x: 10.0, y: 0.0),
318            (x: 10.0, y: 10.0),
319            (x: 0.0, y: 10.0)
320        ];
321        let result = square.closest_point(&point!(x: 1.0, y: 2.0));
322
323        // the point is within the square, so the closest point should be the point itself.
324        assert_eq!(result, Closest::Intersection(point!(x: 1.0, y: 2.0)));
325    }
326
327    #[test]
328    fn multi_polygon_with_internal_and_external_points() {
329        use crate::{point, polygon};
330
331        let square_1 = polygon![
332            (x: 0.0, y: 0.0),
333            (x: 1.0, y: 0.0),
334            (x: 1.0, y: 1.0),
335            (x: 0.0, y: 1.0)
336        ];
337        use crate::Translate;
338        let square_10 = square_1.translate(10.0, 10.0);
339        let square_50 = square_1.translate(50.0, 50.0);
340
341        let multi_polygon = MultiPolygon::new(vec![square_1, square_10, square_50]);
342        let result = multi_polygon.closest_point(&point!(x: 8.0, y: 8.0));
343        assert_eq!(result, Closest::SinglePoint(point!(x: 10.0, y: 10.0)));
344
345        let result = multi_polygon.closest_point(&point!(x: 10.5, y: 10.5));
346        assert_eq!(result, Closest::Intersection(point!(x: 10.5, y: 10.5)));
347    }
348}