geometry_algorithm/
discrete_hausdorff.rs

1//! `discrete_hausdorff_distance(&l1, &l2)` — sup-sup distance over
2//! vertex sets.
3//!
4//! Mirrors `boost::geometry::discrete_hausdorff_distance` from
5//! `boost/geometry/algorithms/discrete_hausdorff_distance.hpp`. The
6//! Boost overload is symmetric — `max(directed(A,B), directed(B,A))`
7//! where `directed(A,B) = max_{p ∈ A} min_{q ∈ B} dist(p, q)` — the
8//! Rust port matches. `O(m × n)` time, `O(1)` space.
9
10use alloc::vec::Vec;
11
12use geometry_cs::CoordinateSystem;
13use geometry_strategy::distance::DefaultDistance;
14use geometry_strategy::{DefaultDistanceStrategy, DistanceStrategy};
15use geometry_trait::{Geometry, Linestring, Point};
16
17/// Shorthand for the coordinate-system family of a point type.
18type Family<P> = <<P as Point>::Cs as CoordinateSystem>::Family;
19
20/// The scalar output type of the default distance strategy between the
21/// point types of two linestrings.
22type DefaultDistOut<L1, L2> = <DefaultDistanceStrategy<
23    <L1 as Geometry>::Point,
24    <L2 as Geometry>::Point,
25> as DistanceStrategy<<L1 as Geometry>::Point, <L2 as Geometry>::Point>>::Out;
26
27/// Symmetric discrete Hausdorff distance between two linestrings, using
28/// the default distance strategy for their coordinate systems.
29///
30/// Mirrors `boost::geometry::discrete_hausdorff_distance(l1, l2)` from
31/// `boost/geometry/algorithms/discrete_hausdorff_distance.hpp`.
32///
33/// # Panics
34///
35/// Panics on an empty linestring — Boost treats empty input as an
36/// error; the Rust port panics with a clear message.
37#[inline]
38#[must_use]
39pub fn discrete_hausdorff_distance<L1, L2>(l1: &L1, l2: &L2) -> DefaultDistOut<L1, L2>
40where
41    L1: Linestring,
42    L2: Linestring,
43    Family<L1::Point>: DefaultDistance<Family<L2::Point>>,
44    DefaultDistanceStrategy<L1::Point, L2::Point>: DistanceStrategy<L1::Point, L2::Point> + Default,
45{
46    let s = <DefaultDistanceStrategy<L1::Point, L2::Point>>::default();
47    discrete_hausdorff_distance_with(l1, l2, s)
48}
49
50/// Symmetric discrete Hausdorff distance using an explicit distance
51/// strategy `dist`.
52///
53/// Mirrors the strategy-taking `boost::geometry::discrete_hausdorff_distance`
54/// overload from
55/// `boost/geometry/algorithms/discrete_hausdorff_distance.hpp`. The
56/// reverse direction reuses the same strategy with the arguments
57/// swapped, so a single `DistanceStrategy<L1::Point, L2::Point>` serves
58/// both directed suprema.
59///
60/// # Panics
61///
62/// Panics if either linestring is empty.
63#[must_use]
64#[allow(
65    clippy::needless_pass_by_value,
66    reason = "Distance strategies are zero-sized/Copy; taking by value matches `distance_with`."
67)]
68pub fn discrete_hausdorff_distance_with<L1, L2, S>(l1: &L1, l2: &L2, dist: S) -> S::Out
69where
70    L1: Linestring,
71    L2: Linestring,
72    S: DistanceStrategy<L1::Point, L2::Point>,
73    S::Out: PartialOrd + Copy,
74{
75    let seq1: Vec<&L1::Point> = l1.points().collect();
76    let seq2: Vec<&L2::Point> = l2.points().collect();
77    assert!(
78        !seq1.is_empty() && !seq2.is_empty(),
79        "empty linestring in discrete_hausdorff"
80    );
81
82    // directed(A, B): outer over `seq1` (L1 points), inner over `seq2`.
83    let sup_ab = directed_sup(&seq1, &seq2, |p, q| dist.distance(p, q));
84    // directed(B, A): outer over `seq2`, inner over `seq1`. The distance
85    // call keeps its `(L1::Point, L2::Point)` argument order.
86    let sup_ba = directed_sup(&seq2, &seq1, |p, q| dist.distance(q, p));
87
88    if sup_ab > sup_ba { sup_ab } else { sup_ba }
89}
90
91/// Directed supremum `max_{p ∈ a} min_{q ∈ b} f(p, q)`.
92fn directed_sup<PA, PB, O, F>(a: &[&PA], b: &[&PB], mut f: F) -> O
93where
94    O: PartialOrd + Copy,
95    F: FnMut(&PA, &PB) -> O,
96{
97    let mut sup: Option<O> = None;
98    for p in a {
99        let mut inf: Option<O> = None;
100        for q in b {
101            let d = f(p, q);
102            inf = Some(match inf {
103                None => d,
104                Some(cur) => {
105                    if d < cur {
106                        d
107                    } else {
108                        cur
109                    }
110                }
111            });
112        }
113        let inf = inf.expect("non-empty inner sequence");
114        sup = Some(match sup {
115            None => inf,
116            Some(cur) => {
117                if inf > cur {
118                    inf
119                } else {
120                    cur
121                }
122            }
123        });
124    }
125    sup.expect("non-empty outer sequence")
126}
127
128#[cfg(test)]
129#[allow(
130    clippy::float_cmp,
131    reason = "Hausdorff reference values compared with an epsilon."
132)]
133mod tests {
134    //! Reference values from
135    //! `boost/geometry/test/algorithms/similarity/discrete_hausdorff_distance.cpp`.
136
137    use super::discrete_hausdorff_distance;
138    use geometry_cs::Cartesian;
139    use geometry_model::{Linestring, Point2D, linestring};
140
141    type Pt = Point2D<f64, Cartesian>;
142
143    #[test]
144    fn identical_linestrings_distance_zero() {
145        let ls: Linestring<Pt> = linestring![(0., 0.), (1., 0.), (2., 0.)];
146        assert!(discrete_hausdorff_distance(&ls, &ls) < 1e-12);
147    }
148
149    /// Subset case: dropping the last vertex — the Hausdorff distance is
150    /// the distance from the dropped vertex `(2,0)` to its nearest kept
151    /// neighbour `(1,0)`, i.e. 1.
152    #[test]
153    fn subset_drops_last_vertex() {
154        let a: Linestring<Pt> = linestring![(0., 0.), (1., 0.), (2., 0.)];
155        let b: Linestring<Pt> = linestring![(0., 0.), (1., 0.)];
156        assert!((discrete_hausdorff_distance(&a, &b) - 1.0).abs() < 1e-9);
157    }
158
159    /// Parallel tracks 1 unit apart with aligned endpoints → 1.
160    #[test]
161    fn parallel_lines() {
162        let a: Linestring<Pt> = linestring![(0., 0.), (5., 0.)];
163        let b: Linestring<Pt> = linestring![(0., 1.), (5., 1.)];
164        assert!((discrete_hausdorff_distance(&a, &b) - 1.0).abs() < 1e-9);
165    }
166
167    /// The symmetric property: `H(a, b) == H(b, a)`.
168    #[test]
169    fn symmetric() {
170        let a: Linestring<Pt> = linestring![(0., 0.), (1., 0.), (2., 0.)];
171        let b: Linestring<Pt> = linestring![(0., 0.), (1., 0.)];
172        assert_eq!(
173            discrete_hausdorff_distance(&a, &b).to_bits(),
174            discrete_hausdorff_distance(&b, &a).to_bits(),
175        );
176    }
177}
geometry_algorithm/discrete_hausdorff.rs

geometry_algorithm/
discrete_hausdorff.rs
