ogeom_algo/thickness/shape.rs
1//! Wall thickness against a shape's own surfaces.
2
3use std::sync::OnceLock;
4
5use ogeom_core::{OgeomResult, parallel};
6use ogeom_math::{Point, Vector};
7use ogeom_topo::{Model, Shape};
8
9use super::{Method, Thickness, inward, meets, spread, through};
10use crate::classify::{Containment, SolidBoundary, exact_ring_chord};
11use crate::distance::{Passed, Target};
12
13/// The wall thickness at each point of `shape`, a solid or a shell,
14/// measured into the material against the point's outward normal on the
15/// faces' own surfaces; in the points' order, `None` where it is thicker
16/// than [`Thickness::max`] or the inward ray leaves the shape without
17/// meeting a face.
18///
19/// A ray meets a face where it crosses the face's surface inside its trim,
20/// or within the band of the trimming rings drawn at ten thousand times the
21/// confusion distance, as [`classify_in_solid_exact`](crate::classify_in_solid_exact)
22/// draws them. A sphere's nearest surface is
23/// [`to_shape`](crate::distance::to_shape)'s nearest point, among the
24/// faces, edges and vertices not within the facing limit of edge-on to the
25/// point's normal. See the [module](super) for the rest.
26///
27/// A point that is not finite, or whose normal has no direction, answers
28/// `None`.
29///
30/// # Errors
31///
32/// As [`Thickness::check`];
33/// [`OgeomError::Dimension`](ogeom_core::OgeomError::Dimension) if there is
34/// not one normal per point;
35/// [`OgeomError::Construction`](ogeom_core::OgeomError::Construction) if
36/// the shape is not a solid, a shell or a compound of solids; and as the
37/// intersections and projections report.
38pub fn thickness_of_shape(
39 model: &Model,
40 shape: &Shape,
41 points: &[Point],
42 normals: &[Vector],
43 options: &Thickness,
44) -> OgeomResult<Vec<Option<f64>>> {
45 let sin = options.checked(points.len(), normals.len())?;
46 let tol = options.tol;
47 let boundary = SolidBoundary::of_closed(model, shape, exact_ring_chord(tol), tol)?;
48 let walls = Walls {
49 model,
50 shape,
51 boundary: &boundary,
52 target: OnceLock::new(),
53 sin,
54 options,
55 };
56 let answers = parallel::map_ordered(points, |i, &p| match inward(p, normals[i]) {
57 Some(u) => walls.at(p, u),
58 None => Ok(None),
59 });
60 answers.into_iter().collect()
61}
62
63/// A shape measured against, with what every point asks of it.
64struct Walls<'a> {
65 model: &'a Model,
66 shape: &'a Shape,
67 boundary: &'a SolidBoundary,
68 /// The shape's elements, for nearest points, gathered when first asked.
69 target: OnceLock<OgeomResult<Target>>,
70 sin: f64,
71 options: &'a Thickness,
72}
73
74impl Walls<'_> {
75 fn target(&self) -> OgeomResult<&Target> {
76 let tol = self.options.tol;
77 match self
78 .target
79 .get_or_init(|| Target::of(self.model, self.shape, false, tol))
80 {
81 Ok(target) => Ok(target),
82 Err(e) => Err(e.clone()),
83 }
84 }
85
86 /// The faces `p` lies on: those that hold it, and where one holds it
87 /// in the band of its trim, those that hold its nearest point on the
88 /// shape. A point in the band is on the face only to within the band,
89 /// and as much on the faces meeting it there, which hold the point of
90 /// the shape it stands for even where it stands off their surfaces by
91 /// more than their tolerance.
92 fn own(&self, p: Point) -> OgeomResult<Vec<Shape>> {
93 let tol = self.options.tol;
94 let held = self.boundary.faces_holding(self.model, p, tol)?;
95 let banded = held.iter().any(|(_, at)| *at == Containment::On);
96 let mut on: Vec<Shape> = held.into_iter().map(|(face, _)| face).collect();
97 if banded
98 && let Some(found) =
99 self.target()?
100 .nearest(self.model, p, None, None, f64::INFINITY, tol)?
101 {
102 for face in self.boundary.faces_at(self.model, found.point, tol)? {
103 if !on.iter().any(|held| held.is_same(&face)) {
104 on.push(face);
105 }
106 }
107 }
108 Ok(on)
109 }
110
111 /// The wall at `p`, inward along the unit `u`.
112 fn at(&self, p: Point, u: Vector) -> OgeomResult<Option<f64>> {
113 let on = self.own(p)?;
114 let max = self.options.max;
115 match self.options.method {
116 Method::Ray { cone, rays } => {
117 let mut least: Option<f64> = None;
118 for k in 0..rays {
119 let v = spread(u, cone, rays, k);
120 if let Some(length) = self.ray(p, v, max, &on)?
121 && least.is_none_or(|held| length < held)
122 {
123 least = Some(length);
124 }
125 }
126 Ok(least)
127 }
128 Method::Sphere => {
129 let ray = self.ray(p, u, max, &on)?;
130 self.sphere(p, u, ray, &on)
131 }
132 }
133 }
134
135 fn ray(&self, p: Point, v: Vector, max: f64, on: &[Shape]) -> OgeomResult<Option<f64>> {
136 self.boundary
137 .first_crossing(self.model, (p, v), max, on, self.sin, self.options.tol)
138 }
139
140 /// The diameter of the medial ball at `p` along the unit `u`, shrunk
141 /// from half the inward ray's reading `ray`, or from half the thickest
142 /// wall where the ray met nothing within it.
143 ///
144 /// Each step finds the surface nearest the centre within the ball and
145 /// shrinks to the ball through it, until the nearest surface is within
146 /// the confusion distance of the ball. Of the faces `p` lies on, `on`,
147 /// and the edges and vertices they alone hold, one whose nearest point
148 /// lies on the ball's radius from `p` is passed over: there the ray
149 /// along it passed over it too, as `p`'s own at its start, or edge-on
150 /// where it holds the centre.
151 fn sphere(
152 &self,
153 p: Point,
154 u: Vector,
155 ray: Option<f64>,
156 on: &[Shape],
157 ) -> OgeomResult<Option<f64>> {
158 let target = self.target()?;
159 let tol = self.options.tol;
160 let mut radius = match ray {
161 Some(length) => length / 2.0,
162 None if self.options.max.is_finite() => self.options.max / 2.0,
163 None => return Ok(None),
164 };
165 let mut bounded = ray.is_some();
166 let faces = |n: Vector| meets(n, u, self.sin);
167 let passed = Passed { faces: on, from: p };
168 loop {
169 let centre = p + u * radius;
170 let Some(found) =
171 target.nearest(self.model, centre, Some(&faces), Some(&passed), radius, tol)?
172 else {
173 break;
174 };
175 bounded = true;
176 if radius - found.distance <= tol.confusion() {
177 break;
178 }
179 match through(p, u, found.point) {
180 Some(r) if r < radius => radius = r,
181 _ => break,
182 }
183 }
184 Ok(bounded.then_some(2.0 * radius))
185 }
186}