geo-repair 0.13.0

Fix invalid GIS geometries in parallel using rewritten GEOS/JTS repair module
Documentation 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121

use geo::{Area, BooleanOps, MultiPolygon, OpType, Polygon};

pub(crate) fn subtract_holes(shell: &Polygon<f64>, holes: &[Polygon<f64>]) -> Option<Polygon<f64>> {
    if holes.is_empty() {
        return Some(shell.clone());
    }
    if holes.len() == 1 {
        let diff = shell.boolean_op(&holes[0], OpType::Difference);
        return largest_polygon(diff);
    }

    // Batch: convert to MultiPolygon, single boolean difference
    // This processes all holes in one overlay operation (GEOS-style)
    let shell_mp = MultiPolygon::new(vec![shell.clone()]);
    let holes_mp = MultiPolygon::new(holes.to_vec());
    let diff = shell_mp.boolean_op(&holes_mp, OpType::Difference);
    if diff.0.is_empty() {
        return None;
    }
    largest_polygon(diff)
}

fn largest_polygon(mp: MultiPolygon<f64>) -> Option<Polygon<f64>> {
    mp.0.into_iter().max_by(|a, b| {
        a.unsigned_area()
            .partial_cmp(&b.unsigned_area())
            .unwrap_or(std::cmp::Ordering::Equal)
    })
}

#[cfg(test)]
mod tests {
    use super::*;
    use geo::{Coord, LineString};

    fn make_shell() -> Polygon<f64> {
        Polygon::new(
            LineString::new(vec![
                Coord { x: 0.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 },
                Coord { x: 0.0, y: 0.0 },
            ]),
            Vec::new(),
        )
    }

    fn make_hole() -> Polygon<f64> {
        Polygon::new(
            LineString::new(vec![
                Coord { x: 2.0, y: 2.0 },
                Coord { x: 8.0, y: 2.0 },
                Coord { x: 8.0, y: 8.0 },
                Coord { x: 2.0, y: 8.0 },
                Coord { x: 2.0, y: 2.0 },
            ]),
            Vec::new(),
        )
    }

    #[test]
    fn test_subtract_holes_empty() {
        let shell = make_shell();
        let result = subtract_holes(&shell, &[]);
        assert!(result.is_some());
    }

    #[test]
    fn test_subtract_holes_single() {
        let shell = make_shell();
        let hole = make_hole();
        let result = subtract_holes(&shell, &[hole]);
        assert!(result.is_some());
        let r = result.unwrap();
        assert_eq!(r.interiors().len(), 1);
    }

    #[test]
    fn test_subtract_holes_multiple() {
        let shell = make_shell();
        let hole1 = make_hole();
        let hole2 = Polygon::new(
            LineString::new(vec![
                Coord { x: 2.0, y: 2.0 },
                Coord { x: 8.0, y: 2.0 },
                Coord { x: 8.0, y: 8.0 },
                Coord { x: 2.0, y: 8.0 },
                Coord { x: 2.0, y: 2.0 },
            ]),
            Vec::new(),
        );
        let result = subtract_holes(&shell, &[hole1, hole2]);
        assert!(result.is_some());
    }

    #[test]
    fn test_subtract_holes_hole_equals_shell() {
        let shell = make_shell();
        let hole = make_shell();
        let result = subtract_holes(&shell, &[hole]);
        assert!(result.is_none());
    }

    #[test]
    fn test_subtract_holes_hole_outside_shell() {
        let shell = make_shell();
        let hole = Polygon::new(
            LineString::new(vec![
                Coord { x: 20.0, y: 20.0 },
                Coord { x: 25.0, y: 20.0 },
                Coord { x: 25.0, y: 25.0 },
                Coord { x: 20.0, y: 25.0 },
                Coord { x: 20.0, y: 20.0 },
            ]),
            Vec::new(),
        );
        let result = subtract_holes(&shell, &[hole]);
        assert!(result.is_some());
    }
}