Skip to main content

symtropy_math/
hyperplane.rs

1// Copyright (C) 2024-2026 Tristan Stoltz / Luminous Dynamics
2// SPDX-License-Identifier: Apache-2.0 OR MIT
3// Commercial licensing: see COMMERCIAL_LICENSE.md at repository root
4use crate::point::Point;
5use nalgebra::SVector;
6
7/// D-dimensional hyperplane: { x : normal ยท x = offset }.
8/// Used for BSP partitioning and 4D cross-section slicing.
9#[derive(Clone, Copy, Debug)]
10pub struct Hyperplane<const D: usize> {
11    pub normal: SVector<f64, D>,
12    pub offset: f64,
13}
14
15impl<const D: usize> Hyperplane<D> {
16    pub fn from_normal_and_point(normal: SVector<f64, D>, point: &Point<D>) -> Self {
17        let n = normal.normalize();
18        Self {
19            normal: n,
20            offset: n.dot(&point.0),
21        }
22    }
23
24    pub fn new(normal: SVector<f64, D>, offset: f64) -> Self {
25        Self { normal, offset }
26    }
27
28    #[inline]
29    pub fn signed_distance(&self, point: &Point<D>) -> f64 {
30        self.normal.dot(&point.0) - self.offset
31    }
32
33    #[inline]
34    pub fn classify(&self, point: &Point<D>) -> Side {
35        let d = self.signed_distance(point);
36        if d > 1e-10 {
37            Side::Front
38        } else if d < -1e-10 {
39            Side::Back
40        } else {
41            Side::On
42        }
43    }
44
45    pub fn project(&self, point: &Point<D>) -> Point<D> {
46        let d = self.signed_distance(point);
47        Point(point.0 - self.normal * d)
48    }
49
50    pub fn intersect_segment(&self, a: &Point<D>, b: &Point<D>) -> Option<f64> {
51        let da = self.signed_distance(a);
52        let db = self.signed_distance(b);
53        if da * db > 0.0 {
54            return None;
55        }
56        let denom = da - db;
57        if denom.abs() < 1e-15 {
58            return None;
59        }
60        Some(da / denom)
61    }
62}
63
64#[derive(Clone, Copy, Debug, PartialEq, Eq)]
65pub enum Side {
66    Front,
67    Back,
68    On,
69}
70
71#[cfg(test)]
72mod tests {
73    use super::*;
74
75    #[test]
76    fn classify_3d() {
77        let hp = Hyperplane::new(SVector::from([0.0, 1.0, 0.0]), 0.0);
78        assert_eq!(hp.classify(&Point::new([0.0, 1.0, 0.0])), Side::Front);
79        assert_eq!(hp.classify(&Point::new([0.0, -1.0, 0.0])), Side::Back);
80        assert_eq!(hp.classify(&Point::new([5.0, 0.0, 3.0])), Side::On);
81    }
82
83    #[test]
84    fn signed_distance() {
85        let hp = Hyperplane::new(SVector::from([1.0, 0.0]), 5.0);
86        let p = Point::new([8.0, 3.0]);
87        assert!((hp.signed_distance(&p) - 3.0).abs() < 1e-12);
88    }
89
90    #[test]
91    fn project() {
92        let hp = Hyperplane::new(SVector::from([0.0, 1.0, 0.0]), 0.0);
93        let p = Point::new([3.0, 7.0, 2.0]);
94        let proj = hp.project(&p);
95        assert!((proj.coord(0) - 3.0).abs() < 1e-12);
96        assert!(proj.coord(1).abs() < 1e-12);
97    }
98
99    #[test]
100    fn segment_intersection() {
101        let hp = Hyperplane::new(SVector::from([0.0, 1.0, 0.0]), 0.0);
102        let a = Point::new([0.0, -1.0, 0.0]);
103        let b = Point::new([0.0, 1.0, 0.0]);
104        let t = hp.intersect_segment(&a, &b).unwrap();
105        assert!((t - 0.5).abs() < 1e-12);
106    }
107
108    #[test]
109    fn no_intersection_same_side() {
110        let hp = Hyperplane::new(SVector::from([0.0, 1.0]), 0.0);
111        assert!(hp
112            .intersect_segment(&Point::new([0.0, 1.0]), &Point::new([0.0, 2.0]))
113            .is_none());
114    }
115
116    #[test]
117    fn hyperplane_4d_slicing() {
118        let hp = Hyperplane::new(SVector::from([0.0, 0.0, 0.0, 1.0]), 0.5);
119        assert_eq!(hp.classify(&Point::new([1.0, 2.0, 3.0, 0.5])), Side::On);
120        assert_eq!(hp.classify(&Point::new([1.0, 2.0, 3.0, 1.0])), Side::Front);
121        assert_eq!(hp.classify(&Point::new([1.0, 2.0, 3.0, 0.0])), Side::Back);
122    }
123}