diff --git a/crates/engine/src/cmd/mod.rs b/crates/engine/src/cmd/mod.rs index 96a06f68..02bcbe59 100644 --- a/crates/engine/src/cmd/mod.rs +++ b/crates/engine/src/cmd/mod.rs @@ -20,6 +20,7 @@ mod parameters; mod sketch; mod sketch_constraints; mod sketch_create; +mod sketch_freeform; mod sketch_import; mod sketch_modify; mod sketch_project; diff --git a/crates/engine/src/cmd/sketch_freeform.rs b/crates/engine/src/cmd/sketch_freeform.rs new file mode 100644 index 00000000..b15c1825 --- /dev/null +++ b/crates/engine/src/cmd/sketch_freeform.rs @@ -0,0 +1,242 @@ +//! Trim and Break on free-form curves (ellipses, splines, conics). The curve is cut into exact +//! pieces (cubic Béziers become single-span control splines, conic arcs conics), so the result +//! stays exact; the pieces replace the original curve. + +use serde_json::{Value, json}; +use solvecraft_geom::{Seg2, Vec2}; +use solvecraft_sketch::{ConstraintKind, CurveKind, Sketch}; + +use super::sketch::add_c; +use crate::Result; +use crate::params::bad; + +/// Exact pieces of a free-form curve (polyline lines for the few kinds without them). +fn pieces(sk: &Sketch, ci: usize) -> Vec { + sk.exact_segs(ci).unwrap_or_else(|| sk.polyline(ci).windows(2).map(|w| Seg2::Line { a: w[0], b: w[1] }).collect()) +} + +/// Position along the pieces (piece index + local parameter) closest to `p`, and the distance. +fn locate(segs: &[Seg2], p: Vec2) -> (f64, f64) { + segs.iter() + .enumerate() + .map(|(i, s)| { + let t = s.closest_param(p); + (i as f64 + t, s.point_at(t).dist(p)) + }) + .min_by(|a, b| a.1.total_cmp(&b.1)) + .unwrap_or((0.0, f64::INFINITY)) +} + +fn point_at(segs: &[Seg2], u: f64) -> Vec2 { + let n = segs.len(); + let i = (u.floor() as usize).min(n.saturating_sub(1)); + segs.get(i).map(|s| s.point_at((u - i as f64).clamp(0.0, 1.0))).unwrap_or_default() +} + +/// Closest point on curve `ci` (exact pieces, or the analytic shape). +fn project(sk: &Sketch, ci: usize, p: Vec2) -> Vec2 { + if sk.curves.get(ci).is_some_and(|c| c.kind.is_freeform()) { + let segs = pieces(sk, ci); + let (u, _) = locate(&segs, p); + return point_at(&segs, u); + } + sk.shape(ci).map(|s| s.project(p)).unwrap_or(p) +} + +/// Where other curves meet free-form curve `ci`: (position along it, the point, the other curve). +fn cuts(sk: &Sketch, ci: usize, segs: &[Seg2]) -> Vec<(f64, Vec2, usize)> { + let mine = sk.polyline(ci); + let tol = 1e-6; + let mut out: Vec<(f64, Vec2, usize)> = Vec::new(); + let seg_x = |a: Vec2, b: Vec2, c: Vec2, d: Vec2| -> Option { + let (r, s) = (b - a, d - c); + let den = r.cross(s); + if den.abs() < 1e-15 { + return None; + } + let (t, u) = ((c - a).cross(s) / den, (c - a).cross(r) / den); + ((-1e-9..=1.0 + 1e-9).contains(&t) && (-1e-9..=1.0 + 1e-9).contains(&u)).then(|| a + r * t) + }; + for (oi, _) in sk.curves.iter().enumerate() { + if oi == ci { + continue; + } + let other = sk.polyline(oi); + let mut found: Vec = Vec::new(); + for w in mine.windows(2) { + for v in other.windows(2) { + if let Some(x) = seg_x(w[0], w[1], v[0], v[1]) { + found.push(x); + } + } + } + // End points of the other curve lying on this one. + if let Some((a, b)) = sk.curves.get(oi).and_then(|c| c.kind.ends()) { + for q in [a, b].iter().filter_map(|q| sk.point(*q)) { + if locate(segs, q).1 < 1e-6 { + found.push(q); + } + } + } + for x in found { + // Refine onto both exact curves by alternating projection. + let mut q = x; + for _ in 0..40 { + let n = project(sk, ci, project(sk, oi, q)); + if n.dist(q) < 1e-13 { + q = n; + break; + } + q = n; + } + let (u, d) = locate(segs, q); + if d < tol && !out.iter().any(|(_, y, _)| y.dist(q) < tol * 10.0) { + out.push((u, q, oi)); + } + } + } + out.sort_by(|a, b| a.0.total_cmp(&b.0)); + out +} + +/// Exact sub-pieces between positions `ua` < `ub` (positions may run past the end on closed +/// curves; they wrap). +fn span(segs: &[Seg2], ua: f64, ub: f64) -> Vec { + let n = segs.len() as f64; + let mut out = Vec::new(); + let mut u = ua; + while u < ub - 1e-12 && out.len() < 10_000 { + let i = u.floor(); + let next = (i + 1.0).min(ub); + let k = (i.rem_euclid(n)) as usize; + let Some(s) = segs.get(k).copied() else { break }; + let (t0, t1) = (u - i, next - i); + let mut piece = s; + if t1 < 1.0 - 1e-12 { + piece = piece.split_at(t1).0; + } + if t0 > 1e-12 { + // Split at t0 measured on the (possibly shortened) piece. + let tt = if t1 < 1.0 - 1e-12 { piece.closest_param(s.point_at(t0)) } else { t0 }; + piece = piece.split_at(tt).1; + } + out.push(piece); + u = next; + } + out +} + +/// Add pieces as curves, chained from point `start` (an existing point index or a new one) to +/// `end`. Returns the new curve indices. +fn add_pieces(sk: &mut Sketch, segs: &[Seg2], start: Option, end: Option, construction: bool) -> Result> { + let mut made = Vec::new(); + let n = segs.len(); + let mut prev = match start { + Some(i) => i, + None => match segs.first() { + Some(s) => sk.add_point(s.start(), None)?, + None => return Ok(made), + }, + }; + for (k, s) in segs.iter().enumerate() { + let to = match (k + 1 == n, end) { + (true, Some(e)) => e, + _ => sk.add_point(s.end(), None)?, + }; + let kind = match *s { + Seg2::Line { .. } => CurveKind::Line { a: prev, b: to }, + Seg2::Cubic { p1, p2, .. } => { + let (c1, c2) = (sk.add_point(p1, None)?, sk.add_point(p2, None)?); + CurveKind::Spline { pts: vec![prev, c1, c2, to], control: true, degree: 3 } + } + Seg2::Conic { apex, w, .. } => { + let x = sk.add_point(apex, None)?; + CurveKind::Conic { a: prev, b: to, apex: x, rho: w / (1.0 + w) } + } + Seg2::Arc { .. } => return Err(bad("trim", "unexpected arc piece")), + }; + let c = sk.add_curve(kind, None)?; + if let Some(cu) = sk.curves.get_mut(c) { + cu.construction = construction; + } + made.push(c); + prev = to; + } + Ok(made) +} + +/// Trim (`remove`) or break a free-form curve at the crossings around `at`. +pub(super) fn trim_or_break(sk: &mut Sketch, ci: usize, at: Vec2, remove: bool, cmd: &str) -> Result { + let segs = pieces(sk, ci); + if segs.is_empty() { + return Err(bad(cmd, "the curve has no shape")); + } + let n = segs.len() as f64; + let closed = segs.first().zip(segs.last()).is_some_and(|(a, b)| a.start().dist(b.end()) < 1e-9); + let cs: Vec<(f64, Vec2, usize)> = cuts(sk, ci, &segs).into_iter().filter(|c| closed || (c.0 > 1e-9 && c.0 < n - 1e-9)).collect(); + let (uq, _) = locate(&segs, at); + let ends = sk.curves.get(ci).and_then(|c| c.kind.ends()); + let construction = sk.curves.get(ci).is_some_and(|c| c.construction); + // Intervals between consecutive cuts (closed curves wrap round). + let mut bounds: Vec<(f64, Option)> = cs.iter().map(|c| (c.0, Some(c.2))).collect(); + if bounds.is_empty() { + if remove { + sk.remove_curves(&[ci]); + return Ok(json!({"deleted": true})); + } + return Err(bad(cmd, "nothing to break at there")); + } + let intervals: Vec<(f64, f64, Option, Option)> = if closed { + if bounds.len() < 2 && remove { + return Err(bad(cmd, "a closed curve needs two crossings to trim")); + } + let m = bounds.len(); + (0..m) + .filter_map(|i| { + let (a, ca) = *bounds.get(i)?; + let (b, cb) = *bounds.get((i + 1) % m)?; + let b = if b <= a + 1e-12 { b + n } else { b }; + Some((a, b, ca, cb)) + }) + .collect() + } else { + bounds.insert(0, (0.0, None)); + bounds.push((n, None)); + bounds.windows(2).map(|w| (w[0].0, w[1].0, w[0].1, w[1].1)).collect() + }; + let inside = |(a, b, _, _): &(f64, f64, Option, Option)| { + let u = if closed && uq < *a { uq + n } else { uq }; + u >= *a && u <= *b + }; + let gone = if remove { intervals.iter().position(inside) } else { None }; + // Points at the cuts are shared by neighbouring intervals. + let mut cut_pt: Vec<(f64, usize)> = Vec::new(); + for (u, q, _) in &cs { + cut_pt.push((*u, sk.add_point(*q, None)?)); + } + let point_for = + |u: f64, cut_pt: &[(f64, usize)]| cut_pt.iter().find(|(x, _)| (x - u.rem_euclid(n)).abs() < 1e-9 || (x - u).abs() < 1e-9).map(|p| p.1); + let mut made = Vec::new(); + for (k, iv) in intervals.iter().enumerate() { + if Some(k) == gone { + continue; + } + let (a, b, ..) = *iv; + let start = if !closed && a <= 1e-12 { ends.map(|e| e.0) } else { point_for(a, &cut_pt) }; + let end = if !closed && b >= n - 1e-12 { ends.map(|e| e.1) } else { point_for(b, &cut_pt) }; + let sub = span(&segs, a, b); + made.extend(add_pieces(sk, &sub, start, end, construction)?); + } + // Cut points stay on the curves that cut there. + for ((_, _, oi), (_, pi)) in cs.iter().zip(&cut_pt) { + let used = sk.curves.iter().any(|c| c.kind.uses(*pi)); + if used { + add_c(sk, ConstraintKind::PointOnCurve { p: *pi, c: *oi })?; + } + } + let ids: Vec = made.iter().filter_map(|c| sk.curves.get(*c).map(|c| c.id.clone())).collect(); + sk.remove_curves(&[ci]); + let unused: Vec = cut_pt.iter().map(|x| x.1).filter(|pi| !sk.curves.iter().any(|c| c.kind.uses(*pi))).collect(); + sk.remove_points(&unused); + Ok(json!({"curves": ids})) +} diff --git a/crates/engine/src/cmd/sketch_modify.rs b/crates/engine/src/cmd/sketch_modify.rs index 94bdc4d9..55be056e 100644 --- a/crates/engine/src/cmd/sketch_modify.rs +++ b/crates/engine/src/cmd/sketch_modify.rs @@ -258,6 +258,11 @@ fn trim(s: &mut Session, p: &Value) -> Result { let (out, info) = edit(s, p, cmd, false, |sk, doc| { let ci = curve_arg(sk, p, cmd)?; let id = sk.curves.get(ci).map(|c| c.id.clone()).unwrap_or_default(); + if sk.curves.get(ci).is_some_and(|c| c.kind.is_freeform()) { + let r = super::sketch_freeform::trim_or_break(sk, ci, at, true, cmd)?; + doc.prune_model_params(); + return Ok(json!({"trimmed": id, "pieces": r["curves"]})); + } let sh = sk.shape(ci).ok_or_else(|| bad(cmd, "curve"))?; let tq = sh.param(sh.project(at)); let (cuts, full) = bracket(sk, ci, tq); @@ -311,6 +316,11 @@ fn break_cmd(s: &mut Session, p: &Value) -> Result { let cmd = "BreakSketchCmd"; let (out, info) = edit(s, p, cmd, false, |sk, _| { let ci = curve_arg(sk, p, cmd)?; + if sk.curves.get(ci).is_some_and(|c| c.kind.is_freeform()) { + let at = req_vec2(cmd, p, "at")?; + let r = super::sketch_freeform::trim_or_break(sk, ci, at, false, cmd)?; + return Ok(r["curves"].as_array().map(|a| a.iter().filter_map(|x| x.as_str().map(str::to_string)).collect()).unwrap_or_default()); + } let sh = sk.shape(ci).ok_or_else(|| bad(cmd, "curve"))?; let (ts, cutters): (Vec, Vec>) = if let Some(list) = p.get("points").and_then(Value::as_array) { let mut v: Vec = list.iter().take(1000).filter_map(vec2).map(|q| sh.param(sh.project(q))).collect(); diff --git a/crates/engine/src/cmd/sketch_modify_tests.rs b/crates/engine/src/cmd/sketch_modify_tests.rs index eab71904..b52481ec 100644 --- a/crates/engine/src/cmd/sketch_modify_tests.rs +++ b/crates/engine/src/cmd/sketch_modify_tests.rs @@ -342,3 +342,43 @@ fn free_form_profiles_extrude_exactly() { let m = run(&mut s, "MeasureCommand", json!({})); assert!(m["bodies"][0]["faces"].as_u64().unwrap() <= 5, "{m}"); } + +#[test] +fn trim_and_break_free_form_curves_exactly() { + let pi = std::f64::consts::PI; + let mut s = new_sketch(); + // An ellipse cut by a line through its centre: trimming the lower half leaves the upper + // half, which closes with the line into half the ellipse. + let e = ids(&run(&mut s, "CircleElipse", json!({"center": [0, 0], "major": [10, 0], "minor_radius": 4}))["curves"])[0].clone(); + let l = ids(&run(&mut s, "DrawPolyline", json!({"points": [[-15, 0], [15, 0]]}))["curves"])[0].clone(); + let r = run(&mut s, "TrimSketchCmd", json!({"curve": e, "at": [0, -4]})); + assert!(!r["result"]["pieces"].as_array().unwrap().is_empty(), "{r}"); + let a = profiles(&s); + assert_eq!(a.len(), 1, "{a:?}"); + assert!((a[0] - pi * 20.0).abs() < 1e-6, "{a:?}"); + let _ = l; + // A fit spline broken where a line crosses it: two pieces that still trace the curve. + let mut s = new_sketch(); + let sp = ids(&run(&mut s, "DrawSpline", json!({"points": [[0, 0], [10, 10], [20, 0], [30, 5]]}))["curves"])[0].clone(); + let before = sketch(&s).polyline(sketch(&s).curve_index(&sp).unwrap()); + run(&mut s, "DrawPolyline", json!({"points": [[15, -10], [15, 20]]})); + let r = run(&mut s, "BreakSketchCmd", json!({"curve": sp, "at": [15, 5]})); + let pieces = ids(&r["curves"]); + assert_eq!(pieces.len(), 4, "three Béziers, the middle one split: {r}"); + let sk = sketch(&s); + // Every original sample lies on one of the pieces. + let polys: Vec> = pieces.iter().map(|id| sk.polyline(sk.curve_index(id).unwrap())).collect(); + for q in before.iter().step_by(5) { + let d = polys + .iter() + .flat_map(|p| { + p.windows(2).map(|w| { + let (a, b) = (w[0], w[1]); + let t = ((*q - a).dot(b - a) / (b - a).len2()).clamp(0.0, 1.0); + q.dist(a + (b - a) * t) + }) + }) + .fold(f64::INFINITY, f64::min); + assert!(d < 0.05, "{q:?} {d}"); + } +} diff --git a/crates/geom/src/profile.rs b/crates/geom/src/profile.rs index 7197c90c..7fc235e1 100644 --- a/crates/geom/src/profile.rs +++ b/crates/geom/src/profile.rs @@ -132,6 +132,63 @@ impl Seg2 { Seg2::Cubic { .. } | Seg2::Conic { .. } => self.integrate(|s, t| s.derivative(t).len()), } } + /// The segment split at parameter `t` in (0, 1) (same kind, exact; conic halves are + /// re-normalised so their end weights stay 1). + pub fn split_at(&self, t: f64) -> (Seg2, Seg2) { + let t = t.clamp(1e-9, 1.0 - 1e-9); + match *self { + Seg2::Line { a, b } => { + let m = a.lerp(b, t); + (Seg2::Line { a, b: m }, Seg2::Line { a: m, b }) + } + Seg2::Arc { center, radius, start, sweep } => ( + Seg2::Arc { center, radius, start, sweep: sweep * t }, + Seg2::Arc { center, radius, start: start + sweep * t, sweep: sweep * (1.0 - t) }, + ), + Seg2::Cubic { p0, p1, p2, p3 } => { + let (a, b, c) = (p0.lerp(p1, t), p1.lerp(p2, t), p2.lerp(p3, t)); + let (d, e) = (a.lerp(b, t), b.lerp(c, t)); + let m = d.lerp(e, t); + (Seg2::Cubic { p0, p1: a, p2: d, p3: m }, Seg2::Cubic { p0: m, p1: e, p2: c, p3 }) + } + Seg2::Conic { a, apex, b, w } => { + // Homogeneous de Casteljau: (a, 1), (w·apex, w), (b, 1). + let h = |p: Vec2, wt: f64| (p * wt, wt); + let lerp = |x: (Vec2, f64), y: (Vec2, f64)| (x.0 * (1.0 - t) + y.0 * t, x.1 * (1.0 - t) + y.1 * t); + let (p0, p1, p2) = (h(a, 1.0), h(apex, w), h(b, 1.0)); + let (q0, q1) = (lerp(p0, p1), lerp(p1, p2)); + let r = lerp(q0, q1); + let pt = |x: (Vec2, f64)| x.0 / x.1.max(1e-300); + let m = pt(r); + ( + Seg2::Conic { a, apex: pt(q0), b: m, w: q0.1 / r.1.max(1e-300).sqrt() }, + Seg2::Conic { a: m, apex: pt(q1), b, w: q1.1 / r.1.max(1e-300).sqrt() }, + ) + } + } + } + /// Parameter of the point of the segment closest to `p` (sampling, then Newton). + pub fn closest_param(&self, p: Vec2) -> f64 { + let n = 64; + let mut t = + (0..=n).map(|i| i as f64 / n as f64).min_by(|x, y| self.point_at(*x).dist(p).total_cmp(&self.point_at(*y).dist(p))).unwrap_or(0.0); + for _ in 0..30 { + let (q, d) = (self.point_at(t), self.derivative(t)); + let (d2, h) = ((self.derivative((t + 1e-6).min(1.0)) - self.derivative((t - 1e-6).max(0.0))) / 2e-6, q - p); + let f = h.dot(d); + let df = d.dot(d) + h.dot(d2); + if df.abs() < 1e-300 { + break; + } + let nt = (t - f / df).clamp(0.0, 1.0); + if (nt - t).abs() < 1e-14 { + t = nt; + break; + } + t = nt; + } + t + } /// The two halves of the segment (same kind, exact). pub fn split_half(&self) -> (Seg2, Seg2) { match *self { diff --git a/crates/geom/src/tests.rs b/crates/geom/src/tests.rs index a8be2b49..5fe4e25b 100644 --- a/crates/geom/src/tests.rs +++ b/crates/geom/src/tests.rs @@ -138,3 +138,24 @@ fn cubic_and_conic_segments_measure_exactly() { assert!((l.area_term() + r.area_term() - seg.area_term()).abs() < 1e-12); } } + +#[test] +fn split_at_and_closest_param() { + use crate::{Seg2, Vec2}; + let w = std::f64::consts::FRAC_1_SQRT_2; + let quarter = Seg2::Conic { a: Vec2::new(10.0, 0.0), apex: Vec2::new(10.0, 10.0), b: Vec2::new(0.0, 10.0), w }; + for t in [0.2, 0.5, 0.8] { + let (l, r) = quarter.split_at(t); + assert!(l.end().dist(quarter.point_at(t)) < 1e-12 && r.start().dist(l.end()) < 1e-12); + for s in [0.3, 0.7] { + assert!((l.point_at(s).len() - 10.0).abs() < 1e-12 && (r.point_at(s).len() - 10.0).abs() < 1e-12); + } + assert!((l.area_term() + r.area_term() - quarter.area_term()).abs() < 1e-12); + } + let c = Seg2::Cubic { p0: Vec2::new(0.0, 0.0), p1: Vec2::new(1.0, 2.0), p2: Vec2::new(3.0, 2.0), p3: Vec2::new(4.0, 0.0) }; + let q = c.point_at(0.37); + assert!((c.closest_param(q + Vec2::new(0.0, 0.0)) - 0.37).abs() < 1e-9); + let (l, r) = c.split_at(0.37); + assert!(l.end().dist(q) < 1e-12 && r.start().dist(q) < 1e-12); + assert!((quarter.closest_param(Vec2::new(20.0, 20.0)) - 0.5).abs() < 1e-9); +}