Another project
1use bone_types::{
2 Aabb3, Angle, AngleTolerance, ChordHeightTolerance, Length, Parameter, Plane3, Point3,
3 Tolerance, UnitVec3, Vec3,
4};
5use core::f64::consts::TAU;
6use uom::si::angle::radian;
7use uom::si::length::millimeter;
8
9use crate::KernelError;
10use crate::angles;
11use crate::circular3;
12use crate::mesh::{MeshVertex, TriMesh};
13use crate::surface3::Surface3;
14
15const MIN_SEGMENTS: u32 = 3;
16
17#[derive(Copy, Clone, Debug, PartialEq)]
18pub struct CylinderSurface {
19 plane: Plane3,
20 radius: Length,
21 start_angle: Angle,
22 sweep_angle: Angle,
23 height: Length,
24}
25
26impl CylinderSurface {
27 pub fn new(
28 plane: Plane3,
29 radius: Length,
30 start_angle: Angle,
31 sweep_angle: Angle,
32 height: Length,
33 tolerance: Tolerance,
34 ) -> Result<Self, KernelError> {
35 let r = radius.get::<millimeter>();
36 let h = height.get::<millimeter>();
37 if !r.is_finite() || !h.is_finite() || r < tolerance.value() || h < tolerance.value() {
38 return Err(KernelError::DegenerateCylinder);
39 }
40 let sweep = sweep_angle.get::<radian>().abs();
41 let angle_eps = AngleTolerance::from_arc_length(tolerance, radius).radians();
42 if !sweep.is_finite() || sweep < angle_eps || sweep > TAU + angle_eps {
43 return Err(KernelError::DegenerateCylinder);
44 }
45 Ok(Self {
46 plane,
47 radius,
48 start_angle,
49 sweep_angle,
50 height,
51 })
52 }
53
54 #[must_use]
55 pub const fn plane(self) -> Plane3 {
56 self.plane
57 }
58
59 #[must_use]
60 pub fn center(self) -> Point3 {
61 self.plane.origin()
62 }
63
64 #[must_use]
65 pub fn axis(self) -> UnitVec3 {
66 self.plane.normal()
67 }
68
69 #[must_use]
70 pub const fn radius(self) -> Length {
71 self.radius
72 }
73
74 #[must_use]
75 pub const fn start_angle(self) -> Angle {
76 self.start_angle
77 }
78
79 #[must_use]
80 pub const fn sweep_angle(self) -> Angle {
81 self.sweep_angle
82 }
83
84 #[must_use]
85 pub const fn height(self) -> Length {
86 self.height
87 }
88
89 #[must_use]
90 pub fn radius_mm(self) -> f64 {
91 self.radius.get::<millimeter>()
92 }
93
94 #[must_use]
95 pub fn height_mm(self) -> f64 {
96 self.height.get::<millimeter>()
97 }
98
99 #[must_use]
100 pub fn start_rad(self) -> f64 {
101 self.start_angle.get::<radian>()
102 }
103
104 #[must_use]
105 pub fn sweep_rad(self) -> f64 {
106 self.sweep_angle.get::<radian>()
107 }
108
109 #[must_use]
110 fn angle_at(self, u: f64) -> Angle {
111 Angle::new::<radian>(self.start_rad() + self.sweep_rad() * u)
112 }
113
114 #[must_use]
115 fn point_at(self, u: f64, v: f64) -> Point3 {
116 let base = circular3::point_at(self.plane, self.radius, self.angle_at(u));
117 base + self.plane.normal().into_vec(self.height * v)
118 }
119
120 #[must_use]
121 fn normal_at(self, u: f64) -> UnitVec3 {
122 let theta = self.start_rad() + self.sweep_rad() * u;
123 let (xx, xy, xz) = self.plane.x_axis().components();
124 let (yx, yy, yz) = self.plane.y_axis().components();
125 let (c, s) = (theta.cos(), theta.sin());
126 let radial = UnitVec3::new_unchecked(c * xx + s * yx, c * xy + s * yy, c * xz + s * yz);
127 if self.sweep_rad() < 0.0 {
128 radial.reversed()
129 } else {
130 radial
131 }
132 }
133}
134
135impl Surface3 for CylinderSurface {
136 fn evaluate(&self, u: Parameter, v: Parameter) -> Point3 {
137 self.point_at(u.value(), v.value())
138 }
139
140 fn partials(&self, u: Parameter, _v: Parameter) -> (Vec3, Vec3) {
141 let tangent = circular3::tangent_vec(self.plane, self.radius, self.angle_at(u.value()))
142 * self.sweep_rad();
143 let axial = self.plane.normal().into_vec(self.height);
144 (tangent, axial)
145 }
146
147 fn normal(&self, u: Parameter, _v: Parameter) -> UnitVec3 {
148 self.normal_at(u.value())
149 }
150
151 fn bounding_box(&self) -> Aabb3 {
152 let base =
153 circular3::bounding_box(self.plane, self.radius, self.start_angle, self.sweep_angle);
154 let top_origin = self.plane.origin() + self.plane.normal().into_vec(self.height);
155 let top_plane = Plane3::new_unchecked(top_origin, self.plane.x_axis(), self.plane.y_axis());
156 let top =
157 circular3::bounding_box(top_plane, self.radius, self.start_angle, self.sweep_angle);
158 base.union(top)
159 }
160
161 fn contains_point(&self, p: Point3, tolerance: Tolerance) -> bool {
162 let (px, py, pn) = circular3::local_coords(self.plane, p);
163 let radial = (px * px + py * py).sqrt();
164 if (radial - self.radius_mm()).abs() > tolerance.value() {
165 return false;
166 }
167 if pn < -tolerance.value() || pn > self.height_mm() + tolerance.value() {
168 return false;
169 }
170 let angle_eps = AngleTolerance::from_arc_length(tolerance, self.radius).radians();
171 angles::contains(py.atan2(px), self.start_rad(), self.sweep_rad(), angle_eps)
172 }
173
174 fn tessellate(&self, chord: ChordHeightTolerance, angle: AngleTolerance) -> TriMesh {
175 let n = circular3::sample_count(self.radius, self.sweep_angle, chord, MIN_SEGMENTS).max(
176 circular3::angular_sample_count(self.sweep_angle, angle, MIN_SEGMENTS),
177 );
178 let inv_n = 1.0 / f64::from(n);
179 TriMesh::from_grid(n + 1, 2, |i, j| {
180 let u = f64::from(i) * inv_n;
181 MeshVertex::new(self.point_at(u, f64::from(j)), self.normal_at(u))
182 })
183 }
184}
185
186impl core::fmt::Display for CylinderSurface {
187 fn fmt(&self, f: &mut core::fmt::Formatter<'_>) -> core::fmt::Result {
188 write!(
189 f,
190 "cylinder_surface{{ c={}, r={} mm, start={} rad, sweep={} rad, h={} mm, axis={} }}",
191 self.center(),
192 self.radius_mm(),
193 self.start_rad(),
194 self.sweep_rad(),
195 self.height_mm(),
196 self.axis(),
197 )
198 }
199}