geometry_kernel/
overlay.rs1use crate::canonicalize::{canonicalize_multi_polygon, canonicalize_polygon};
2use crate::error::{GeometryError, Result};
3use crate::precision::PrecisionModel;
4use crate::predicates::{is_convex_ring, is_ring_ccw, point_in_polygon, PointLocation};
5use crate::types::{BBox, Coord, LinearRing, MultiPolygon, Polygon};
6
7pub fn intersection(
8 subject: &MultiPolygon,
9 clip: &MultiPolygon,
10 precision: PrecisionModel,
11) -> Result<MultiPolygon> {
12 let mut out = Vec::new();
13
14 for subject_polygon in &subject.polygons {
15 for clip_polygon in &clip.polygons {
16 if let Some(polygon) = polygon_intersection(subject_polygon, clip_polygon, precision)? {
17 out.push(polygon);
18 }
19 }
20 }
21
22 Ok(canonicalize_multi_polygon(
23 &MultiPolygon::new(out),
24 precision,
25 ))
26}
27
28pub fn difference(
29 subject: &MultiPolygon,
30 clip: &MultiPolygon,
31 precision: PrecisionModel,
32) -> Result<MultiPolygon> {
33 let mut out = Vec::new();
34
35 for subject_polygon in &subject.polygons {
36 let subject_bbox = polygon_bbox(subject_polygon);
37 let mut changed = false;
38 let mut fully_covered = false;
39
40 for clip_polygon in &clip.polygons {
41 let Some(subject_bbox) = subject_bbox else {
42 continue;
43 };
44 let Some(clip_bbox) = polygon_bbox(clip_polygon) else {
45 continue;
46 };
47
48 if !subject_bbox.intersects(clip_bbox) {
49 continue;
50 }
51
52 changed = true;
53 if polygon_vertices_inside(subject_polygon, clip_polygon, precision) {
54 fully_covered = true;
55 break;
56 }
57 }
58
59 if !changed {
60 out.push(subject_polygon.clone());
61 } else if !fully_covered {
62 return Err(GeometryError::Unsupported(
63 "pure Rust difference currently supports disjoint or fully covered polygons"
64 .to_owned(),
65 ));
66 }
67 }
68
69 Ok(canonicalize_multi_polygon(
70 &MultiPolygon::new(out),
71 precision,
72 ))
73}
74
75fn polygon_intersection(
76 subject: &Polygon,
77 clip: &Polygon,
78 precision: PrecisionModel,
79) -> Result<Option<Polygon>> {
80 if subject.is_empty() || clip.is_empty() {
81 return Ok(None);
82 }
83 if !subject.holes.is_empty() || !clip.holes.is_empty() {
84 return Err(GeometryError::Unsupported(
85 "pure Rust intersection currently supports polygons without holes".to_owned(),
86 ));
87 }
88 if !is_convex_ring(&clip.exterior, precision) {
89 return Err(GeometryError::Unsupported(
90 "pure Rust intersection currently requires a convex clip polygon".to_owned(),
91 ));
92 }
93
94 let mut output = subject.exterior.coords[..subject.exterior.coords.len() - 1].to_vec();
95 let clip_coords = &clip.exterior.coords;
96 let clip_ccw = is_ring_ccw(&clip.exterior);
97
98 for edge in clip_coords.windows(2) {
99 if output.is_empty() {
100 break;
101 }
102
103 let input = output;
104 output = Vec::new();
105 let mut previous = *input.last().expect("non-empty input");
106 let mut previous_inside = inside_clip_edge(previous, edge[0], edge[1], clip_ccw, precision);
107
108 for current in input {
109 let current_inside = inside_clip_edge(current, edge[0], edge[1], clip_ccw, precision);
110 match (previous_inside, current_inside) {
111 (true, true) => output.push(current),
112 (true, false) => {
113 output.push(line_intersection(
114 previous, current, edge[0], edge[1], precision,
115 ));
116 }
117 (false, true) => {
118 output.push(line_intersection(
119 previous, current, edge[0], edge[1], precision,
120 ));
121 output.push(current);
122 }
123 (false, false) => {}
124 }
125 previous = current;
126 previous_inside = current_inside;
127 }
128 }
129
130 dedup_ring_coords(&mut output, precision);
131 if output.len() < 3 {
132 return Ok(None);
133 }
134
135 output.push(output[0]);
136 let polygon = canonicalize_polygon(
137 &Polygon::new(LinearRing::new(output), Vec::new()),
138 precision,
139 );
140 if polygon.is_empty() {
141 Ok(None)
142 } else {
143 Ok(Some(polygon))
144 }
145}
146
147fn inside_clip_edge(
148 point: Coord,
149 edge_start: Coord,
150 edge_end: Coord,
151 clip_ccw: bool,
152 precision: PrecisionModel,
153) -> bool {
154 let cross = (edge_end.x - edge_start.x) * (point.y - edge_start.y)
155 - (edge_end.y - edge_start.y) * (point.x - edge_start.x);
156 if clip_ccw {
157 cross >= -precision.epsilon()
158 } else {
159 cross <= precision.epsilon()
160 }
161}
162
163fn line_intersection(
164 a1: Coord,
165 a2: Coord,
166 b1: Coord,
167 b2: Coord,
168 precision: PrecisionModel,
169) -> Coord {
170 let denom = (a1.x - a2.x) * (b1.y - b2.y) - (a1.y - a2.y) * (b1.x - b2.x);
171 if denom.abs() <= precision.epsilon() {
172 return precision.snap_coord(a2);
173 }
174
175 let a_cross = a1.x * a2.y - a1.y * a2.x;
176 let b_cross = b1.x * b2.y - b1.y * b2.x;
177 precision.snap_coord(Coord::new(
178 (a_cross * (b1.x - b2.x) - (a1.x - a2.x) * b_cross) / denom,
179 (a_cross * (b1.y - b2.y) - (a1.y - a2.y) * b_cross) / denom,
180 ))
181}
182
183fn polygon_vertices_inside(subject: &Polygon, clip: &Polygon, precision: PrecisionModel) -> bool {
184 subject.exterior.coords[..subject.exterior.coords.len().saturating_sub(1)]
185 .iter()
186 .all(|point| {
187 matches!(
188 point_in_polygon(*point, clip, precision),
189 PointLocation::Interior | PointLocation::Boundary
190 )
191 })
192}
193
194fn polygon_bbox(polygon: &Polygon) -> Option<BBox> {
195 BBox::from_coords(&polygon.exterior.coords)
196}
197
198fn dedup_ring_coords(coords: &mut Vec<Coord>, precision: PrecisionModel) {
199 let mut deduped = Vec::with_capacity(coords.len());
200 for coord in coords.iter().copied() {
201 if deduped
202 .last()
203 .is_none_or(|last| !precision.same_coord(*last, coord))
204 {
205 deduped.push(coord);
206 }
207 }
208 if deduped.len() > 1 && precision.same_coord(deduped[0], *deduped.last().unwrap()) {
209 deduped.pop();
210 }
211 *coords = deduped;
212}