mirror of
https://github.com/bherbruck/solidcraft.git
synced 2026-11-03 10:35:58 +00:00
brep: a seam from pole to pole gets pcurves a period apart
A whole sphere bounded only by its seam (SMLib's sphere: one meridian used twice) came out with both uses on the same u line, a loop that encloses nothing. Two causes in the per-loop closure: - The pole-junction search took any positive signed area as the right turn. With the two uses on one line the area is rounding, and in most frames that rounding is a tiny positive number, so "no shift" won over the period shift. A turn now needs more than 1e-9 of the loop's period box. - Junctions counted as poles only where the surface collapses exactly. Ends within the tolerance of a pole (as a kernel's records have them) were closed to the nearest period instead. A junction is now free where moving its periodic parameter by half a period stays within the tolerance, and fit samples within a quarter of it take their free parameter from a neighbour, so the pcurve runs straight into the pole. Tests: the sphere bounded by its seam in four frames, three meridian angles, circle and NURBS meridians, both loop orders, ends at or a tenth of the tolerance off the poles; and a torus bounded by its two seams. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01QWhMwdqJYQKQq8RCLSL1LT
This commit is contained in:
@@ -101,7 +101,7 @@ pub fn pcurves(surface_record: &Surface, loops: &[Vec<EdgeUse<'_>>], tolerance:
|
||||
};
|
||||
fitted.push((pc, u.forward));
|
||||
}
|
||||
on.close(&mut fitted)?;
|
||||
on.close(&mut fitted, tolerance)?;
|
||||
out.push(fitted.into_iter().map(|(c, _)| c).collect());
|
||||
}
|
||||
Ok(out)
|
||||
@@ -254,10 +254,29 @@ impl On {
|
||||
Some(q)
|
||||
}
|
||||
|
||||
fn sample(&self, curve: &CurveEval, t: f64, near: Option<[f64; 2]>, accept: f64) -> Result<Sample> {
|
||||
fn sample(&self, curve: &CurveEval, t: f64, near: Option<[f64; 2]>, accept: f64, reach: f64) -> Result<Sample> {
|
||||
let p = curve.point(t).ok_or_else(|| Error::Failed(format!("the curve has no point at {t}")))?;
|
||||
let uv = self.project(p, near, accept).ok_or_else(|| self.off(p, t, near))?;
|
||||
Ok(Sample { t, p, uv, pole: self.surface.pole(uv), mid: None })
|
||||
Ok(Sample { t, p, uv, pole: self.pole(uv, reach), mid: None })
|
||||
}
|
||||
|
||||
/// Whether `uv` is at a pole along each direction, or within `reach` of one: moving that
|
||||
/// parameter by half its period (or across its range) moves the point by no more than `reach`
|
||||
/// at its speed there. There the parameter is as free as at the pole itself.
|
||||
fn pole(&self, uv: [f64; 2], reach: f64) -> [bool; 2] {
|
||||
let mut out = self.surface.pole(uv);
|
||||
let Some((_, su, sv)) = self.surface.derivs(uv[0], uv[1]) else { return out };
|
||||
for (d, speed) in [su.len(), sv.len()].into_iter().enumerate() {
|
||||
let span = match (self.dom.period.get(d).copied().flatten(), self.dom.range.get(d).copied().flatten()) {
|
||||
(Some(period), _) => 0.5 * period,
|
||||
(None, Some((a, b))) => (b - a).abs(),
|
||||
(None, None) => continue,
|
||||
};
|
||||
if let Some(flag) = out.get_mut(d) {
|
||||
*flag |= speed * span <= reach;
|
||||
}
|
||||
}
|
||||
out
|
||||
}
|
||||
|
||||
/// The error for a curve point `p` (at `t`) farther than the fit allows from the surface.
|
||||
@@ -272,6 +291,9 @@ impl On {
|
||||
/// sampled twice as densely until none does.
|
||||
fn fit(&self, record: &Curve, curve: &CurveEval, [t0, t1]: [f64; 2], tol: f64) -> Result<Curve2> {
|
||||
let (target, accept) = (0.5 * tol, 0.5 * tol);
|
||||
// Samples this near a pole take their free parameter from a neighbour, which moves
|
||||
// their lift by at most this (well inside the target).
|
||||
let reach = 0.25 * tol;
|
||||
// Interior knots: where the curve's own parameterisation is less smooth (a pcurve is no
|
||||
// smoother than its curve), else even.
|
||||
let mut inner = breaks(record, t0, t1);
|
||||
@@ -289,7 +311,7 @@ impl On {
|
||||
let (a, b) = (w.first().copied().unwrap_or(t0), w.get(1).copied().unwrap_or(t1));
|
||||
for k in usize::from(!samples.is_empty())..=PER_SPAN {
|
||||
let t = if k == PER_SPAN { b } else { a + (b - a) * k as f64 / PER_SPAN as f64 };
|
||||
let s = self.sample(curve, t, near, accept)?;
|
||||
let s = self.sample(curve, t, near, accept, reach)?;
|
||||
near = Some(s.uv);
|
||||
samples.push(s);
|
||||
}
|
||||
@@ -340,7 +362,7 @@ impl On {
|
||||
let t = 0.5 * (s.t + n.t);
|
||||
let uv = self.project(m, Some(s.uv), accept).ok_or_else(|| self.off(m, t, Some(s.uv)))?;
|
||||
dense.push(Sample { mid: None, ..s });
|
||||
dense.push(Sample { t, p: m, uv, pole: self.surface.pole(uv), mid: None });
|
||||
dense.push(Sample { t, p: m, uv, pole: self.pole(uv, reach), mid: None });
|
||||
}
|
||||
_ => dense.push(s),
|
||||
}
|
||||
@@ -352,10 +374,11 @@ impl On {
|
||||
|
||||
/// The loop's uses moved by whole periods so each starts where the one before it ends (in
|
||||
/// loop order and direction), then all together into the surface's base period where the
|
||||
/// loop's lowest point allows. Where uses meet at a pole the periodic parameter is free:
|
||||
/// there the offsets that turn the loop counter-clockwise in (u, v) (a face boundary seen
|
||||
/// from the natural normal) are taken, else the nearest.
|
||||
fn close(&self, uses: &mut [(Curve2, bool)]) -> Result<()> {
|
||||
/// loop's lowest point allows. Where uses meet at a pole (within `tol`) the periodic
|
||||
/// parameter is free: there the offsets that turn the loop counter-clockwise in (u, v) (a
|
||||
/// face boundary seen from the natural normal) by more than rounding are taken, else the
|
||||
/// nearest. A seam used twice from pole to pole thus encloses its face, never nothing.
|
||||
fn close(&self, uses: &mut [(Curve2, bool)], tol: f64) -> Result<()> {
|
||||
let n = uses.len();
|
||||
let mut ends = Vec::with_capacity(n);
|
||||
for (c, fwd) in uses.iter() {
|
||||
@@ -367,7 +390,7 @@ impl On {
|
||||
let Some(period) = self.dom.period.get(d).copied().flatten().filter(|p| *p > 0.0) else { continue };
|
||||
// Start after a pole junction if there is one: from there continuity fixes the
|
||||
// offsets up to the next pole.
|
||||
let at_pole: Vec<bool> = ends.iter().map(|(s, _)| self.surface.pole(*s).get(d).copied().unwrap_or(false)).collect();
|
||||
let at_pole: Vec<bool> = ends.iter().map(|(s, _)| self.pole(*s, tol).get(d).copied().unwrap_or(false)).collect();
|
||||
let first = at_pole.iter().position(|&p| p).unwrap_or(0);
|
||||
let order: Vec<usize> = (0..n).map(|k| (first + k) % n).collect();
|
||||
let mut shift = vec![0.0; n];
|
||||
@@ -384,6 +407,11 @@ impl On {
|
||||
// After each later pole one period less, none or more: the right turn, then the
|
||||
// fewest moves.
|
||||
if !free.is_empty() && free.len() <= 6 {
|
||||
// Less than this is rounding: a loop that encloses nothing (a seam's two uses
|
||||
// on one line) has no turn, whatever the sign of its area.
|
||||
let other = uses.iter().flat_map(|(c, _)| points(c).into_iter().filter_map(|q| q.get(1 - d).copied()));
|
||||
let (lo, hi) = other.fold((f64::INFINITY, f64::NEG_INFINITY), |(a, b), x| (a.min(x), b.max(x)));
|
||||
let flat = if hi >= lo { 1e-9 * period * (hi - lo) } else { 0.0 };
|
||||
let mut best: Option<(bool, usize, Vec<f64>)> = None;
|
||||
for combo in 0..3usize.pow(free.len() as u32) {
|
||||
let (mut code, mut moved, mut add) = (combo, 0, 0.0);
|
||||
@@ -399,7 +427,7 @@ impl On {
|
||||
*slot += add;
|
||||
}
|
||||
}
|
||||
let right = area(uses, &extra, d) > 0.0;
|
||||
let right = area(uses, &extra, d) > flat;
|
||||
if best.as_ref().is_none_or(|(r, m, _)| (right && !r) || (right == *r && moved < *m)) {
|
||||
best = Some((right, moved, extra));
|
||||
}
|
||||
|
||||
@@ -234,6 +234,93 @@ fn meridians_run_straight_into_the_poles() {
|
||||
}
|
||||
}
|
||||
|
||||
/// A whole sphere bounded only by its seam (one meridian from pole to pole, used twice), in
|
||||
/// frames where rounding leaves the loop with the uses on one line a tiny positive area: the
|
||||
/// two uses come out a period apart and the loop turns counter-clockwise around the whole
|
||||
/// (u, v) box, in either loop order, for a circle or a NURBS meridian, with its ends at the
|
||||
/// poles or within a tenth of the tolerance of them.
|
||||
#[test]
|
||||
fn a_seam_from_pole_to_pole_encloses_the_whole_sphere() {
|
||||
let w = std::f64::consts::FRAC_1_SQRT_2;
|
||||
let tol = 1e-7;
|
||||
let frames = [
|
||||
about_z([0.; 3], false),
|
||||
about_z([1., -2., 3.], true),
|
||||
frame([1., 2., 3.], [0., 0.6, 0.8], [0., -0.8, 0.6]),
|
||||
frame([-4., 0.5, 2.], [0.36, 0.48, 0.8], [-0.8, 0.6, 0.]),
|
||||
];
|
||||
for sf in frames {
|
||||
let sphere = Surface::Sphere { frame: sf.clone(), radius: 6. };
|
||||
let ax = crate::brep::eval::Axes::new(&sf);
|
||||
for turn in [0f64, 2.1, -0.7] {
|
||||
let x = ax.x * turn.cos() + ax.y * turn.sin();
|
||||
let p = |v: Vec3| [v.x, v.y, v.z];
|
||||
let (o, r, z) = (ax.o, x * 6., ax.z * 6.);
|
||||
let circle = Curve::Circle { frame: Frame { origin: p(o), x_axis: p(x), y_axis: p(ax.z) }, radius: 6. };
|
||||
for off in [0., 0.1 * tol] {
|
||||
// South pole to north pole through the equator, as two rational quarter arcs,
|
||||
// the ends `off` to the side of the poles.
|
||||
let h = |v: Vec3, wt: f64| [v.x * wt, v.y * wt, v.z * wt, wt];
|
||||
let side = ax.y * off;
|
||||
let nurbs = Curve::RationalSpline {
|
||||
knots: vec![0., 0., 0., 0.5, 0.5, 1., 1., 1.],
|
||||
homogeneous: vec![h(o - z + side, 1.), h(o - z + r, w), h(o + r, 1.), h(o + z + r, w), h(o + z - side, 1.)],
|
||||
};
|
||||
for (seam, iv) in [(&circle, [-FRAC_PI_2, FRAC_PI_2]), (&nurbs, [0., 1.])] {
|
||||
for order in [[true, false], [false, true]] {
|
||||
let l: Vec<EdgeUse> = order.iter().map(|&forward| EdgeUse { edge: 0, curve: seam, interval: iv, forward }).collect();
|
||||
let pcs = &pcurves(&sphere, std::slice::from_ref(&l), tol).unwrap()[0];
|
||||
let what = format!("{sf:?} {turn} {off} {seam:?} {order:?}");
|
||||
for pc in pcs {
|
||||
assert!(lift_error(&sphere, seam, iv, pc) <= tol, "{what}");
|
||||
}
|
||||
let mid = 0.5 * (iv[0] + iv[1]);
|
||||
let (a, b) = (pcs[0].eval(mid).unwrap(), pcs[1].eval(mid).unwrap());
|
||||
assert!(((a[0] - b[0]).abs() - TAU).abs() < 1e-6 && (a[1] - b[1]).abs() < 1e-12, "a period apart: {a:?} {b:?} {what}");
|
||||
let (ends, area) = walk(&l, pcs);
|
||||
for k in 0..2 {
|
||||
assert!((ends[k].1[1] - ends[(k + 1) % 2].0[1]).abs() < 1e-6, "closed in v: {ends:?} {what}");
|
||||
}
|
||||
// The box less the corners `walk` cuts (it joins each use's 16th
|
||||
// sample to the next use's start): about 15/16 of it.
|
||||
assert!((area / (TAU * PI * 15. / 16.) - 1.).abs() < 1e-2, "counter-clockwise around the box: {area} {what}");
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
/// A whole torus bounded only by its two seams (a circle about the axis and a circle around
|
||||
/// the tube, each used twice): no pole, continuity alone puts each seam's uses a period apart
|
||||
/// and the loop around the whole (u, v) box.
|
||||
#[test]
|
||||
fn a_torus_bounded_by_its_seams_spans_both_periods() {
|
||||
let frames = [about_z([0.; 3], false), frame([1., 2., 3.], [0., 0.6, 0.8], [0., -0.8, 0.6])];
|
||||
for sf in frames {
|
||||
let torus = Surface::Torus { frame: sf.clone(), major_radius: 10., minor_radius: 2. };
|
||||
let ax = crate::brep::eval::Axes::new(&sf);
|
||||
let p = |v: Vec3| [v.x, v.y, v.z];
|
||||
// The outer equator (v = 0, u from 0) and the tube circle at u = 0 (v from 0, upwards).
|
||||
let around = Curve::Circle { frame: sf.clone(), radius: 12. };
|
||||
let tube = Curve::Circle { frame: Frame { origin: p(ax.o + ax.x * 10.), x_axis: p(ax.x), y_axis: p(ax.z) }, radius: 2. };
|
||||
let full = [0., TAU];
|
||||
let l = vec![
|
||||
EdgeUse { edge: 0, curve: &around, interval: full, forward: true },
|
||||
EdgeUse { edge: 1, curve: &tube, interval: full, forward: true },
|
||||
EdgeUse { edge: 0, curve: &around, interval: full, forward: false },
|
||||
EdgeUse { edge: 1, curve: &tube, interval: full, forward: false },
|
||||
];
|
||||
let pcs = &pcurves(&torus, std::slice::from_ref(&l), 1e-9).unwrap()[0];
|
||||
let (ends, area) = walk(&l, pcs);
|
||||
for k in 0..4 {
|
||||
let (end, start) = (ends[k].1, ends[(k + 1) % 4].0);
|
||||
assert!((end[0] - start[0]).abs() < 1e-12 && (end[1] - start[1]).abs() < 1e-12, "{k}: {end:?} {start:?}");
|
||||
}
|
||||
assert!((area / (TAU * TAU) - 1.).abs() < 1e-12, "counter-clockwise around the box: {area}");
|
||||
}
|
||||
}
|
||||
|
||||
/// An intersection record evaluates to the point of both surfaces in the leader's normal
|
||||
/// plane, and its pcurves on either surface fit.
|
||||
#[test]
|
||||
|
||||
Reference in New Issue
Block a user