mirror of
https://github.com/bherbruck/solidcraft.git
synced 2026-11-06 03:51:08 +00:00
kernel: plane_fit fits and certifies planes on the schema 4 record
plane_through gives a plane record (its natural normal turning with the loop); EdgeSampler samples a record edge through CurveEval in the direction a loop runs along it (an intersection curve traced as the record holds it, with the evaluator's own hint); planar_face certifies a loop of record edges from their records (lines by their ends, circles and ellipses by their frames, polylines and splines by their points with positive weights, intersection curves by a plane among their surfaces) and returns the plane, the caller building the face on the record. (not built: base in progress) Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_017itGh5fpMEousnZ91cYcVt
This commit is contained in:
jyviko
@@ -1,40 +1,26 @@
|
||||
//! Planes fitted to points around a loop (Newell's method), for planar faces of polygons read
|
||||
//! from files and for patches filling a planar loop of edges, and points sampled along edges
|
||||
//! (shared with the fill of a loop that is not planar).
|
||||
//! from files and for patches filling a planar loop of edges, and points sampled along record
|
||||
//! edges (shared with the fill of a loop that is not planar).
|
||||
|
||||
use mt::{BoundedCurve, EuclideanSpace, InnerSpace, MetricSpace, ParametricCurve};
|
||||
use solidcraft_geom::Vec3;
|
||||
use truck_modeling as mt;
|
||||
use solidcraft_geom::brep::eval::CurveEval;
|
||||
use solidcraft_geom::brep::{Brep, Curve, Edge, Surface};
|
||||
|
||||
use crate::body::{from_p3, p3};
|
||||
use crate::{KernelError, Result};
|
||||
use crate::Result;
|
||||
|
||||
/// The plane through points round a loop: its normal by Newell's method (turning with the
|
||||
/// loop), through their centroid; `None` for fewer than three points or no area.
|
||||
pub(crate) fn plane_through(pts: &[mt::Point3]) -> Option<mt::Plane> {
|
||||
let n = pts.len();
|
||||
if n < 3 {
|
||||
/// loop, the plane's natural normal), through their centroid; `None` for fewer than three points
|
||||
/// or no area.
|
||||
pub(crate) fn plane_through(pts: &[Vec3]) -> Option<Surface> {
|
||||
if pts.len() < 3 {
|
||||
return None;
|
||||
}
|
||||
let mut nrm = mt::Vector3::new(0.0, 0.0, 0.0);
|
||||
let mut c = mt::Vector3::new(0.0, 0.0, 0.0);
|
||||
for i in 0..n {
|
||||
let (a, b) = (pts.get(i)?, pts.get((i + 1) % n)?);
|
||||
nrm.x += (a.y - b.y) * (a.z + b.z);
|
||||
nrm.y += (a.z - b.z) * (a.x + b.x);
|
||||
nrm.z += (a.x - b.x) * (a.y + b.y);
|
||||
c += a.to_vec();
|
||||
}
|
||||
let len = nrm.magnitude();
|
||||
if !(len > 1e-12) {
|
||||
return None;
|
||||
}
|
||||
let z = nrm / len;
|
||||
let o = mt::Point3::from_vec(c / n as f64);
|
||||
let a = if z.x.abs() < 0.9 { mt::Vector3::unit_x() } else { mt::Vector3::unit_y() };
|
||||
let x = (a - z * a.dot(z)).normalize();
|
||||
let z = loop_area(pts).normalized()?;
|
||||
let o = pts.iter().fold(Vec3::ZERO, |a, p| a + *p) / pts.len() as f64;
|
||||
let a = if z.x.abs() < 0.9 { Vec3::X } else { Vec3::Y };
|
||||
let x = (a - z * a.dot(z)).normalized()?;
|
||||
let y = z.cross(x);
|
||||
Some(mt::Plane::new(o, o + x, o + y))
|
||||
Some(Surface::Plane { origin: [o.x, o.y, o.z], u_axis: [x.x, x.y, x.z], v_axis: [y.x, y.y, y.z] })
|
||||
}
|
||||
|
||||
/// The area vector of a closed polygon (Newell's method): normal to the loop, as long as the area
|
||||
@@ -48,97 +34,88 @@ pub(crate) fn loop_area(points: &[Vec3]) -> Vec3 {
|
||||
n * 0.5
|
||||
}
|
||||
|
||||
/// The point of `c` at `t` (`None` where an intersection curve cannot find it).
|
||||
fn curve_point(c: &mt::Curve, t: f64) -> Option<Vec3> {
|
||||
let p = match c {
|
||||
mt::Curve::IntersectionCurve(ic) => ic.search_triple(t)?.0,
|
||||
other => other.subs(t),
|
||||
};
|
||||
let p = from_p3(p);
|
||||
p.is_finite().then_some(p)
|
||||
}
|
||||
|
||||
/// Points along an edge, in its own direction, at `s ∈ [0, 1]` of its parameter range. An
|
||||
/// intersection curve's point is searched from the surface parameters of the point found before
|
||||
/// when that one is near (an eighth of the range), and kept only when it lies near the curve's
|
||||
/// leader (a seeded search can slide to another branch of the surfaces' intersection); else the
|
||||
/// unseeded search decides. `None` where no point is found.
|
||||
/// Points along a record edge, in the direction a loop runs along it, at `s ∈ [0, 1]` of its
|
||||
/// interval: the record curve's own points (an intersection curve traced as the record holds it).
|
||||
pub(crate) struct EdgeSampler {
|
||||
curve: mt::Curve,
|
||||
curve: CurveEval,
|
||||
range: (f64, f64),
|
||||
hint: Option<(f64, mt::Point2, mt::Point2)>,
|
||||
}
|
||||
|
||||
impl EdgeSampler {
|
||||
pub(crate) fn new(edge: &mt::Edge) -> Option<EdgeSampler> {
|
||||
let curve = edge.oriented_curve();
|
||||
let range = curve.range_tuple();
|
||||
(range.0.is_finite() && range.1.is_finite() && range.0 < range.1).then_some(EdgeSampler { curve, range, hint: None })
|
||||
/// The sampler of `edge` run from its start (`forward`) or from its end.
|
||||
pub(crate) fn new(edge: &Edge, forward: bool) -> Option<EdgeSampler> {
|
||||
let [t0, t1] = edge.interval;
|
||||
let curve = CurveEval::new(&edge.curve).ok()?;
|
||||
(t0.is_finite() && t1.is_finite() && t0 < t1).then_some(EdgeSampler { curve, range: if forward { (t0, t1) } else { (t1, t0) } })
|
||||
}
|
||||
|
||||
pub(crate) fn at(&mut self, s: f64) -> Option<Vec3> {
|
||||
let (t0, t1) = self.range;
|
||||
let t = if s >= 1.0 { t1 } else { t0 + (t1 - t0) * s };
|
||||
let mt::Curve::IntersectionCurve(ic) = &self.curve else { return curve_point(&self.curve, t) };
|
||||
let leader = ic.leader();
|
||||
let (p, d) = (leader.subs(t), leader.der(t));
|
||||
let near = |q: mt::Point3| q.distance(p) <= 4.0 * ic.tolerance() + 1e-6;
|
||||
let seeded = self
|
||||
.hint
|
||||
.filter(|(tp, _, _)| (t - tp).abs() <= (t1 - t0) / 8.0)
|
||||
.and_then(|(_, a, b)| mt::double_projection(ic.surface0(), Some((a.x, a.y)), ic.surface1(), Some((b.x, b.y)), p, d, 100))
|
||||
.filter(|f| near(f.0));
|
||||
let found = seeded.or_else(|| mt::double_projection(ic.surface0(), None, ic.surface1(), None, p, d, 100))?;
|
||||
self.hint = Some((t, found.1, found.2));
|
||||
let q = from_p3(found.0);
|
||||
q.is_finite().then_some(q)
|
||||
self.curve.point(t).filter(|p| p.is_finite())
|
||||
}
|
||||
}
|
||||
|
||||
/// The planar face bounded by the closed loop `chain` (each edge oriented along it) when its
|
||||
/// edges are certainly in one plane, `None` otherwise (also when an edge cannot be sampled).
|
||||
/// The plane is fitted to points along the edges ([`plane_through`], the normal turning with the
|
||||
/// loop); the edges are certified within `tol` of it conservatively: a line by its ends, a
|
||||
/// B-spline or NURBS curve with positive weights by its control points (the curve lies in their
|
||||
/// convex hull), an intersection curve only when one of its two surfaces is a plane that the
|
||||
/// loop's points lie on and the fitted plane agrees with over the loop (the curve lies in that
|
||||
/// plane). Samples alone are never trusted: a spline can bulge out of the plane between them.
|
||||
pub(crate) fn planar_face(chain: &[mt::Edge], tol: f64) -> Result<Option<mt::Face>> {
|
||||
/// The plane of the closed loop `chain` of edges of `b` (each with the way the loop runs along
|
||||
/// it) when its edges are certainly in one plane, `None` otherwise (also when an edge cannot be
|
||||
/// sampled). The plane is fitted to points along the edges ([`plane_through`], its natural normal
|
||||
/// turning with the loop); the edges are certified within `tol` of it from their records: a line
|
||||
/// by its ends, a circle or ellipse by its frame, a polyline, B-spline or NURBS curve with
|
||||
/// positive weights by its points (the curve lies in their convex hull), an intersection curve
|
||||
/// only when one of its two surfaces is a plane that the loop's points lie on (the curve lies in
|
||||
/// that plane). Samples alone are never trusted: a spline can bulge out of the plane between them.
|
||||
pub(crate) fn planar_face(b: &Brep, chain: &[(u32, bool)], tol: f64) -> Result<Option<Surface>> {
|
||||
const SAMPLES: usize = 16;
|
||||
let edges: Vec<&Edge> = match chain.iter().map(|(e, _)| b.edges.get(*e as usize)).collect::<Option<_>>() {
|
||||
Some(e) => e,
|
||||
None => return Ok(None),
|
||||
};
|
||||
let mut pts = Vec::with_capacity(chain.len() * SAMPLES);
|
||||
for e in chain {
|
||||
let Some(mut sampler) = EdgeSampler::new(e) else { return Ok(None) };
|
||||
for (e, (_, forward)) in edges.iter().zip(chain) {
|
||||
let Some(mut sampler) = EdgeSampler::new(e, *forward) else { return Ok(None) };
|
||||
for k in 0..SAMPLES {
|
||||
let Some(p) = sampler.at(k as f64 / SAMPLES as f64) else { return Ok(None) };
|
||||
pts.push(p3(p));
|
||||
pts.push(p);
|
||||
}
|
||||
}
|
||||
let Some(plane) = plane_through(&pts) else { return Ok(None) };
|
||||
let (o, n) = (plane.origin(), plane.normal());
|
||||
let off = |p: mt::Point3, o: mt::Point3, n: mt::Vector3| (p - o).dot(n).abs();
|
||||
let on_plane = |p: mt::Point3| off(p, o, n) <= tol;
|
||||
let Some((o, n)) = plane_frame(&plane) else { return Ok(None) };
|
||||
let on = |p: Vec3, o: Vec3, n: Vec3| (p - o).dot(n).abs() <= tol;
|
||||
let on_plane = |p: Vec3| on(p, o, n);
|
||||
if !pts.iter().all(|p| on_plane(*p)) {
|
||||
return Ok(None);
|
||||
}
|
||||
// A support plane `q` of an intersection curve holds the loop when the loop's points lie on
|
||||
// it and their feet on the fitted plane do too.
|
||||
let holds = |q: &mt::Plane| {
|
||||
let (qo, qn) = (q.origin(), q.normal());
|
||||
pts.iter().all(|p| off(*p, qo, qn) <= tol && off(*p - n * (*p - o).dot(n), qo, qn) <= tol)
|
||||
};
|
||||
for e in chain {
|
||||
let certified = match e.curve() {
|
||||
mt::Curve::Line(l) => on_plane(l.0) && on_plane(l.1),
|
||||
mt::Curve::BSplineCurve(b) => b.control_points().iter().all(|p| on_plane(*p)),
|
||||
mt::Curve::NurbsCurve(c) => {
|
||||
c.non_rationalized().control_points().iter().all(|h| h.w > 0.0 && on_plane(mt::Point3::new(h.x / h.w, h.y / h.w, h.z / h.w)))
|
||||
for e in &edges {
|
||||
let certified = match &e.curve {
|
||||
Curve::Line { start, end } => on_plane(Vec3::from(*start)) && on_plane(Vec3::from(*end)),
|
||||
Curve::Circle { frame, .. } | Curve::Ellipse { frame, .. } => {
|
||||
let z = Vec3::from(frame.x_axis).cross(Vec3::from(frame.y_axis));
|
||||
on_plane(Vec3::from(frame.origin)) && z.cross(n).len() <= 1e-9
|
||||
}
|
||||
Curve::Polyline { points } | Curve::Spline { points, .. } => points.iter().all(|p| on_plane(Vec3::from(*p))),
|
||||
Curve::RationalSpline { homogeneous, .. } => {
|
||||
homogeneous.iter().all(|h| h[3] > 0.0 && on_plane(Vec3::new(h[0] / h[3], h[1] / h[3], h[2] / h[3])))
|
||||
}
|
||||
// A support plane holds the loop when the loop's points lie on it.
|
||||
Curve::SurfaceIntersection { surfaces, .. } => {
|
||||
surfaces.iter().filter_map(|s| plane_frame(s)).any(|(qo, qn)| pts.iter().all(|p| on(*p, qo, qn)))
|
||||
}
|
||||
mt::Curve::IntersectionCurve(ic) => [ic.surface0(), ic.surface1()].into_iter().any(|s| matches!(s, mt::Surface::Plane(q) if holds(q))),
|
||||
};
|
||||
if !certified {
|
||||
return Ok(None);
|
||||
}
|
||||
}
|
||||
let wire: mt::Wire = chain.to_vec().into();
|
||||
Ok(Some(mt::Face::try_new(vec![wire], mt::Surface::Plane(plane)).map_err(|e| KernelError::Failed(format!("patch: {e}")))?))
|
||||
Ok(Some(plane))
|
||||
}
|
||||
|
||||
/// A plane record's origin and unit normal (`None` for any other surface).
|
||||
fn plane_frame(s: &Surface) -> Option<(Vec3, Vec3)> {
|
||||
match s {
|
||||
Surface::Plane { origin, u_axis, v_axis } => Some((Vec3::from(*origin), Vec3::from(*u_axis).cross(Vec3::from(*v_axis)).normalized()?)),
|
||||
Surface::PlaneThrough { origin, u_point, v_point } => {
|
||||
let o = Vec3::from(*origin);
|
||||
Some((o, (Vec3::from(*u_point) - o).cross(Vec3::from(*v_point) - o).normalized()?))
|
||||
}
|
||||
_ => None,
|
||||
}
|
||||
}
|
||||
|
||||
Reference in New Issue
Block a user