symtropy_math/
hyperplane.rs1use crate::point::Point;
5use nalgebra::SVector;
6
7#[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}