geo/algorithm/
line_locate_point.rs

1use crate::{
2    CoordFloat, Line, LineString, Point,
3    {euclidean_distance::EuclideanDistance, euclidean_length::EuclideanLength},
4};
5use std::ops::AddAssign;
6
7/// Returns a (option of the) fraction of the line's total length
8/// representing the location of the closest point on the line to
9/// the given point.
10///
11/// If the line has zero length the fraction returned is zero.
12///
13/// If either the point's coordinates or any coordinates of the line
14/// are not finite, returns `None`.
15///
16/// # Examples
17///
18/// ```
19/// use geo::{LineString, point};
20/// use geo::LineLocatePoint;
21///
22/// let linestring: LineString = vec![
23///     [-1.0, 0.0],
24///     [0.0, 0.0],
25///     [0.0, 1.0]
26/// ].into();
27///
28/// assert_eq!(linestring.line_locate_point(&point!(x: -1.0, y: 0.0)), Some(0.0));
29/// assert_eq!(linestring.line_locate_point(&point!(x: -0.5, y: 0.0)), Some(0.25));
30/// assert_eq!(linestring.line_locate_point(&point!(x: 0.0, y: 0.0)), Some(0.5));
31/// assert_eq!(linestring.line_locate_point(&point!(x: 0.0, y: 0.5)), Some(0.75));
32/// assert_eq!(linestring.line_locate_point(&point!(x: 0.0, y: 1.0)), Some(1.0));
33/// ```
34pub trait LineLocatePoint<T, Rhs> {
35    type Output;
36    type Rhs;
37
38    fn line_locate_point(&self, p: &Rhs) -> Self::Output;
39}
40
41impl<T> LineLocatePoint<T, Point<T>> for Line<T>
42where
43    T: CoordFloat,
44{
45    type Output = Option<T>;
46    type Rhs = Point<T>;
47
48    fn line_locate_point(&self, p: &Self::Rhs) -> Self::Output {
49        // let $s$ be the starting point of the line, and $v$ its
50        // direction vector. We want to find $l$ such that
51        // $(p - (s + lv)) \cdot v = 0$, i.e. the vector from
52        // $l$ along the line to $p$ is perpendicular to $v$.a
53
54        // vector $p - s$
55        let sp: Point<_> = *p - self.start_point();
56
57        // direction vector of line, $v$
58        let v: Point<_> = (self.end - self.start).into();
59
60        // $v \cdot v$
61        let v_sq = v.dot(v);
62        if v_sq == T::zero() {
63            // The line has zero length, return zero
64            Some(T::zero())
65        } else {
66            // $v \cdot (p - s)$
67            let v_dot_sp = v.dot(sp);
68            let l = v_dot_sp / v_sq;
69            if l.is_finite() {
70                Some(l.max(T::zero()).min(T::one()))
71            } else {
72                None
73            }
74        }
75    }
76}
77
78impl<T> LineLocatePoint<T, Point<T>> for LineString<T>
79where
80    T: CoordFloat + AddAssign,
81    Line<T>: EuclideanDistance<T, Point<T>> + EuclideanLength<T>,
82    LineString<T>: EuclideanLength<T>,
83{
84    type Output = Option<T>;
85    type Rhs = Point<T>;
86
87    fn line_locate_point(&self, p: &Self::Rhs) -> Self::Output {
88        let total_length = (*self).euclidean_length();
89        if total_length == T::zero() {
90            return Some(T::zero());
91        }
92        let mut cum_length = T::zero();
93        let mut closest_dist_to_point = T::infinity();
94        let mut fraction = T::zero();
95        for segment in self.lines() {
96            let segment_distance_to_point = segment.euclidean_distance(p);
97            let segment_length = segment.euclidean_length();
98            let segment_fraction = segment.line_locate_point(p)?; // if any segment has a None fraction, return None
99            if segment_distance_to_point < closest_dist_to_point {
100                closest_dist_to_point = segment_distance_to_point;
101                fraction = (cum_length + segment_fraction * segment_length) / total_length;
102            }
103            cum_length += segment_length;
104        }
105        Some(fraction)
106    }
107}
108
109#[cfg(test)]
110mod test {
111    use super::*;
112    use crate::geo_types::coord;
113    use crate::point;
114    use num_traits::Float;
115
116    #[test]
117    fn test_line_locate_point_line() {
118        // Some finite examples
119        let line = Line::new(coord! { x: -1.0, y: 0.0 }, coord! { x: 1.0, y: 0.0 });
120        let point = Point::new(0.0, 1.0);
121        assert_eq!(line.line_locate_point(&point), Some(0.5));
122
123        let point = Point::new(1.0, 1.0);
124        assert_eq!(line.line_locate_point(&point), Some(1.0));
125
126        let point = Point::new(2.0, 1.0);
127        assert_eq!(line.line_locate_point(&point), Some(1.0));
128
129        let point = Point::new(-1.0, 1.0);
130        assert_eq!(line.line_locate_point(&point), Some(0.0));
131
132        let point = Point::new(-2.0, 1.0);
133        assert_eq!(line.line_locate_point(&point), Some(0.0));
134
135        // point contains inf or nan
136        let point = Point::new(Float::nan(), 1.0);
137        assert_eq!(line.line_locate_point(&point), None);
138
139        let point = Point::new(Float::infinity(), 1.0);
140        assert_eq!(line.line_locate_point(&point), None);
141
142        let point = Point::new(Float::neg_infinity(), 1.0);
143        assert_eq!(line.line_locate_point(&point), None);
144
145        // line contains inf or nan
146        let line = Line::new(
147            coord! { x: 0.0, y: 0.0 },
148            coord! {
149                x: Float::infinity(),
150                y: 0.0,
151            },
152        );
153        let point = Point::new(1000.0, 1000.0);
154        assert_eq!(line.line_locate_point(&point), None);
155
156        let line = Line::new(
157            coord! { x: 0.0, y: 0.0 },
158            coord! {
159                x: Float::neg_infinity(),
160                y: 0.0,
161            },
162        );
163        let point = Point::new(1000.0, 1000.0);
164        assert_eq!(line.line_locate_point(&point), None);
165
166        let line = Line::new(
167            coord! { x: 0.0, y: 0.0 },
168            coord! {
169                x: Float::nan(),
170                y: 0.0,
171            },
172        );
173        let point = Point::new(1000.0, 1000.0);
174        assert_eq!(line.line_locate_point(&point), None);
175
176        // zero length line
177        let line: Line = Line::new(coord! { x: 1.0, y: 1.0 }, coord! { x: 1.0, y: 1.0 });
178        let pt = point!(x: 2.0, y: 2.0);
179        assert_eq!(line.line_locate_point(&pt), Some(0.0));
180
181        // another concrete example
182        let line: Line = Line::new(coord! { x: 0.0, y: 0.0 }, coord! { x: 10.0, y: 0.0 });
183        let pt = Point::new(555.0, 555.0);
184        assert_eq!(line.line_locate_point(&pt), Some(1.0));
185        let pt = Point::new(10.0000001, 0.0);
186        assert_eq!(line.line_locate_point(&pt), Some(1.0));
187        let pt = Point::new(9.0, 0.001);
188        assert_eq!(line.line_locate_point(&pt), Some(0.9));
189    }
190
191    #[test]
192    fn test_line_locate_point_linestring() {
193        // finite example using the ring
194        let ring: LineString = geo_test_fixtures::ring::<f64>();
195        let pt = point!(x: 10.0, y: 1.0);
196        assert_eq!(ring.line_locate_point(&pt), Some(0.0));
197
198        let pt = point!(x: 10.0, y: 1.0000000000000742);
199        assert_eq!(ring.line_locate_point(&pt), Some(0.9999999999999988));
200
201        let pt = point!(x: 10.0, y: 1.0);
202        assert_eq!(ring.line_locate_point(&pt), Some(0.0));
203
204        // point contains inf or nan
205        let pt = point!(x: Float::nan(), y: 1.0);
206        assert_eq!(ring.line_locate_point(&pt), None);
207
208        let pt = point!(x: Float::infinity(), y: 1.0);
209        assert_eq!(ring.line_locate_point(&pt), None);
210
211        let pt = point!(x: Float::neg_infinity(), y: 1.0);
212        assert_eq!(ring.line_locate_point(&pt), None);
213
214        // point is equidistant to two line segments - return the fraction from the first closest
215        let line: LineString = LineString::new(vec![
216            (0.0, 0.0).into(),
217            (1.0, 0.0).into(),
218            (1.0, 1.0).into(),
219            (0.0, 1.0).into(),
220        ]);
221        let pt = point!(x: 0.0, y: 0.5);
222        assert_eq!(line.line_locate_point(&pt), Some(0.0));
223
224        let line: LineString = LineString::new(vec![
225            (1.0, 1.0).into(),
226            (1.0, 1.0).into(),
227            (1.0, 1.0).into(),
228        ]);
229        let pt = point!(x: 2.0, y: 2.0);
230        assert_eq!(line.line_locate_point(&pt), Some(0.0));
231
232        // line contains inf or nan
233        let line: LineString = LineString::new(vec![
234            coord! { x: 1.0, y: 1.0 },
235            coord! {
236                x: Float::nan(),
237                y: 1.0,
238            },
239            coord! { x: 0.0, y: 0.0 },
240        ]);
241        let pt = point!(x: 2.0, y: 2.0);
242        assert_eq!(line.line_locate_point(&pt), None);
243
244        let line: LineString = LineString::new(vec![
245            coord! { x: 1.0, y: 1.0 },
246            coord! {
247                x: Float::infinity(),
248                y: 1.0,
249            },
250            coord! { x: 0.0, y: 0.0 },
251        ]);
252        let pt = point!(x: 2.0, y: 2.0);
253        assert_eq!(line.line_locate_point(&pt), None);
254        let line: LineString = LineString::new(vec![
255            coord! { x: 1.0, y: 1.0 },
256            coord! {
257                x: Float::neg_infinity(),
258                y: 1.0,
259            },
260            coord! { x: 0.0, y: 0.0 },
261        ]);
262        let pt = point!(x: 2.0, y: 2.0);
263        assert_eq!(line.line_locate_point(&pt), None);
264    }
265}
geo/algorithm/line_locate_point.rs

geo/algorithm/
line_locate_point.rs
