osm_pbf2json/
geo.rs

1use geo::prelude::*;
2use geo::Closest;
3use geo_types::{Coordinate, Geometry, Line, LineString, MultiPoint, MultiPolygon, Point, Polygon};
4use serde::{Deserialize, Serialize};
5use std::convert::{TryFrom, TryInto};
6
7const EQ_PRECISION: f64 = 1.0e-5;
8
9#[derive(Serialize, Deserialize, Debug)]
10pub struct Location {
11    pub lat: f64,
12    pub lon: f64,
13}
14
15impl From<(f64, f64)> for Location {
16    fn from(tuple: (f64, f64)) -> Location {
17        Location {
18            lon: tuple.0,
19            lat: tuple.1,
20        }
21    }
22}
23
24#[derive(Clone, Debug)]
25pub struct SegmentGeometry {
26    bounding_box: BoundingBox,
27    line_string: LineString<f64>,
28}
29
30#[derive(Clone, Debug)]
31pub struct BoundaryGeometry {
32    bounding_box: BoundingBox,
33    multi_polygon: MultiPolygon<f64>,
34}
35
36impl BoundaryGeometry {
37    pub fn new(multi_polygon: MultiPolygon<f64>) -> Result<Self, &'static str> {
38        let bounding_box = (&multi_polygon).try_into()?;
39        Ok(BoundaryGeometry {
40            bounding_box,
41            multi_polygon,
42        })
43    }
44
45    pub fn coordinates(&self) -> Vec<Vec<Vec<(f64, f64)>>> {
46        self.multi_polygon
47            .clone()
48            .into_iter()
49            .map(|polygon| {
50                let (exterior, interiours) = polygon.into_inner();
51                let mut rings = vec![exterior];
52                rings.extend(interiours);
53                rings
54            })
55            .map(|line_strings| {
56                line_strings
57                    .iter()
58                    .map(|ls| ls.points_iter().map(|p| (p.x(), p.y())).collect())
59                    .collect()
60            })
61            .collect()
62    }
63
64    pub fn sw_ne(&self) -> ([f64; 2], [f64; 2]) {
65        (self.bounding_box.sw, self.bounding_box.ne)
66    }
67
68    pub fn intersects(&self, geometry: &SegmentGeometry) -> bool {
69        self.multi_polygon
70            .0
71            .iter()
72            .any(|polygon| polygon.intersects(&geometry.line_string))
73    }
74
75    pub fn owns(&self, geometry: &SegmentGeometry) -> bool {
76        if let Some(centroid) = geometry.line_string.centroid() {
77            self.multi_polygon.contains(&centroid)
78        } else {
79            false
80        }
81    }
82}
83
84#[derive(Clone, Debug)]
85struct BoundingBox {
86    sw: [f64; 2],
87    ne: [f64; 2],
88}
89
90impl PartialEq<Location> for Location {
91    fn eq(&self, other: &Self) -> bool {
92        let self_point = Point::new(self.lon, self.lat);
93        let other_point = Point::new(other.lon, other.lat);
94        let distance = self_point.euclidean_distance(&other_point);
95        distance < EQ_PRECISION
96    }
97}
98
99impl PartialEq<Bounds> for Bounds {
100    fn eq(&self, other: &Self) -> bool {
101        let (self_ne, self_sw) = self.into();
102        let (other_ne, other_sw) = other.into();
103        self_ne == other_ne && self_sw == other_sw
104    }
105}
106
107impl BoundingBox {
108    pub fn pad(&self, distance: f64) -> BoundingBox {
109        let sw: Point<f64> = self.sw.into();
110        let ne: Point<f64> = self.ne.into();
111        let padding: Point<f64> = (distance, distance).into();
112        let sw_padded = sw - padding;
113        let ne_padded = ne + padding;
114        BoundingBox {
115            sw: [sw_padded.lng(), sw_padded.lat()],
116            ne: [ne_padded.lng(), ne_padded.lat()],
117        }
118    }
119}
120
121impl TryFrom<&MultiPolygon<f64>> for BoundingBox {
122    type Error = &'static str;
123
124    fn try_from(multi_polygon: &MultiPolygon<f64>) -> Result<Self, Self::Error> {
125        let rect = multi_polygon
126            .bounding_rect()
127            .ok_or("cannot get bounding box for the given set of coordinates")?;
128        let sw = [rect.min().x, rect.min().y];
129        let ne = [rect.max().x, rect.max().y];
130        Ok(BoundingBox { sw, ne })
131    }
132}
133
134impl TryFrom<&LineString<f64>> for BoundingBox {
135    type Error = &'static str;
136
137    fn try_from(line_string: &LineString<f64>) -> Result<Self, Self::Error> {
138        let rect = line_string
139            .bounding_rect()
140            .ok_or("cannot get bounding box for the given set of coordinates")?;
141        let sw = [rect.min().x, rect.min().y];
142        let ne = [rect.max().x, rect.max().y];
143        Ok(BoundingBox { sw, ne })
144    }
145}
146
147impl SegmentGeometry {
148    pub fn new(coordinates: Vec<(f64, f64)>) -> Result<Self, &'static str> {
149        let line_string: LineString<f64> = coordinates.into();
150        let bounding_box: BoundingBox = (&line_string).try_into()?;
151        let geometry = SegmentGeometry {
152            bounding_box,
153            line_string,
154        };
155        Ok(geometry)
156    }
157
158    pub fn len(&self) -> usize {
159        self.line_string.points_iter().count()
160    }
161
162    pub fn sw_ne(&self) -> ([f64; 2], [f64; 2]) {
163        (self.bounding_box.sw, self.bounding_box.ne)
164    }
165
166    pub fn padded_sw_ne(&self, distance: f64) -> ([f64; 2], [f64; 2]) {
167        let BoundingBox { sw, ne } = self.bounding_box.pad(distance);
168        (sw, ne)
169    }
170}
171
172pub trait Length {
173    fn length(&self) -> f64;
174}
175
176impl Length for SegmentGeometry {
177    fn length(&self) -> f64 {
178        let sw: Coordinate<f64> = self.bounding_box.sw.into();
179        let ne: Coordinate<f64> = self.bounding_box.ne.into();
180        let line = Line::new(sw, ne);
181        line.euclidean_length()
182    }
183}
184
185impl Length for Vec<&SegmentGeometry> {
186    fn length(&self) -> f64 {
187        self.iter().map(|segment| segment.length()).sum()
188    }
189}
190
191impl From<&SegmentGeometry> for Vec<(f64, f64)> {
192    fn from(geometry: &SegmentGeometry) -> Vec<(f64, f64)> {
193        geometry
194            .line_string
195            .points_iter()
196            .map(|c| (c.x(), c.y()))
197            .collect()
198    }
199}
200
201impl From<SegmentGeometry> for Vec<(f64, f64)> {
202    fn from(geometry: SegmentGeometry) -> Vec<(f64, f64)> {
203        geometry
204            .line_string
205            .points_iter()
206            .map(|c| (c.x(), c.y()))
207            .collect()
208    }
209}
210
211impl From<&(f64, f64)> for Location {
212    fn from(coordinates: &(f64, f64)) -> Self {
213        Location {
214            lon: coordinates.0,
215            lat: coordinates.1,
216        }
217    }
218}
219
220impl From<Location> for [f64; 2] {
221    fn from(loc: Location) -> Self {
222        [loc.lon, loc.lat]
223    }
224}
225
226impl From<Point<f64>> for Location {
227    fn from(point: Point<f64>) -> Self {
228        Location {
229            lat: point.lat(),
230            lon: point.lng(),
231        }
232    }
233}
234
235impl From<Coordinate<f64>> for Location {
236    fn from(coordinate: Coordinate<f64>) -> Self {
237        Location {
238            lat: coordinate.y,
239            lon: coordinate.x,
240        }
241    }
242}
243
244#[derive(Serialize, Deserialize, Debug)]
245pub struct Bounds {
246    e: f64,
247    n: f64,
248    s: f64,
249    w: f64,
250}
251
252pub trait Midpoint {
253    fn midpoint(&self) -> Option<(f64, f64)>;
254}
255
256impl Midpoint for Vec<&SegmentGeometry> {
257    fn midpoint(&self) -> Option<(f64, f64)> {
258        let flattened: Vec<_> = self
259            .iter()
260            .flat_map(|geometry| {
261                let coordinates: Vec<(f64, f64)> = (*geometry).into();
262                coordinates
263            })
264            .collect();
265        let multi_points: MultiPoint<f64> = flattened.into();
266        let centroid = multi_points.centroid()?;
267        let closest = multi_points.closest_point(&centroid);
268        match closest {
269            Closest::Intersection(p) => Some((p.lng(), p.lat())),
270            Closest::SinglePoint(p) => Some((p.lng(), p.lat())),
271            _ => None,
272        }
273    }
274}
275
276impl From<&Bounds> for (Location, Location) {
277    fn from(bounds: &Bounds) -> Self {
278        let ne = Location {
279            lon: bounds.n,
280            lat: bounds.e,
281        };
282        let sw = Location {
283            lon: bounds.s,
284            lat: bounds.w,
285        };
286
287        (ne, sw)
288    }
289}
290
291fn get_geometry(coordinates: &[(f64, f64)]) -> Option<Geometry<f64>> {
292    let line_string: LineString<f64> = coordinates.to_vec().into();
293    let first = line_string.points_iter().next()?;
294    let last = line_string.points_iter().last()?;
295    if first == last {
296        let polygon = Polygon::new(line_string, vec![]);
297        Some(Geometry::Polygon(polygon))
298    } else {
299        Some(Geometry::LineString(line_string))
300    }
301}
302
303fn get_bounds(geometry: &Geometry<f64>) -> Option<Bounds> {
304    let rect = match geometry {
305        Geometry::LineString(ls) => ls.bounding_rect(),
306        Geometry::Polygon(p) => p.bounding_rect(),
307        _ => None,
308    }?;
309    Some(Bounds {
310        e: rect.max().x,
311        n: rect.max().y,
312        s: rect.min().y,
313        w: rect.min().x,
314    })
315}
316
317pub trait Centerable {
318    fn get_centroid(&self) -> Option<Location>;
319}
320
321impl Centerable for Vec<(f64, f64)> {
322    fn get_centroid(&self) -> Option<Location> {
323        let geometry = get_geometry(self)?;
324        geometry.get_centroid()
325    }
326}
327
328impl Centerable for Geometry<f64> {
329    fn get_centroid(&self) -> Option<Location> {
330        let point = match self {
331            Geometry::LineString(ls) => ls.centroid(),
332            Geometry::Polygon(p) => p.centroid(),
333            _ => None,
334        }?;
335        Some(point.into())
336    }
337}
338
339pub fn get_geo_info(coordinates: &[(f64, f64)]) -> (Option<Location>, Option<Bounds>) {
340    if let Some(geo) = get_geometry(coordinates) {
341        let centroid = geo.get_centroid();
342        let bounds = get_bounds(&geo);
343        return (centroid, bounds);
344    }
345    (None, None)
346}
347
348pub fn get_compound_coordinates(coordinates: Vec<(f64, f64)>) -> Vec<(f64, f64)> {
349    let multi_points: MultiPoint<_> = coordinates.into();
350    let convex_hull = multi_points.convex_hull();
351    convex_hull
352        .exterior()
353        .points_iter()
354        .map(|p| (p.lng(), p.lat()))
355        .collect()
356}
357
358#[cfg(test)]
359mod tests {
360    use super::*;
361    use approx::*;
362
363    fn approx_eq<T: Into<[f64; 2]>>(a: [f64; 2], o: Option<T>) {
364        let b: [f64; 2] = o.unwrap().into();
365        assert_relative_eq!(a[0], b[0], epsilon = f64::EPSILON);
366        assert_relative_eq!(a[1], b[1], epsilon = f64::EPSILON);
367    }
368
369    #[test]
370    fn get_centroid_for_line() {
371        let coordinates = vec![(9., 50.), (9., 51.), (10., 51.)];
372        // 1     2
373        //  c
374        //
375        // 0
376        approx_eq([9.25, 50.75], coordinates.get_centroid());
377    }
378
379    #[test]
380    fn midpoint() {
381        let coordinates = vec![(9., 50.), (9., 51.), (10., 51.)];
382        let geometry_1 = SegmentGeometry::new(coordinates).unwrap();
383        let coordinates = vec![(12., 51.), (12., 50.)];
384        let geometry_2 = SegmentGeometry::new(coordinates).unwrap();
385        // 1.1   1.2        2.0
386        //
387        //
388        // 1.0              2.1
389        let midpoint = vec![&geometry_1, &geometry_2]
390            .midpoint()
391            .map(|(lng, lat)| [lng, lat]);
392        approx_eq([10., 51.], midpoint);
393    }
394
395    #[test]
396    fn get_geo_info_open() {
397        let coordinates = vec![(5., 49.), (6., 50.), (7., 49.)];
398        let (centroid, bounds) = get_geo_info(&coordinates);
399        let reference_loc = Location { lat: 49.5, lon: 6. };
400        assert_eq!(centroid.unwrap(), reference_loc);
401        let reference_bounds = Bounds {
402            e: 7.,
403            n: 50.,
404            s: 49.,
405            w: 5.,
406        };
407        assert_eq!(bounds.unwrap(), reference_bounds);
408    }
409
410    #[test]
411    fn get_geo_info_closed() {
412        let coordinates = vec![(5., 49.), (6., 50.), (7., 49.), (5., 49.)];
413        let (centroid, bounds) = get_geo_info(&coordinates);
414        let reference_loc = Location {
415            lat: 49.333_333,
416            lon: 6.,
417        };
418        assert_eq!(centroid.unwrap(), reference_loc);
419        let reference_bounds = Bounds {
420            e: 7.,
421            n: 50.,
422            s: 49.,
423            w: 5.,
424        };
425        assert_eq!(bounds.unwrap(), reference_bounds);
426    }
427}