Skip to main content

dualis_core/
vector.rs

1//! Vector arithmetic that no single domain owns.
2//!
3//! An orthonormal basis about a direction is needed to sample a hemisphere, to
4//! orient a scatter lobe and to build a local frame for a contact normal.
5//! Reflection about a plane is a light ray off a mirror and a ball off a wall, in
6//! the same three lines. Neither belongs to optics, so neither lives there.
7//!
8//! Refraction does belong to optics โ€” Snell's law is about refractive indices โ€”
9//! and is in `dualis-optics` accordingly.
10
11use glam::DVec3;
12
13/// Two unit vectors completing an orthonormal frame with `dir`.
14///
15/// `dir` must be normalised. The helper axis is chosen by which component of
16/// `dir` is small, which keeps the cross product well conditioned โ€” picking a
17/// fixed helper would give a near-zero cross product whenever `dir` happened to
18/// be parallel to it, and the resulting basis would be numerical noise.
19pub fn basis_for(dir: DVec3) -> (DVec3, DVec3) {
20    let helper = if dir.x.abs() < 0.9 {
21        DVec3::X
22    } else {
23        DVec3::Y
24    };
25    let u = helper.cross(dir).normalize();
26    let v = dir.cross(u);
27    (u, v)
28}
29
30/// Flip `n` so it opposes the incident direction (`d ยท n < 0`).
31///
32/// A surface has two sides and a stored normal points at one of them. Whether a
33/// ray is arriving from the front or the back is a fact about the ray, not the
34/// surface, so the normal is oriented per hit rather than per surface.
35pub fn oriented_against(n: DVec3, d: DVec3) -> DVec3 {
36    if d.dot(n) > 0.0 {
37        -n
38    } else {
39        n
40    }
41}
42
43/// Mirror `d` about the plane with normal `n`. Both should be normalised.
44pub fn reflect(d: DVec3, n: DVec3) -> DVec3 {
45    (d - n * (2.0 * d.dot(n))).normalize()
46}
47
48#[cfg(test)]
49mod tests {
50    use super::*;
51
52    /// The basis is orthonormal and right-handed for every direction, including
53    /// the axis-aligned ones where a badly chosen helper vector would collapse it.
54    #[test]
55    fn the_basis_is_orthonormal_everywhere() {
56        let directions = [
57            DVec3::X,
58            DVec3::Y,
59            DVec3::Z,
60            -DVec3::X,
61            -DVec3::Z,
62            DVec3::new(1.0, 1.0, 1.0).normalize(),
63            DVec3::new(0.9999, 0.01, 0.0).normalize(),
64            DVec3::new(-0.7, 0.2, 0.68).normalize(),
65        ];
66        for d in directions {
67            let (u, v) = basis_for(d);
68            assert!((u.length() - 1.0).abs() < 1e-12, "u not unit for {d}");
69            assert!((v.length() - 1.0).abs() < 1e-12, "v not unit for {d}");
70            assert!(u.dot(d).abs() < 1e-12, "u not perpendicular to {d}");
71            assert!(v.dot(d).abs() < 1e-12, "v not perpendicular to {d}");
72            assert!(u.dot(v).abs() < 1e-12, "u and v not perpendicular for {d}");
73            // Right-handed: u x v = d.
74            assert!((u.cross(v) - d).length() < 1e-12, "left-handed for {d}");
75        }
76    }
77
78    /// Reflection reverses the normal component and keeps the tangential one, so
79    /// the angle in equals the angle out and the reflection of a reflection is the
80    /// original.
81    #[test]
82    fn reflection_preserves_the_angle_and_is_its_own_inverse() {
83        let n = DVec3::new(0.0, 1.0, 0.0);
84        let d = DVec3::new(1.0, -1.0, 0.0).normalize();
85        let r = reflect(d, n);
86        assert!((r - DVec3::new(1.0, 1.0, 0.0).normalize()).length() < 1e-12);
87        // Equal angles either side.
88        assert!((d.dot(n).abs() - r.dot(n).abs()).abs() < 1e-12);
89        // Reflecting the reflection comes back, up to direction.
90        assert!((reflect(r, n) - d).length() < 1e-12);
91        // Normal incidence turns straight around.
92        assert!((reflect(-n, n) - n).length() < 1e-12);
93    }
94
95    #[test]
96    fn a_normal_is_oriented_per_hit_not_per_surface() {
97        let n = DVec3::Z;
98        // Arriving from the front, the normal is left alone.
99        assert_eq!(oriented_against(n, -DVec3::Z), n);
100        // Arriving from behind, it is flipped to face the ray.
101        assert_eq!(oriented_against(n, DVec3::Z), -n);
102        // And either way the result opposes the ray, which is the contract.
103        for d in [DVec3::Z, -DVec3::Z, DVec3::new(0.3, 0.1, 0.9).normalize()] {
104            assert!(oriented_against(n, d).dot(d) <= 0.0);
105        }
106    }
107}