axiolid_nurbs/
curve_analysis.rs

1//! Differential properties of regular parametric curves.
2
3use axiolid_contracts::{GeomError, GeomResult};
4use axiolid_core::{Point2, Point3, Scalar, Tolerance, Vec2, Vec3};
5use axiolid_curve::{Curve2, Curve3};
6use axiolid_evaluate::curve::{jet2, jet3};
7
8/// Parameter-invariant differential properties of a planar curve.
9#[derive(Debug, Clone, Copy, PartialEq)]
10pub struct CurveDifferential2 {
11    /// Evaluated point.
12    pub point: Point2,
13    /// Unit tangent in increasing-parameter direction.
14    pub unit_tangent: Vec2,
15    /// Unsigned parameter-invariant curvature.
16    pub curvature: Scalar,
17    /// Signed planar curvature; positive bends toward the tangent's left normal.
18    pub signed_curvature: Scalar,
19    /// Left unit normal scaled by signed curvature.
20    pub curvature_vector: Vec2,
21}
22
23/// Parameter-invariant differential properties of a spatial curve.
24#[derive(Debug, Clone, Copy, PartialEq)]
25pub struct CurveDifferential3 {
26    /// Evaluated point.
27    pub point: Point3,
28    /// Unit tangent in increasing-parameter direction.
29    pub unit_tangent: Vec3,
30    /// Unsigned parameter-invariant curvature.
31    pub curvature: Scalar,
32    /// Normal-acceleration direction scaled by curvature.
33    pub curvature_vector: Vec3,
34}
35
36/// Analyze a regular planar curve at `t`.
37pub fn analyze_curve2(
38    curve: &Curve2,
39    t: Scalar,
40    tolerance: Tolerance,
41) -> GeomResult<CurveDifferential2> {
42    reject_polyline_corner2(curve, t)?;
43    let j = jet2(curve, t)?;
44    let speed = j.first.length();
45    regular(speed, tolerance)?;
46    let signed = j.first.perp_dot(j.second) / speed.powi(3);
47    let tangent = j.first / speed;
48    let normal = Vec2::new(-tangent.y, tangent.x);
49    finite(signed)?;
50    Ok(CurveDifferential2 {
51        point: j.point,
52        unit_tangent: tangent,
53        curvature: signed.abs(),
54        signed_curvature: signed,
55        curvature_vector: normal * signed,
56    })
57}
58
59/// Analyze a regular spatial curve at `t`.
60pub fn analyze_curve3(
61    curve: &Curve3,
62    t: Scalar,
63    tolerance: Tolerance,
64) -> GeomResult<CurveDifferential3> {
65    reject_polyline_corner3(curve, t)?;
66    let j = jet3(curve, t)?;
67    let speed = j.first.length();
68    regular(speed, tolerance)?;
69    let tangent = j.first / speed;
70    let normal_acceleration = j.second - tangent * j.second.dot(tangent);
71    let curvature_vector = normal_acceleration / speed.powi(2);
72    let curvature = curvature_vector.length();
73    finite(curvature)?;
74    Ok(CurveDifferential3 {
75        point: j.point,
76        unit_tangent: tangent,
77        curvature,
78        curvature_vector,
79    })
80}
81
82fn reject_polyline_corner2(curve: &Curve2, t: Scalar) -> GeomResult<()> {
83    if let Curve2::Polyline(polyline) = curve {
84        reject_corner(polyline.points.len(), polyline.closed, t)?;
85    }
86    Ok(())
87}
88
89fn reject_polyline_corner3(curve: &Curve3, t: Scalar) -> GeomResult<()> {
90    if let Curve3::Polyline(polyline) = curve {
91        reject_corner(polyline.points.len(), polyline.closed, t)?;
92    }
93    Ok(())
94}
95
96fn reject_corner(point_count: usize, closed: bool, t: Scalar) -> GeomResult<()> {
97    if !t.is_finite() {
98        return Ok(());
99    }
100    let segment_count = if closed {
101        point_count
102    } else {
103        point_count.saturating_sub(1)
104    };
105    if segment_count == 0 {
106        return Ok(());
107    }
108    let end = segment_count as Scalar;
109    let parameter = t.clamp(0.0, end);
110    let is_vertex = parameter.fract() == 0.0;
111    let is_two_sided = closed || (parameter > 0.0 && parameter < end);
112    if is_vertex && is_two_sided {
113        Err(GeomError::Degenerate(
114            "polyline curvature is undefined at a non-smooth vertex".to_owned(),
115        ))
116    } else {
117        Ok(())
118    }
119}
120
121fn regular(speed: Scalar, tolerance: Tolerance) -> GeomResult<()> {
122    if !speed.is_finite() || speed <= tolerance.linear() {
123        Err(GeomError::Degenerate(format!(
124            "curve speed {speed} does not exceed the linear tolerance"
125        )))
126    } else {
127        Ok(())
128    }
129}
130
131fn finite(value: Scalar) -> GeomResult<()> {
132    if value.is_finite() {
133        Ok(())
134    } else {
135        Err(GeomError::Degenerate(
136            "curve differential is non-finite".to_owned(),
137        ))
138    }
139}