geop_core_sketch/space/
geometry.rs

1//! The geometry of a 3-D sketch's arcs, generic over [`Number`] so the same
2//! formulas serve the residuals (differentiated as
3//! [`geop_core_math::dual::Dual`]s) and building the NURBS.
4//!
5//! An arc is given by three points: from `s` through `m` to `e`. Everything
6//! about it follows from the two legs `a = s - m` and `b = e - m` without a
7//! single angle: the cosine of half its sweep is `-a.b / (|a| |b|)` and the
8//! sine `|a x b| / (|a| |b|)` — the angle the arc makes at `m` is the
9//! inscribed angle over the chord `se`.
10
11use geop_core_geometry::nurb_curve::cos_sin;
12use geop_core_math::{
13    geop_error::GeopResult,
14    scalars::{Number, Scalar},
15    vector::Vector3,
16};
17
18/// The helix round the axis from `base` to `top` from `start`, as far up
19/// the axis as `top` is from `base`, in `turns` turns (see
20/// [`super::CurveKind3d::Helix`]).
21#[derive(Clone, Copy, Debug)]
22pub struct Helix3<T: Number> {
23    pub start: Vector3<T>,
24    pub base: Vector3<T>,
25    pub top: Vector3<T>,
26    pub turns: Scalar,
27    pub left: bool,
28}
29
30impl<T: Number> Helix3<T> {
31    /// The unit axis, and the way from the axis to its start, square to it.
32    pub fn axes(&self) -> GeopResult<(Vector3<T>, Vector3<T>)> {
33        let w = self.top.sub(&self.base).normalize()?;
34        let off = self.start.sub(&self.base);
35        Ok((w, off.sub(&w.prod_scalar(off.prod_dot(&w)))))
36    }
37
38    /// Where it arrives: its start, up the axis, turned about it. The turn
39    /// is design data, so its cosine and sine are constants. A whole number
40    /// of turns comes back exactly to where it started, seen down the axis.
41    pub fn end(&self) -> GeopResult<Vector3<T>> {
42        let (w, out) = self.axes()?;
43        let degrees = 360.0 * self.turns.to_f64() * if self.left { -1.0 } else { 1.0 };
44        let (cos, sin) = cos_sin(degrees);
45        Ok(self
46            .start
47            .add(&self.top.sub(&self.base))
48            .add(&out.prod_scalar(T::from_f64(cos - 1.0)))
49            .add(&w.prod_cross(&out).prod_scalar(T::from_f64(sin))))
50    }
51}
52
53/// The circular arc from `s` through `m` to `e`.
54#[derive(Clone, Copy, Debug)]
55pub struct Arc3<T: Number> {
56    pub s: Vector3<T>,
57    pub m: Vector3<T>,
58    pub e: Vector3<T>,
59}
60
61impl<T: Number> Arc3<T> {
62    fn legs(&self) -> (Vector3<T>, Vector3<T>) {
63        (self.s.sub(&self.m), self.e.sub(&self.m))
64    }
65
66    /// The unit normal of its plane, the way it turns: from `s` to `e` it
67    /// runs counter-clockwise about it.
68    pub fn normal(&self) -> GeopResult<Vector3<T>> {
69        let (a, b) = self.legs();
70        b.prod_cross(&a).normalize()
71    }
72
73    /// `cos` and `sin` of half its sweep.
74    pub fn half_sweep(&self) -> GeopResult<(T, T)> {
75        let (a, b) = self.legs();
76        let ab = a.try_norm()?.mul(b.try_norm()?);
77        Ok((
78            a.prod_dot(&b).neg().div(ab)?,
79            a.prod_cross(&b).try_norm()?.div(ab)?,
80        ))
81    }
82
83    /// Its radius: `|a| |b| |e - s| / (2 |a x b|)`.
84    pub fn radius(&self) -> GeopResult<T> {
85        let (a, b) = self.legs();
86        a.try_norm()?
87            .mul(b.try_norm()?)
88            .mul(self.e.sub(&self.s).try_norm()?)
89            .div(T::TWO.mul(a.prod_cross(&b).try_norm()?))
90    }
91
92    /// Its center: `m + (|a|^2 b - |b|^2 a) x (a x b) / (2 |a x b|^2)`.
93    pub fn center(&self) -> GeopResult<Vector3<T>> {
94        let (a, b) = self.legs();
95        let n = a.prod_cross(&b);
96        let k = T::ONE.div(T::TWO.mul(n.norm_sq()))?;
97        let v = b.prod_scalar(a.norm_sq()).sub(&a.prod_scalar(b.norm_sq()));
98        Ok(self.m.add(&v.prod_cross(&n).prod_scalar(k)))
99    }
100
101    /// `cos` and `sin` of half the sweep from `p` to `q` along its circle,
102    /// the way it turns — two points of the arc, `q` after `p`: the part of
103    /// the arc between two of its points. The sine is the chord over the
104    /// diameter, `|q - p| / (2 r)`; the cosine is the sine of the whole
105    /// sweep over twice the sine of half of it,
106    /// `n . ((p - c) x (q - c)) / (r |q - p|)`, which carries the sign of a
107    /// sweep past a half turn without a case.
108    pub fn half_sweep_between(&self, p: &Vector3<T>, q: &Vector3<T>) -> GeopResult<(T, T)> {
109        let (c, r, n) = (self.center()?, self.radius()?, self.normal()?);
110        let chord = q.sub(p).try_norm()?;
111        let sin = chord.div(T::TWO.mul(r))?;
112        let turned = p.sub(&c).prod_cross(&q.sub(&c));
113        let cos = n.prod_dot(&turned).div(r.mul(chord))?;
114        Ok((cos, sin))
115    }
116
117    /// Its unit tangent at `end`, in its own direction — from `s` to `e`:
118    /// the chord's direction turned back (at the start) or on (at the end)
119    /// by half the sweep, in the arc's plane.
120    pub fn tangent(&self, at_end: bool) -> GeopResult<Vector3<T>> {
121        let chord = self.e.sub(&self.s).normalize()?;
122        let side = self.normal()?.prod_cross(&chord);
123        let (cos, sin) = self.half_sweep()?;
124        let sin = if at_end { sin } else { sin.neg() };
125        Ok(chord.prod_scalar(cos).add(&side.prod_scalar(sin)))
126    }
127}
128
129/// `(u, v)`: two unit vectors that span, with `direction`, all of space —
130/// to measure how far another vector is from parallel to it along. Which
131/// two is a free choice (see [`Vector3::across`]): worked out in plain
132/// numbers, sharp.
133pub fn across(direction: &Vector3) -> GeopResult<[Vector3; 2]> {
134    let sharp = |v: &Vector3| v.map(|c| Scalar::from_f64(c.to_f64()));
135    let [u, v] = sharp(direction).normalize()?.across()?;
136    Ok([sharp(&u), sharp(&v)])
137}
geop_core_sketch/space/geometry.rs

geop_core_sketch/space/
geometry.rs
