geo/algorithm/relate/geomgraph/intersection_matrix.rs
1use crate::{
2 ContainsProperly, GeoNum, GeometryCow, HasDimensions, coordinate_position::CoordPos,
3 dimensions::Dimensions,
4};
5
6use crate::geometry_cow::GeometryCow::Point;
7use crate::relate::geomgraph::intersection_matrix::dimension_matcher::DimensionMatcher;
8use std::str::FromStr;
9
10/// Output from [`Relate::relate`](trait.Relate.html) which models a *Dimensionally Extended Nine-Intersection Model (DE-9IM)* matrix.
11///
12/// Represented as 3x3 matrices, they express the topological relationships between the Interior, Boundary, and Exterior of two geometries.
13///
14/// Consider these partially overlapping geometries:
15/// ```rust
16/// # use geo::wkt;
17/// #
18/// let a = wkt!(POLYGON((0. 0.,0. 4.,4. 4.,4. 0., 0. 0.)));
19/// let b = wkt!(LINESTRING(2. 2.,6. 2.));
20/// ```
21///
22/// <svg xmlns="http://www.w3.org/2000/svg" viewBox="-1 -1 8 6" width="100%" height="200" font-size="1px">
23/// <!-- blue square -->
24/// <rect x="0" y="0" width="4" height="4" fill="#8888FF" />
25/// <text x="0.25" y="1">A</text>
26/// <text x="5" y="1.8">B</text>
27/// <!-- yellow line sticking into the square -->
28/// <line x1="2" y1="2" x2="6" y2="2" stroke="#EEEE00" stroke-width="0.2" />
29/// </svg>
30///
31/// Their intersection matrix looks like this:
32///
33/// ```rust
34/// # use geo::wkt;
35/// # use geo::relate::{Relate, IntersectionMatrix};
36/// # use std::str::FromStr;
37/// #
38/// # let a = wkt!(POLYGON((0. 0.,0. 4.,4. 4.,4. 0., 0. 0.)));
39/// # let b = wkt!(LINESTRING(2. 2.,6. 2.));
40/// let intersection_matrix = a.relate(&b);
41/// ```
42///
43/// | I.M. | Interior B | Boundary B | Exterior B |
44/// |----------------|------------|------------|------------|
45/// | **Interior A** | 1 (line) | 0 (point) | 2 (area) |
46/// | **Boundary A** | 0 (point) | empty | 1 (line) |
47/// | **Exterior A** | 1 (line) | 0 (point) | 2 (area) |
48///
49/// The indices of the matrix are expressed as the enum cases [CoordPos::Inside, CoordPos::OnBoundary, CoordPos::Outside](CoordPos).
50/// The entries represent the [Dimensions](enum.Dimension.html) of the intersection (if any) at those positions.
51///
52/// # Examples
53///
54/// <svg xmlns="http://www.w3.org/2000/svg" viewBox="-1 -1 8 6" width="100%" height="200" font-size="1px">
55/// <!-- blue square -->
56/// <rect x="0" y="0" width="4" height="4" fill="#8888FF" />
57/// <text x="0.25" y="1">A</text>
58/// <text x="5" y="1.8">B</text>
59/// <!-- yellow line sticking into the square -->
60/// <line x1="2" y1="2" x2="6" y2="2" stroke="#FFFF11" stroke-width="0.2" />
61///
62/// <!-- green line showing overlap -->
63/// <line x1="2" y1="2" x2="4" y2="2" stroke="#006600" stroke-width="0.2" />
64/// </svg>
65///
66/// ```rust
67/// # use geo::wkt;
68/// # use geo::relate::{Relate, IntersectionMatrix};
69/// # use geo::coordinate_position::CoordPos;
70/// # use geo::dimensions::Dimensions;
71/// #
72/// # let a = wkt!(POLYGON((0. 0.,0. 4.,4. 4.,4. 0., 0. 0.)));
73/// # let b = wkt!(LINESTRING(2. 2.,6. 2.));
74/// # let intersection_matrix = a.relate(&b);
75/// let intersection = intersection_matrix.get(
76/// // The interior of A...
77/// CoordPos::Inside,
78/// // ...and the interior of B...
79/// CoordPos::Inside
80/// );
81/// // ...intersect along a line.
82/// assert_eq!(intersection, Dimensions::OneDimensional); // lines have one dimension
83/// ```
84///
85/// <svg xmlns="http://www.w3.org/2000/svg" viewBox="-1 -1 8 6" width="100%" height="200" font-size="1px">
86/// <!-- blue square -->
87/// <rect x="0" y="0" width="4" height="4" fill="#8888FF" />
88/// <text x="0.25" y="1">A</text>
89/// <text x="5" y="1.8">B</text>
90/// <!-- yellow line sticking into the square -->
91/// <line x1="2" y1="2" x2="6" y2="2" stroke="#EEEE00" stroke-width="0.2" />
92///
93/// <!-- green dot showing overlap -->
94/// <circle cx="4" cy="2" r="0.3" stroke-width="0" fill="#006600" />
95/// </svg>
96///
97/// ```rust
98/// # use geo::wkt;
99/// # use geo::relate::{Relate, IntersectionMatrix};
100/// # use geo::coordinate_position::CoordPos;
101/// # use geo::dimensions::Dimensions;
102/// #
103/// # let a = wkt!(POLYGON((0. 0.,0. 4.,4. 4.,4. 0., 0. 0.)));
104/// # let b = wkt!(LINESTRING(2. 2.,6. 2.));
105/// # let intersection_matrix = a.relate(&b);
106/// let intersection = intersection_matrix.get(
107/// // The boundary of A...
108/// CoordPos::OnBoundary,
109/// // ...and the interior of B...
110/// CoordPos::Inside
111/// );
112/// // ...intersect at a point.
113/// assert_eq!(intersection, Dimensions::ZeroDimensional); // points have zero dimension
114/// ```
115///
116/// <svg xmlns="http://www.w3.org/2000/svg" viewBox="-1 -1 8 6" width="100%" height="200" font-size="1px">
117/// <!-- blue square -->
118/// <rect x="0" y="0" width="4" height="4" fill="#8888FF" stroke-width="0.1" stroke="black" stroke-dasharray="0.15" />
119/// <text x="0.25" y="1">A</text>
120/// <text x="5" y="1.8">B</text>
121/// <!-- yellow line sticking into the square -->
122/// <line x1="2" y1="2" x2="6" y2="2" stroke="#EEEE00" stroke-width="0.2" />
123/// <circle cx="2" cy="2" r="0.3" stroke-width="0.1" stroke="black" stroke-dasharray="0.15" fill="none" />
124/// <circle cx="6" cy="2" r="0.3" stroke-width="0.1" stroke="black" stroke-dasharray="0.15" fill="none" />
125/// </svg>
126///
127/// ```rust
128/// # use geo::wkt;
129/// # use geo::relate::{Relate, IntersectionMatrix};
130/// # use geo::coordinate_position::CoordPos;
131/// # use geo::dimensions::Dimensions;
132/// #
133/// # let a = wkt!(POLYGON((0. 0.,0. 4.,4. 4.,4. 0., 0. 0.)));
134/// # let b = wkt!(LINESTRING(2. 2.,6. 2.));
135/// # let intersection_matrix = a.relate(&b);
136/// let intersection = intersection_matrix.get(
137/// // The boundary of A...
138/// CoordPos::OnBoundary,
139/// // ...and the boundary of B...
140/// CoordPos::OnBoundary
141/// );
142/// // ...do not intersect.
143/// assert_eq!(intersection, Dimensions::Empty)
144/// ```
145///
146/// ## Predicates
147///
148/// Computing an `IntersectionMatrix` can be expensive, but once you have it, checking predicates is essentially free.
149///
150/// ```rust
151/// use geo::wkt;
152/// use geo::relate::Relate;
153///
154/// let a = wkt!(POLYGON((0. 0.,0. 4.,4. 4.,4. 0., 0. 0.)));
155/// let b = wkt!(LINESTRING(2. 2.,6. 2.));
156/// // potentially expensive
157/// let intersection_matrix = a.relate(&b);
158/// // fast!
159/// assert!(intersection_matrix.is_intersects());
160/// // fast!
161/// assert!(intersection_matrix.is_crosses());
162/// // fast!
163/// assert!(!intersection_matrix.is_contains());
164/// ```
165///
166/// ## String representation
167///
168/// DE-9IM matrix values (such as `1020F1102`) are a concise way to specify the topological relationship between two geometries.
169///
170/// ```rust
171/// use geo::wkt;
172/// use geo::relate::{Relate, IntersectionMatrix};
173/// use std::str::FromStr;
174///
175/// let a = wkt!(POLYGON((0. 0.,0. 4.,4. 4.,4. 0., 0. 0.)));
176/// let b = wkt!(LINESTRING(2. 2.,6. 2.));
177/// let intersection_matrix = a.relate(&b);
178/// let expected = IntersectionMatrix::from_str("1020F1102").expect("valid DE-9IM specification");
179/// assert_eq!(intersection_matrix, expected)
180/// ```
181///
182/// # References
183///
184/// For a description of the DE-9IM and the spatial predicates derived from it,
185/// see the following references:
186/// - [OGC 99-049 OpenGIS Simple Features Specification for SQL](http://portal.opengeospatial.org/files/?artifact_id=829), Section 2.1.13
187/// - [OGC 06-103r4 OpenGIS Implementation Standard for Geographic information - Simple feature access - Part 1: Common architecture](http://portal.opengeospatial.org/files/?artifact_id=25355), Section 6.1.15 (which provides some further details on certain predicate specifications).
188/// - Wikipedia article on [DE-9IM](https://en.wikipedia.org/wiki/DE-9IM)
189///
190/// This implementation is heavily based on that from the [JTS project](https://github.com/locationtech/jts/blob/master/modules/core/src/main/java/org/locationtech/jts/geom/IntersectionMatrix.java).
191#[derive(PartialEq, Eq, Clone)]
192pub struct IntersectionMatrix(LocationArray<LocationArray<Dimensions>>);
193
194/// Helper struct so we can index IntersectionMatrix by CoordPos
195///
196/// CoordPos enum members are ordered: OnBoundary, Inside, Outside
197/// DE-9IM matrices are ordered: Inside, Boundary, Exterior
198///
199/// So we can't simply `CoordPos as usize` without losing the conventional ordering
200/// of elements, which is useful for debug / interop.
201#[derive(PartialEq, Eq, Clone, Copy)]
202struct LocationArray<T>([T; 3]);
203
204impl<T> LocationArray<T> {
205 fn iter(&self) -> impl Iterator<Item = &T> {
206 self.0.iter()
207 }
208}
209
210impl<T> std::ops::Index<CoordPos> for LocationArray<T> {
211 type Output = T;
212
213 fn index(&self, index: CoordPos) -> &Self::Output {
214 match index {
215 CoordPos::Inside => &self.0[0],
216 CoordPos::OnBoundary => &self.0[1],
217 CoordPos::Outside => &self.0[2],
218 }
219 }
220}
221
222impl<T> std::ops::IndexMut<CoordPos> for LocationArray<T> {
223 fn index_mut(&mut self, index: CoordPos) -> &mut Self::Output {
224 match index {
225 CoordPos::Inside => &mut self.0[0],
226 CoordPos::OnBoundary => &mut self.0[1],
227 CoordPos::Outside => &mut self.0[2],
228 }
229 }
230}
231
232#[derive(Debug)]
233pub struct InvalidInputError {
234 message: String,
235}
236
237impl InvalidInputError {
238 fn new(message: String) -> Self {
239 Self { message }
240 }
241}
242
243impl std::error::Error for InvalidInputError {}
244impl std::fmt::Display for InvalidInputError {
245 fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
246 write!(f, "invalid input: {}", self.message)
247 }
248}
249
250impl std::fmt::Debug for IntersectionMatrix {
251 fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
252 fn char_for_dim(dim: &Dimensions) -> &'static str {
253 match dim {
254 Dimensions::Empty => "F",
255 Dimensions::ZeroDimensional => "0",
256 Dimensions::OneDimensional => "1",
257 Dimensions::TwoDimensional => "2",
258 }
259 }
260 let text = self
261 .0
262 .iter()
263 .flat_map(|r| r.iter().map(char_for_dim))
264 .collect::<Vec<&str>>()
265 .join("");
266
267 write!(f, "IntersectionMatrix({})", &text)
268 }
269}
270
271impl IntersectionMatrix {
272 pub const fn empty() -> Self {
273 IntersectionMatrix(LocationArray([LocationArray([Dimensions::Empty; 3]); 3]))
274 }
275
276 pub(crate) const fn empty_disjoint() -> Self {
277 IntersectionMatrix(LocationArray([
278 LocationArray([Dimensions::Empty, Dimensions::Empty, Dimensions::Empty]),
279 LocationArray([Dimensions::Empty, Dimensions::Empty, Dimensions::Empty]),
280 // since Geometries are finite and embedded in a 2-D space,
281 // the `(Outside, Outside)` element must always be 2-D
282 LocationArray([
283 Dimensions::Empty,
284 Dimensions::Empty,
285 Dimensions::TwoDimensional,
286 ]),
287 ]))
288 }
289
290 /// If the Geometries are disjoint, we need to enter their dimension and boundary dimension in
291 /// the `Outside` rows in the IM
292 pub(crate) fn compute_disjoint(
293 &mut self,
294 geometry_a: &(impl HasDimensions + ?Sized),
295 geometry_b: &(impl HasDimensions + ?Sized),
296 ) {
297 {
298 let dimensions = geometry_a.dimensions();
299 if dimensions != Dimensions::Empty {
300 self.set(CoordPos::Inside, CoordPos::Outside, dimensions);
301
302 let boundary_dimensions = geometry_a.boundary_dimensions();
303 if boundary_dimensions != Dimensions::Empty {
304 self.set(CoordPos::OnBoundary, CoordPos::Outside, boundary_dimensions);
305 }
306 }
307 }
308
309 {
310 let dimensions = geometry_b.dimensions();
311 if dimensions != Dimensions::Empty {
312 self.set(CoordPos::Outside, CoordPos::Inside, dimensions);
313
314 let boundary_dimensions = geometry_b.boundary_dimensions();
315 if boundary_dimensions != Dimensions::Empty {
316 self.set(CoordPos::Outside, CoordPos::OnBoundary, boundary_dimensions);
317 }
318 }
319 }
320 }
321
322 /// Set `dimensions` of the cell specified by the positions.
323 ///
324 /// `position_a`: which position `dimensions` applies to within the first geometry
325 /// `position_b`: which position `dimensions` applies to within the second geometry
326 /// `dimensions`: the dimension of the incident
327 pub(crate) fn set(
328 &mut self,
329 position_a: CoordPos,
330 position_b: CoordPos,
331 dimensions: Dimensions,
332 ) {
333 self.0[position_a][position_b] = dimensions;
334 }
335
336 /// Reports an incident of `dimensions`, which updates the IntersectionMatrix if it's greater
337 /// than what has been reported so far.
338 ///
339 /// `position_a`: which position `minimum_dimensions` applies to within the first geometry
340 /// `position_b`: which position `minimum_dimensions` applies to within the second geometry
341 /// `minimum_dimensions`: the dimension of the incident
342 pub(crate) fn set_at_least(
343 &mut self,
344 position_a: CoordPos,
345 position_b: CoordPos,
346 minimum_dimensions: Dimensions,
347 ) {
348 if self.0[position_a][position_b] < minimum_dimensions {
349 self.0[position_a][position_b] = minimum_dimensions;
350 }
351 }
352
353 /// If both geometries have `Some` position, then changes the specified element to at
354 /// least `minimum_dimensions`.
355 ///
356 /// Else, if either is none, do nothing.
357 ///
358 /// `position_a`: which position `minimum_dimensions` applies to within the first geometry, or
359 /// `None` if the dimension was not incident with the first geometry.
360 /// `position_b`: which position `minimum_dimensions` applies to within the second geometry, or
361 /// `None` if the dimension was not incident with the second geometry.
362 /// `minimum_dimensions`: the dimension of the incident
363 pub(crate) fn set_at_least_if_in_both(
364 &mut self,
365 position_a: Option<CoordPos>,
366 position_b: Option<CoordPos>,
367 minimum_dimensions: Dimensions,
368 ) {
369 if let (Some(position_a), Some(position_b)) = (position_a, position_b) {
370 self.set_at_least(position_a, position_b, minimum_dimensions);
371 }
372 }
373
374 pub(crate) fn set_at_least_from_string(
375 &mut self,
376 dimensions: &str,
377 ) -> Result<(), InvalidInputError> {
378 if dimensions.len() != 9 {
379 let message = format!("Expected dimensions length 9, found: {}", dimensions.len());
380 return Err(InvalidInputError::new(message));
381 }
382
383 let mut chars = dimensions.chars();
384 for a in &[CoordPos::Inside, CoordPos::OnBoundary, CoordPos::Outside] {
385 for b in &[CoordPos::Inside, CoordPos::OnBoundary, CoordPos::Outside] {
386 match chars.next().expect("already validated length is 9") {
387 '0' => self.0[*a][*b] = self.0[*a][*b].max(Dimensions::ZeroDimensional),
388 '1' => self.0[*a][*b] = self.0[*a][*b].max(Dimensions::OneDimensional),
389 '2' => self.0[*a][*b] = self.0[*a][*b].max(Dimensions::TwoDimensional),
390 'F' => {}
391 other => {
392 let message = format!("expected '0', '1', '2', or 'F'. Found: {other}");
393 return Err(InvalidInputError::new(message));
394 }
395 }
396 }
397 }
398
399 Ok(())
400 }
401
402 // NOTE for implementers
403 // See https://en.wikipedia.org/wiki/DE-9IM#Spatial_predicates for a mapping between predicates and matrices
404 // The number of constraints in your relation function MUST match the number of NON-MASK (T or F) matrix entries
405
406 // Indexes of the IntersectionMatrix map to indexes of a DE-9IM specification string as follows:
407 // ==================================================================
408 // self.0[CoordPos::Inside][CoordPos::Inside]: 0
409 // self.0[CoordPos::Inside][CoordPos::OnBoundary]: 1
410 // self.0[CoordPos::Inside][CoordPos::Outside]: 2
411
412 // self.0[CoordPos::OnBoundary][CoordPos::Inside]: 3
413 // self.0[CoordPos::OnBoundary][CoordPos::OnBoundary]: 4
414 // self.0[CoordPos::OnBoundary][CoordPos::Outside]: 5
415
416 // self.0[CoordPos::Outside][CoordPos::Inside]: 6
417 // self.0[CoordPos::Outside][CoordPos::OnBoundary]: 7
418 // self.0[CoordPos::Outside][CoordPos::Outside]: 8
419 // ==================================================================
420
421 // Relationship between matrix entry and Dimensions
422 // ==================================================================
423 // A `T` entry translates to `!= Dimensions::Empty`
424 // An `F` entry translates to `== Dimensions::Empty`
425 // A `*` (mask) entry is OMITTED
426 // ==================================================================
427
428 // Examples
429 // ==================================================================
430 // `[T********]` -> `self.0[CoordPos::Inside][CoordPos::Inside] != Dimensions::Empty`
431 // `[********F]` -> `self.0[CoordPos::Outside][CoordPos::Outside] == Dimensions::Empty`
432 // `[**T****F*]` -> `self.0[CoordPos::Inside][CoordPos::Outside] != Dimensions::Empty
433 // && self.0[CoordPos::Outside][CoordPos::OnBoundary] == Dimensions::Empty`
434 // ==================================================================
435
436 /// Returns `true` if geometries `a` and `b` are disjoint: they have no point in common,
437 /// forming a set of disconnected geometries.
438 ///
439 /// # Notes
440 /// - Matches `[FF*FF****]`
441 /// - This predicate is **anti-reflexive**
442 pub fn is_disjoint(&self) -> bool {
443 self.0[CoordPos::Inside][CoordPos::Inside] == Dimensions::Empty
444 && self.0[CoordPos::Inside][CoordPos::OnBoundary] == Dimensions::Empty
445 && self.0[CoordPos::OnBoundary][CoordPos::Inside] == Dimensions::Empty
446 && self.0[CoordPos::OnBoundary][CoordPos::OnBoundary] == Dimensions::Empty
447 }
448
449 /// Tests if [`IntersectionMatrix::is_disjoint`] returns `false`.
450 ///
451 /// Returns `true` if the two geometries related by this matrix intersect: they have at least one point in common.
452 ///
453 /// # Notes
454 /// - Matches any of `[T********], [*T*******], [***T*****], [****T****]`
455 /// - This predicate is **reflexive and symmetric**
456 pub fn is_intersects(&self) -> bool {
457 !self.is_disjoint()
458 }
459
460 /// Returns `true` if the first geometry is within the second: `a` lies in the interior of `b`.
461 ///
462 ///
463 /// # Notes
464 /// - Also known as **inside**
465 /// - The mask `[T*F**F***`] occurs in the definition of both [`IntersectionMatrix::is_within`] and [`IntersectionMatrix::is_coveredby`]; For **most** situations, [`IntersectionMatrix::is_coveredby`] should be used in preference to [`IntersectionMatrix::is_within`]
466 /// - This predicate is **reflexive and transitive**
467 pub fn is_within(&self) -> bool {
468 self.0[CoordPos::Inside][CoordPos::Inside] != Dimensions::Empty
469 && self.0[CoordPos::Inside][CoordPos::Outside] == Dimensions::Empty
470 && self.0[CoordPos::OnBoundary][CoordPos::Outside] == Dimensions::Empty
471 }
472
473 /// Returns `true` if geometry `a` contains geometry `b`.
474 ///
475 /// # Notes
476 /// - Matches `[T*****FF*]`
477 /// - This predicate is **reflexive and transitive**
478 pub fn is_contains(&self) -> bool {
479 self.0[CoordPos::Inside][CoordPos::Inside] != Dimensions::Empty
480 && self.0[CoordPos::Outside][CoordPos::Inside] == Dimensions::Empty
481 && self.0[CoordPos::Outside][CoordPos::OnBoundary] == Dimensions::Empty
482 }
483
484 /// Returns `true` if the first geometry is *topologically* equal to the second.
485 ///
486 /// # Notes
487 /// - Matches `[T*F**FFF*]`
488 /// - This predicate is **reflexive, symmetric, and transitive**
489 pub fn is_equal_topo(&self) -> bool {
490 if self == &Self::empty_disjoint() {
491 // Any two empty geometries are topologically equal
492 return true;
493 }
494
495 self.0[CoordPos::Inside][CoordPos::Inside] != Dimensions::Empty
496 && self.0[CoordPos::Inside][CoordPos::Outside] == Dimensions::Empty
497 && self.0[CoordPos::Outside][CoordPos::Inside] == Dimensions::Empty
498 && self.0[CoordPos::Outside][CoordPos::OnBoundary] == Dimensions::Empty
499 && self.0[CoordPos::OnBoundary][CoordPos::Outside] == Dimensions::Empty
500 }
501
502 /// Returns true if every point in Geometry `a` lies inside (i.e. intersects the interior or boundary of) Geometry `b`.
503 ///
504 /// Equivalently, tests that no point of `a` lies outside (in the exterior of) `b`:
505 /// - `a` is covered by `b` (extends [`IntersectionMatrix::is_within`]): geometry `a` lies in `b`. OR
506 /// - At least **one** point of `a` lies in `b`, and **no** point of `a` lies in the **exterior** of `b` OR
507 /// - **Every** point of `a` is a point of (the **interior** or **boundary** of) `b`
508 ///
509 /// returns `true` if the first geometry is covered by the second.
510 ///
511 /// ```
512 /// use geo_types::{Polygon, polygon};
513 /// use geo::relate::Relate;
514 ///
515 /// let poly1 = polygon![
516 /// (x: 125., y: 179.),
517 /// (x: 110., y: 150.),
518 /// (x: 160., y: 160.),
519 /// (x: 125., y: 179.),
520 /// ];
521 /// let poly2 = polygon![
522 /// (x: 124., y: 182.),
523 /// (x: 106., y: 146.),
524 /// (x: 162., y: 159.),
525 /// (x: 124., y: 182.),
526 /// ];
527 ///
528 /// let intersection = poly1.relate(&poly2);
529 /// assert_eq!(intersection.is_coveredby(), true);
530 /// ```
531 ///
532 /// # Notes
533 /// - Matches any of `[T*F**F***], [*TF**F***], [**FT*F***], [**F*TF***]`
534 /// - This predicate is **reflexive and transitive**
535 #[allow(clippy::nonminimal_bool)]
536 pub fn is_coveredby(&self) -> bool {
537 // [T*F**F***]
538 self.0[CoordPos::Inside][CoordPos::Inside] != Dimensions::Empty
539 && self.0[CoordPos::Inside][CoordPos::Outside] == Dimensions::Empty
540 && self.0[CoordPos::OnBoundary][CoordPos::Outside] == Dimensions::Empty ||
541 // [*TF**F***]
542 self.0[CoordPos::Inside][CoordPos::OnBoundary] != Dimensions::Empty
543 && self.0[CoordPos::Inside][CoordPos::Outside] == Dimensions::Empty
544 && self.0[CoordPos::OnBoundary][CoordPos::Outside] == Dimensions::Empty ||
545 // [**FT*F***]
546 self.0[CoordPos::Inside][CoordPos::Outside] == Dimensions::Empty
547 && self.0[CoordPos::OnBoundary][CoordPos::Inside] != Dimensions::Empty
548 && self.0[CoordPos::OnBoundary][CoordPos::Outside] == Dimensions::Empty ||
549 // [**F*TF***]
550 self.0[CoordPos::Inside][CoordPos::Outside] == Dimensions::Empty
551 && self.0[CoordPos::OnBoundary][CoordPos::OnBoundary] != Dimensions::Empty
552 && self.0[CoordPos::OnBoundary][CoordPos::Outside] == Dimensions::Empty
553 }
554
555 /// Returns `true` if every point in Geometry `b` lies inside
556 /// (i.e. intersects the interior or boundary of) Geometry `a`. Equivalently,
557 /// tests that no point of `b` lies outside (in the exterior of) `a`.
558 ///
559 /// # Notes
560 /// - Unlike [`IntersectionMatrix::is_contains`], it does **not** distinguish between points in the boundary and in the interior of geometries
561 /// - For **most** situations, [`IntersectionMatrix::is_covers`] should be used in preference to [`IntersectionMatrix::is_contains`]
562 /// - Matches any of `[T*****FF*], [*T****FF*], [***T**FF*], [****T*FF*]`
563 /// - This predicate is **reflexive and transitive**
564 #[allow(clippy::nonminimal_bool)]
565 pub fn is_covers(&self) -> bool {
566 // [T*****FF*]
567 self.0[CoordPos::Inside][CoordPos::Inside] != Dimensions::Empty
568 && self.0[CoordPos::Outside][CoordPos::Inside] == Dimensions::Empty
569 && self.0[CoordPos::Outside][CoordPos::OnBoundary] == Dimensions::Empty ||
570 // [*T****FF*]
571 self.0[CoordPos::Inside][CoordPos::OnBoundary] != Dimensions::Empty
572 && self.0[CoordPos::Outside][CoordPos::Inside] == Dimensions::Empty
573 && self.0[CoordPos::Outside][CoordPos::OnBoundary] == Dimensions::Empty ||
574 // [***T**FF*]
575 self.0[CoordPos::OnBoundary][CoordPos::Inside] != Dimensions::Empty
576 && self.0[CoordPos::Outside][CoordPos::Inside] == Dimensions::Empty
577 && self.0[CoordPos::Outside][CoordPos::OnBoundary] == Dimensions::Empty ||
578 // [****T*FF*]
579 self.0[CoordPos::OnBoundary][CoordPos::OnBoundary] != Dimensions::Empty
580 && self.0[CoordPos::Outside][CoordPos::Inside] == Dimensions::Empty
581 && self.0[CoordPos::Outside][CoordPos::OnBoundary] == Dimensions::Empty
582 }
583
584 /// Returns `true` if every point in Geometry `b` is a point of Geometry `a` interior.
585 ///
586 /// # Notes
587 /// - If Geometry `b` has any interaction with the boundary of Geometry `a`, then the result is `false`.
588 /// - Geometry`a` will never contains_properly itself.
589 /// - Matches `[T**FF*FF*]`
590 /// - This predicate is **transitive** but not **reflexive**
591 ///
592 /// # Performance
593 ///
594 /// [`ContainsProperly`] is faster when checking between `Polygon` and `MultiPolygon` when inputs are smaller than about 650 vertices, otherwise use `IntersectionMatrix::is_contains_properly`
595 #[allow(clippy::nonminimal_bool)]
596 pub fn is_contains_properly(&self) -> bool {
597 // [T**FF*FF*]
598 self.0[CoordPos::Inside][CoordPos::Inside] != Dimensions::Empty
599 && self.0[CoordPos::OnBoundary][CoordPos::Inside] == Dimensions::Empty
600 && self.0[CoordPos::OnBoundary][CoordPos::OnBoundary] == Dimensions::Empty
601 && self.0[CoordPos::Outside][CoordPos::Inside] == Dimensions::Empty
602 && self.0[CoordPos::Outside][CoordPos::OnBoundary] == Dimensions::Empty
603 }
604
605 /// Returns `true` if `a` touches `b`: they have at least one point in common, but their
606 /// interiors do not intersect.
607 ///
608 /// # Notes
609 /// - Matches any of `[FT*******], [F**T*****], [F***T****]`
610 /// - This predicate is **symmetric**
611 #[allow(clippy::nonminimal_bool)]
612 pub fn is_touches(&self) -> bool {
613 // [FT*******]
614 self.0[CoordPos::Inside][CoordPos::Inside] == Dimensions::Empty
615 && self.0[CoordPos::Inside][CoordPos::OnBoundary] != Dimensions::Empty ||
616 // [F**T*****]
617 self.0[CoordPos::Inside][CoordPos::Inside] == Dimensions::Empty
618 && self.0[CoordPos::OnBoundary][CoordPos::Inside] != Dimensions::Empty ||
619 // [F***T****]
620 self.0[CoordPos::Inside][CoordPos::Inside] == Dimensions::Empty
621 && self.0[CoordPos::OnBoundary][CoordPos::OnBoundary] != Dimensions::Empty
622 }
623
624 /// Compares two geometry objects and returns `true` if their intersection "spatially crosses";
625 /// that is, the geometries have some, but not all interior points in common
626 ///
627 /// ```
628 /// use geo_types::{LineString, line_string, polygon};
629 /// use geo::relate::Relate;
630 ///
631 /// let line_string: LineString = line_string![(x: 85.0, y: 194.0), (x: 162.0, y: 135.0)];
632 /// let poly = polygon![
633 /// (x: 125., y: 179.),
634 /// (x: 110., y: 150.),
635 /// (x: 160., y: 160.),
636 /// (x: 125., y: 179.),
637 /// ];
638 ///
639 /// let intersection = line_string.relate(&poly);
640 /// assert_eq!(intersection.is_crosses(), true);
641 /// ```
642 ///
643 /// # Notes
644 /// - If any of the following do not hold, the function will return `false`:
645 /// - The intersection of the interiors of the geometries must be non-empty
646 /// - The intersection must have dimension less than the maximum dimension of the two input geometries (two polygons cannot cross)
647 /// - The intersection of the two geometries must not equal either geometry (two points cannot cross)
648 /// - Matches one of `[T*T******] (a < b)`, `[T*****T**] (a > b)`, `[0********] (dimensions == 1)`
649 /// - This predicate is **symmetric and irreflexive**
650 pub fn is_crosses(&self) -> bool {
651 let dims_a = self.0[CoordPos::Inside][CoordPos::Inside]
652 .max(self.0[CoordPos::Inside][CoordPos::OnBoundary])
653 .max(self.0[CoordPos::Inside][CoordPos::Outside]);
654
655 let dims_b = self.0[CoordPos::Inside][CoordPos::Inside]
656 .max(self.0[CoordPos::OnBoundary][CoordPos::Inside])
657 .max(self.0[CoordPos::Outside][CoordPos::Inside]);
658 match (dims_a, dims_b) {
659 // a < b
660 _ if dims_a < dims_b =>
661 // [T*T******]
662 {
663 self.0[CoordPos::Inside][CoordPos::Inside] != Dimensions::Empty
664 && self.0[CoordPos::Inside][CoordPos::Outside] != Dimensions::Empty
665 }
666 // a > b
667 _ if dims_a > dims_b =>
668 // [T*****T**]
669 {
670 self.0[CoordPos::Inside][CoordPos::Inside] != Dimensions::Empty
671 && self.0[CoordPos::Outside][CoordPos::Inside] != Dimensions::Empty
672 }
673 // a == b, only line / line permitted
674 (Dimensions::OneDimensional, Dimensions::OneDimensional) =>
675 // [0********]
676 {
677 self.0[CoordPos::Inside][CoordPos::Inside] == Dimensions::ZeroDimensional
678 }
679 _ => false,
680 }
681 }
682
683 /// Returns `true` if geometry `a` and `b` "spatially overlap". Two geometries overlap if they have the
684 /// same dimension, their interiors intersect in that dimension, and each has at least one point
685 /// inside the other (or equivalently, neither one covers the other)
686 ///
687 /// ```
688 /// use geo_types::{Polygon, polygon};
689 /// use geo::relate::Relate;
690 ///
691 /// let poly1 = polygon![
692 /// (x: 125., y: 179.),
693 /// (x: 110., y: 150.),
694 /// (x: 160., y: 160.),
695 /// (x: 125., y: 179.),
696 /// ];
697 /// let poly2 = polygon![
698 /// (x: 126., y: 179.),
699 /// (x: 110., y: 150.),
700 /// (x: 161., y: 160.),
701 /// (x: 126., y: 179.),
702 /// ];
703 ///
704 /// let intersection = poly1.relate(&poly2);
705 /// assert_eq!(intersection.is_overlaps(), true);
706 /// ```
707 ///
708 /// # Notes
709 /// - Matches one of `[1*T***T**] (dimensions == 1)`, `[T*T***T**] (dimensions == 0 OR 2)`
710 /// - This predicate is **symmetric**
711 #[allow(clippy::nonminimal_bool)]
712 pub fn is_overlaps(&self) -> bool {
713 // dimensions must be non-empty, equal, and line / line is a special case
714 let dims_a = self.0[CoordPos::Inside][CoordPos::Inside]
715 .max(self.0[CoordPos::Inside][CoordPos::OnBoundary])
716 .max(self.0[CoordPos::Inside][CoordPos::Outside]);
717
718 let dims_b = self.0[CoordPos::Inside][CoordPos::Inside]
719 .max(self.0[CoordPos::OnBoundary][CoordPos::Inside])
720 .max(self.0[CoordPos::Outside][CoordPos::Inside]);
721 match (dims_a, dims_b) {
722 // line / line: [1*T***T**]
723 (Dimensions::OneDimensional, Dimensions::OneDimensional) => {
724 self.0[CoordPos::Inside][CoordPos::Inside] == Dimensions::OneDimensional
725 && self.0[CoordPos::Inside][CoordPos::Outside] != Dimensions::Empty
726 && self.0[CoordPos::Outside][CoordPos::Inside] != Dimensions::Empty
727 }
728 // point / point or polygon / polygon: [T*T***T**]
729 (Dimensions::ZeroDimensional, Dimensions::ZeroDimensional)
730 | (Dimensions::TwoDimensional, Dimensions::TwoDimensional) => {
731 self.0[CoordPos::Inside][CoordPos::Inside] != Dimensions::Empty
732 && self.0[CoordPos::Inside][CoordPos::Outside] != Dimensions::Empty
733 && self.0[CoordPos::Outside][CoordPos::Inside] != Dimensions::Empty
734 }
735 _ => false,
736 }
737 }
738
739 /// Directly accesses this matrix
740 ///
741 /// ```
742 /// use geo_types::{LineString, Rect, line_string};
743 /// use geo::{coordinate_position::CoordPos, dimensions::Dimensions, relate::Relate};
744 ///
745 /// let line_string: LineString = line_string![(x: 0.0, y: 0.0), (x: 10.0, y: 0.0), (x: 5.0, y: 5.0)];
746 /// let rect = Rect::new((0.0, 0.0), (5.0, 5.0));
747 ///
748 /// let intersection = line_string.relate(&rect);
749 ///
750 /// // The intersection of the two interiors is empty, because no part of the string is inside the rect
751 /// assert_eq!(intersection.get(CoordPos::Inside, CoordPos::Inside), Dimensions::Empty);
752 ///
753 /// // The intersection of the line string's interior with the rect's boundary is one-dimensional, because part of the first line overlaps one of the rect's edges
754 /// assert_eq!(intersection.get(CoordPos::Inside, CoordPos::OnBoundary), Dimensions::OneDimensional);
755 ///
756 /// // The intersection of the line string's interior with the rect's exterior is one-dimensional, because part of the string is outside the rect
757 /// assert_eq!(intersection.get(CoordPos::Inside, CoordPos::Outside), Dimensions::OneDimensional);
758 ///
759 /// // The intersection of the line string's boundary with the rect's interior is empty, because neither of its end points are inside the rect
760 /// assert_eq!(intersection.get(CoordPos::OnBoundary, CoordPos::Inside), Dimensions::Empty);
761 ///
762 /// // The intersection of the line string's boundary with the rect's boundary is zero-dimensional, because the string's start and end points are on the rect's edges
763 /// assert_eq!(intersection.get(CoordPos::OnBoundary, CoordPos::OnBoundary), Dimensions::ZeroDimensional);
764 ///
765 /// // The intersection of the line string's boundary with the rect's exterior is empty, because neither of its end points are outside the rect
766 /// assert_eq!(intersection.get(CoordPos::OnBoundary, CoordPos::Outside), Dimensions::Empty);
767 ///
768 /// // The intersection of the the line's exterior with the rect's interior is two-dimensional, because it's simply the rect's interior
769 /// assert_eq!(intersection.get(CoordPos::Outside, CoordPos::Inside), Dimensions::TwoDimensional);
770 ///
771 /// // The intersection of the line's exterior with the rect's boundary is one-dimensional, because it's the rect's edges (minus where the string overlaps it)
772 /// assert_eq!(intersection.get(CoordPos::Outside, CoordPos::OnBoundary), Dimensions::OneDimensional);
773 ///
774 /// // The intersection of the two exteriors is two-dimensional, because it's the whole plane around the two shapes
775 /// assert_eq!(intersection.get(CoordPos::Outside, CoordPos::Outside), Dimensions::TwoDimensional);
776 /// ```
777 pub fn get(&self, lhs: CoordPos, rhs: CoordPos) -> Dimensions {
778 self.0[lhs][rhs]
779 }
780
781 /// Does the intersection matrix match the provided DE-9IM specification string?
782 ///
783 /// A DE-9IM spec string must be 9 characters long, and each character
784 /// must be one of the following:
785 ///
786 /// - 0: matches a 0-dimensional (point) intersection
787 /// - 1: matches a 1-dimensional (line) intersection
788 /// - 2: matches a 2-dimensional (area) intersection
789 /// - f or F: matches only empty dimensions
790 /// - t or T: matches anything non-empty
791 /// - *: matches anything
792 ///
793 /// ```
794 /// use geo::algorithm::Relate;
795 /// use geo::geometry::Polygon;
796 /// use wkt::TryFromWkt;
797 ///
798 /// let a = Polygon::<f64>::try_from_wkt_str("POLYGON((0 0,4 0,4 4,0 4,0 0))").expect("valid WKT");
799 /// let b = Polygon::<f64>::try_from_wkt_str("POLYGON((1 1,4 0,4 4,0 4,1 1))").expect("valid WKT");
800 /// let im = a.relate(&b);
801 /// assert!(im.matches("212F11FF2").expect("valid DE-9IM spec"));
802 /// assert!(im.matches("TTT***FF2").expect("valid DE-9IM spec"));
803 /// assert!(!im.matches("TTT***FFF").expect("valid DE-9IM spec"));
804 /// ```
805 pub fn matches(&self, spec: &str) -> Result<bool, InvalidInputError> {
806 if spec.len() != 9 {
807 return Err(InvalidInputError::new(format!(
808 "DE-9IM specification must be exactly 9 characters. Got {len}",
809 len = spec.len()
810 )));
811 }
812
813 let mut chars = spec.chars();
814 for a in &[CoordPos::Inside, CoordPos::OnBoundary, CoordPos::Outside] {
815 for b in &[CoordPos::Inside, CoordPos::OnBoundary, CoordPos::Outside] {
816 let dim_spec = dimension_matcher::DimensionMatcher::try_from(
817 chars.next().expect("already validated length is 9"),
818 )?;
819 if !dim_spec.matches(self.0[*a][*b]) {
820 return Ok(false);
821 }
822 }
823 }
824 Ok(true)
825 }
826}
827
828/// Build an IntersectionMatrix based on a string specification.
829/// ```
830/// use geo::algorithm::relate::IntersectionMatrix;
831/// use std::str::FromStr;
832///
833/// let intersection_matrix = IntersectionMatrix::from_str("212101212").expect("valid DE-9IM specification");
834/// assert!(intersection_matrix.is_intersects());
835/// assert!(!intersection_matrix.is_contains());
836/// ```
837impl FromStr for IntersectionMatrix {
838 type Err = InvalidInputError;
839 fn from_str(str: &str) -> Result<Self, Self::Err> {
840 let mut im = IntersectionMatrix::empty();
841 im.set_at_least_from_string(str)?;
842 Ok(im)
843 }
844}
845
846pub(crate) mod dimension_matcher {
847 use super::Dimensions;
848 use super::InvalidInputError;
849
850 /// A single letter from a DE-9IM matching specification like "1*T**FFF*"
851 pub(crate) enum DimensionMatcher {
852 Anything,
853 NonEmpty,
854 Exact(Dimensions),
855 }
856
857 impl DimensionMatcher {
858 pub fn matches(&self, dim: Dimensions) -> bool {
859 match (self, dim) {
860 (Self::Anything, _) => true,
861 (DimensionMatcher::NonEmpty, d) => d != Dimensions::Empty,
862 (DimensionMatcher::Exact(a), b) => a == &b,
863 }
864 }
865 }
866
867 impl TryFrom<char> for DimensionMatcher {
868 type Error = InvalidInputError;
869
870 fn try_from(value: char) -> Result<Self, Self::Error> {
871 Ok(match value {
872 '*' => Self::Anything,
873 't' | 'T' => Self::NonEmpty,
874 'f' | 'F' => Self::Exact(Dimensions::Empty),
875 '0' => Self::Exact(Dimensions::ZeroDimensional),
876 '1' => Self::Exact(Dimensions::OneDimensional),
877 '2' => Self::Exact(Dimensions::TwoDimensional),
878 _ => {
879 return Err(InvalidInputError::new(format!(
880 "invalid DE-9IM specification character: {value}"
881 )));
882 }
883 })
884 }
885 }
886}
887
888#[cfg(test)]
889mod tests {
890 use super::*;
891 use crate::algorithm::Relate;
892 use crate::geometry::*;
893 use crate::wkt;
894
895 #[test]
896 fn test_crosses() {
897 // these polygons look like they cross, but two polygons cannot cross
898 let a: Geometry<_> = wkt! { POLYGON ((3.4 15.7, 2.2 11.3, 5.8 11.4, 3.4 15.7)) }.into();
899 let b: Geometry<_> = wkt! { POLYGON ((5.2 13.1, 4.5 10.9, 6.3 11.1, 5.2 13.1)) }.into();
900 // this linestring is a single leg of b: it can cross polygon a
901 let c: Geometry<_> = wkt! { LINESTRING (5.2 13.1, 4.5 10.9) }.into();
902 let relate_ab = a.relate(&b);
903 let relate_ca = c.relate(&a);
904 assert!(!relate_ab.is_crosses());
905 assert!(relate_ca.is_crosses());
906 }
907
908 #[test]
909 fn test_crosses_2() {
910 // two lines can cross
911 // same geometry as test_crosses: single legs of polygons a and b
912 let a: Geometry<_> = wkt! { LINESTRING (5.2 13.1, 4.5 10.9) }.into();
913 let b: Geometry<_> = wkt! { LINESTRING (3.4 15.7, 2.2 11.3, 5.8 11.4) }.into();
914 let relate_ab = a.relate(&b);
915 assert!(relate_ab.is_crosses());
916 }
917
918 mod test_matches {
919 use super::*;
920
921 fn subject() -> IntersectionMatrix {
922 // Topologically, this is a nonsense IM
923 IntersectionMatrix::from_str("F00111222").unwrap()
924 }
925
926 #[test]
927 fn matches_exactly() {
928 assert!(subject().matches("F00111222").unwrap());
929 }
930
931 #[test]
932 fn doesnt_match() {
933 assert!(!subject().matches("222222222").unwrap());
934 }
935
936 #[test]
937 fn matches_truthy() {
938 assert!(subject().matches("FTTTTTTTT").unwrap());
939 }
940
941 #[test]
942 fn matches_wildcard() {
943 assert!(subject().matches("F0011122*").unwrap());
944 }
945 }
946
947 #[test]
948 fn empty_is_equal_topo() {
949 let empty_polygon = Polygon::<f64>::empty();
950 let im = empty_polygon.relate(&empty_polygon);
951 assert!(im.is_equal_topo());
952 }
953}