geop_core_geometry/nurb_curve/
project.rs

1use geop_core_math::{
2    disjoint_set::DisjointSet, geop_error::GeopResult, scalars::Scalar, vector::Vector3,
3};
4
5use super::NurbCurve;
6
7impl<S: Scalar> NurbCurve<S, 4> {
8    /// Every parameter `t` at which this curve passes through `target`
9    /// (plural since a self-intersecting curve can pass through the same
10    /// point at more than one, genuinely distinct, parameter).
11    ///
12    /// Recursively subdivides the curve, discarding any segment whose
13    /// convex hull could not contain `target`. A surviving segment
14    /// converges once its hull's chord length is no longer definitely
15    /// greater than `min_subdivision_size`, contributing the
16    /// [`Scalar::union`] of its own `[t0, t1]` domain as a candidate,
17    /// folded into a [`DisjointSet`] so no two returned solutions ever
18    /// describe the same physical parameter — for a genuine interval
19    /// scalar each is a real "the true parameter is provably within this
20    /// span" guarantee, not an arbitrarily narrowed single point.
21    pub fn project(
22        &self,
23        target: Vector3<S>,
24        max_nodes: usize,
25        min_subdivision_size: S,
26    ) -> GeopResult<Vec<S>> {
27        let mut stack: Vec<NurbCurve<S, 4>> = vec![self.clone()];
28        let mut explored = 0usize;
29        let mut solutions: DisjointSet<S> = DisjointSet::new();
30
31        while let Some(seg) = stack.pop() {
32            if explored >= max_nodes {
33                break;
34            }
35            explored += 1;
36
37            let hull = match seg.convex_hull() {
38                Ok(hull) => hull,
39                // Degenerate segment (zero weight): can't bound or split it — skip it.
40                Err(_) => continue,
41            };
42            if hull.definitely_not_contains(&target) {
43                continue;
44            }
45
46            let (t0, t1) = seg.domain();
47            let chord_len = hull.points[hull.points.len() - 1]
48                .sub(&hull.points[0])
49                .norm();
50            if !chord_len.definitely_greater(min_subdivision_size) {
51                solutions.insert(t0.union(t1));
52                continue;
53            }
54
55            if let Ok((left, right)) = seg.split_mid() {
56                stack.push(left);
57                stack.push(right);
58            }
59            // Cannot split (e.g. midpoint already at multiplicity p+1) and
60            // hasn't converged — nothing more to do with this segment.
61        }
62
63        Ok(solutions.into_vec())
64    }
65}
66
67#[cfg(test)]
68mod tests {
69    use geop_core_math::for_all_scalars;
70    use geop_core_math::{scalars::Scalar, vector::Vector3};
71
72    use super::super::NurbCurve3D;
73
74    const MAX: usize = 500;
75
76    fn line<S: Scalar>() -> NurbCurve3D<S> {
77        let p = |x: f64, y: f64, z: f64| {
78            geop_core_math::vector::Vector4::from_array([
79                S::from_f64(x),
80                S::from_f64(y),
81                S::from_f64(z),
82                S::ONE,
83            ])
84        };
85        NurbCurve3D::try_new(
86            1,
87            vec![p(0.0, 0.0, 0.0), p(1.0, 0.0, 0.0)],
88            vec![S::ZERO, S::ZERO, S::ONE, S::ONE],
89        )
90        .unwrap()
91    }
92
93    // `ScalInFPA64` is fixed-point at 2^-32 (~2.3e-10) resolution — once
94    // subdivision reaches that granularity, `t` can't be refined any
95    // further, so the tolerance here has to be a bit looser than
96    // `min_subdivision_size` itself to accommodate that scalar's precision
97    // floor.
98    const TOL: f64 = 1e-4;
99
100    fn check_project_point_on_line<S: Scalar>() {
101        let curve = line::<S>();
102        let target = Vector3::from_array([S::from_f64(0.42), S::from_f64(0.0), S::from_f64(0.0)]);
103        let solutions = curve.project(target, MAX, S::from_f64(1e-6)).unwrap();
104        assert_eq!(solutions.len(), 1);
105        assert!(
106            solutions[0]
107                .sub(S::from_f64(0.42))
108                .abs()
109                .could_be_less(S::from_f64(TOL))
110        );
111    }
112    #[test]
113    fn project_point_on_line() {
114        for_all_scalars!(check_project_point_on_line);
115    }
116
117    fn check_project_point_off_line_finds_nothing<S: Scalar>() {
118        let curve = line::<S>();
119        // Off the line entirely (y=1) — no segment's hull could ever
120        // contain it.
121        let target = Vector3::from_array([S::from_f64(0.3), S::from_f64(1.0), S::ZERO]);
122        assert!(
123            curve
124                .project(target, MAX, S::from_f64(1e-6))
125                .unwrap()
126                .is_empty()
127        );
128    }
129    #[test]
130    fn project_point_off_line_finds_nothing() {
131        for_all_scalars!(check_project_point_off_line_finds_nothing);
132    }
133
134    fn check_zero_budget_finds_nothing<S: Scalar>() {
135        let curve = line::<S>();
136        let target = Vector3::from_array([S::from_f64(0.5), S::ZERO, S::ZERO]);
137        assert!(
138            curve
139                .project(target, 0, S::from_f64(1e-6))
140                .unwrap()
141                .is_empty()
142        );
143    }
144    #[test]
145    fn zero_budget_finds_nothing() {
146        for_all_scalars!(check_zero_budget_finds_nothing);
147    }
148}
geop_core_geometry/nurb_curve/project.rs

geop_core_geometry/nurb_curve/
project.rs
