From b1ca05d990f94053bd31bafc37fdd393f1819b62 Mon Sep 17 00:00:00 2001 From: aiongg Date: Mon, 26 Oct 2026 20:04:06 +0000 Subject: [PATCH] refactor(three-d): move solid.rs items to solid/{build,inflate,membrane,rings}.rs --- crates/three-d/src/solid.rs | 1331 +------------------------- crates/three-d/src/solid/build.rs | 460 +++++++++ crates/three-d/src/solid/inflate.rs | 305 ++++++ crates/three-d/src/solid/membrane.rs | 400 ++++++++ crates/three-d/src/solid/rings.rs | 163 ++++ 5 files changed, 1343 insertions(+), 1316 deletions(-) create mode 100644 crates/three-d/src/solid/build.rs create mode 100644 crates/three-d/src/solid/inflate.rs create mode 100644 crates/three-d/src/solid/membrane.rs create mode 100644 crates/three-d/src/solid/rings.rs diff --git a/crates/three-d/src/solid.rs b/crates/three-d/src/solid.rs index eb3814b9a..043e34766 100644 --- a/crates/three-d/src/solid.rs +++ b/crates/three-d/src/solid.rs @@ -13,6 +13,21 @@ use vectorcraft_geom::{BezPath, FillRule, ParamCurve, PathData, Point, Rect, Vec use crate::{MAX_FACES, Revolve, V3, eye_distance, perspective_scale}; +mod build; +mod inflate; +mod membrane; +mod rings; + +#[cfg(test)] +use build::{Camera, cuts, edge_distance, newell, section}; +pub use build::{extrude, plane}; +pub use inflate::{inflate, inflate_together}; +#[cfg(test)] +use membrane::{Edges, Membrane, cell_triangles}; +use rings::signed_area; +#[cfg(test)] +use rings::{MAX_SOLID_POINTS, rings}; + /// The view shared by every 3D effect: rotation, perspective and directional lighting. #[derive(Clone, Copy, Debug, PartialEq)] pub struct View { @@ -289,1322 +304,6 @@ impl Inflate { } } -/// Outline points Extrude & Bevel and Inflate sample at most: every point costs walls, bevel -/// bands and distance tests. -const MAX_SOLID_POINTS: usize = 2048; -/// Outline points Rotate samples at most: it only turns the outline, so an outlined paragraph fits. -const MAX_PLANE_POINTS: usize = 65536; -/// Pieces a twisted wall's strip is cut into along one edge at most. -const MAX_PIECES: usize = 64; - -/// How many straight pieces follow `c` within `tolerance` (Wang's bound on a cubic's distance from -/// its chords), 1 to 32. -fn pieces(c: vectorcraft_geom::CubicBez, tolerance: f64) -> u32 { - let d = |a: Point, b: Point, c: Point| (a.to_vec2() - 2.0 * b.to_vec2() + c.to_vec2()).hypot(); - let bend = d(c.p0, c.p1, c.p2).max(d(c.p1, c.p2, c.p3)); - let n = (0.75 * bend / tolerance).sqrt().ceil(); - if n.is_finite() { (n as u32).clamp(1, 32) } else { 32 } -} - -/// The rings of a shape: closed polygons with the filled side on their left, (-dy, dx): outer -/// outlines run one way and holes the other, so nonzero winding fills them as the art is. Curves -/// are followed within a tolerance of 1/2000 of `size` (at least 0.05 pt), so small curves take -/// few points and large ones more; more than `limit` points in all is refused. -/// -/// Under the nonzero `rule` an outline with filled art on both sides is no edge of the shape and -/// is left out: a contour inside another that runs the same way doesn't cut a hole, as it would -/// if nesting alone decided (#1724). -fn rings(path: &PathData, rule: FillRule, size: f64, limit: usize) -> Result>, String> { - let tolerance = (size / 2000.0).max(0.05); - let mut count = 0; - let mut out = vec![]; - for sp in &path.subpaths { - let segments = if sp.closed { sp.anchors.len() } else { sp.anchors.len().saturating_sub(1) }; - if segments > limit { - return Err("The shape is too complex for a 3D effect".into()); - } - let mut ring: Vec = vec![]; - if let Some(first) = sp.anchors.first() { - ring.push(first.p); - } - for i in 0..segments { - let curve = sp.segment(i); - let samples = if sp.segment_is_line(i) { 1 } else { pieces(curve, tolerance) }; - count += samples as usize; - if count > limit { - return Err("The shape is too complex for a 3D effect".into()); - } - for k in 1..=samples { - ring.push(curve.eval(f64::from(k) / f64::from(samples))); - } - } - if ring.iter().any(|p| !p.x.is_finite() || !p.y.is_finite() || p.x.abs() > 1e8 || p.y.abs() > 1e8) { - return Err("The shape contains invalid coordinates".into()); - } - ring.dedup_by(|a, b| a.distance(*b) < 1e-9); - while ring.len() > 1 && ring.first().zip(ring.last()).is_some_and(|(a, b)| a.distance(*b) < 1e-9) { - ring.pop(); - } - if ring.len() >= 3 && signed_area(&ring).abs() > 1e-9 { - out.push(ring); - } - } - // Which rings hold which: each ring against the others, its edges walked only where its - // bounds hold the point, and within [`MAX_CONTAINMENT_EDGES`] edges in all. - let too_complex = || String::from("The shape is too complex for a 3D effect"); - let mut work = 0usize; - let mut within = |out: &[Vec], bounds: &[Rect], i: usize| -> Result, String> { - let p = out.get(i).and_then(|r| r.first()).copied().unwrap_or_default(); - let mut held_by = vec![]; - for (j, (o, b)) in out.iter().zip(bounds).enumerate() { - if j == i || !(b.x0 <= p.x && p.x <= b.x1 && b.y0 <= p.y && p.y <= b.y1) { - continue; - } - work = work.saturating_add(o.len()); - if work > MAX_CONTAINMENT_EDGES { - return Err(too_complex()); - } - if contains(o, p) { - held_by.push(j); - } - } - Ok(held_by) - }; - let bounds_of = |out: &[Vec]| -> Vec { - out.iter().map(|r| r.iter().map(|q| Rect::from_points(*q, *q)).reduce(|a, b| a.union(b)).unwrap_or_default()).collect() - }; - if rule == FillRule::NonZero { - let turn = |r: &[Point]| if signed_area(r) > 0.0 { 1 } else { -1 }; - let bounds = bounds_of(&out); - let mut edge = Vec::with_capacity(out.len()); - for (i, r) in out.iter().enumerate() { - // The winding number just outside the ring, and just inside it. - let outside: i32 = within(&out, &bounds, i)?.into_iter().filter_map(|j| out.get(j)).map(|o| turn(o)).sum(); - edge.push((outside != 0) != (outside + turn(r) != 0)); - } - let mut edge = edge.into_iter(); - out.retain(|_| edge.next().unwrap_or(true)); - } - // A ring inside an odd number of others is a hole. - let bounds = bounds_of(&out); - let holes: Vec = (0..out.len()).map(|i| within(&out, &bounds, i).map(|held| held.len() % 2 == 1)).collect::>()?; - for (r, hole) in out.iter_mut().zip(holes) { - if (signed_area(r) > 0.0) == hole { - r.reverse(); - } - } - Ok(out) -} - -/// Most ring edges [`rings`] walks to tell which rings hold which: every ring is tested against -/// the others, so many rings of many points each would otherwise take billions of steps. -const MAX_CONTAINMENT_EDGES: usize = 50_000_000; - -fn signed_area(r: &[Point]) -> f64 { - r.iter().zip(r.iter().cycle().skip(1)).take(r.len()).map(|(a, b)| a.x * b.y - b.x * a.y).sum::() / 2.0 -} - -fn contains(r: &[Point], p: Point) -> bool { - let mut inside = false; - for (a, b) in r.iter().zip(r.iter().cycle().skip(1)).take(r.len()) { - if (a.y > p.y) != (b.y > p.y) && p.x < a.x + (p.y - a.y) / (b.y - a.y) * (b.x - a.x) { - inside = !inside; - } - } - inside -} - -/// The unit normal on the left of `a`→`b`: towards the filled side. -fn left(a: Point, b: Point) -> (f64, f64) { - let (dx, dy) = (b.x - a.x, b.y - a.y); - let l = dx.hypot(dy); - if l < 1e-12 { (0.0, 0.0) } else { (-dy / l, dx / l) } -} - -/// `ring` moved `w` points towards its filled side (negative: away from it), each corner along -/// its bisector, its mitre limited to four times `w` so sharp corners don't shoot off. Each point -/// stays the partner of its source point, so rings at different insets join face to face. -fn inset(ring: &[Point], w: f64) -> Vec { - if w == 0.0 { - return ring.to_vec(); - } - let n = ring.len(); - (0..n) - .map(|i| { - let p = ring[i]; - let prev = ring[(i + n - 1) % n]; - let next = ring[(i + 1) % n]; - let (ax, ay) = left(prev, p); - let (bx, by) = left(p, next); - let (mut mx, mut my) = (ax + bx, ay + by); - let ml = mx.hypot(my); - if ml < 1e-9 { - (mx, my) = (ax, ay); - } else { - (mx, my) = (mx / ml, my / ml); - } - let cos = (mx * ax + my * ay).max(0.25); - Point::new(p.x + mx * w / cos, p.y + my * w / cos) - }) - .collect() -} - -/// Turns 3D points into document points and shades faces for one view. -struct Camera { - view: Revolve, - origin: Point, - /// The eye's distance in front of the centre; `None` without perspective. - eye: Option, - /// Deeper than any part of the object: faces drawn first go this far behind it. - far: f64, - light: V3, -} - -impl Camera { - fn new(view: View, bounds: Rect, depth: f64) -> Self { - let view = view.revolve(); - let [lx, ly, lz] = view.light_direction(); - let size = (bounds.width() / 2.0).max(bounds.height() / 2.0).max(depth / 2.0).max(1.0); - Self { view, origin: bounds.center(), eye: eye_distance(view.perspective, size), far: 16.0 * size, light: V3(lx, ly, lz) } - } - fn rotate(&self, v: V3) -> V3 { - self.view.rotate(v) - } - /// A point of the art, `z` points in front of its plane, in view space. - fn place(&self, p: Point, z: f64) -> V3 { - self.rotate(V3(p.x - self.origin.x, p.y - self.origin.y, z)) - } - fn project(&self, v: V3) -> Point { - let f = perspective_scale(self.eye, v.2); - Point::new(self.origin.x + v.0 * f, self.origin.y + v.1 * f) - } - /// The cosine between a face's outward unit normal `n` and the direction from `at` (both in - /// view space) to the eye: positive when the face is turned towards the viewer. - fn facing(&self, at: V3, n: V3) -> f64 { - let to_eye = self.eye.map_or(V3(0.0, 0.0, 1.0), |e| V3(0.0, 0.0, e).sub(at).unit()); - n.dot(to_eye) - } - fn brightness(&self, n: V3) -> f64 { - crate::brightness(&self.view, n, self.light) - } -} - -/// Facets being collected for one solid. -struct Builder<'a> { - camera: &'a Camera, - /// A closed solid hides the faces turned away from the viewer; an open one shows their - /// insides. - closed: bool, - facets: Vec, -} - -/// Where a face goes in the order the facets are drawn, back to front. -#[derive(Clone, Copy, PartialEq, Eq)] -enum Order { - /// The face is sorted by the depth of its center, except that a face seen almost edge-on goes - /// first, under everything. - Depth, - /// The face is sorted by the depth of its center even when it is seen almost edge-on, so that - /// a piece of a twisted wall at the solid's outline covers the walls behind it. - Twisted, - /// The face is drawn after every other face when it is turned toward the viewer. A flat end of - /// the solid has the whole solid behind its plane, so it covers whatever it overlaps. Its depth - /// is `f64::MAX`, so the faces drawn last by several solids sorted together keep the order they - /// were added in. - Last, -} - -impl Builder<'_> { - /// Add a face of view-space contours with outward normal `n` (view space), drawn in `order`. - fn face(&mut self, contours: Vec>, n: V3, order: Order) -> Result<(), String> { - let count: usize = contours.iter().map(Vec::len).sum(); - if count < 3 { - return Ok(()); - } - let centre = contours.iter().flatten().fold(V3::default(), |s, v| V3(s.0 + v.0, s.1 + v.1, s.2 + v.2)); - let centre = V3(centre.0 / count as f64, centre.1 / count as f64, centre.2 / count as f64); - let n = n.unit(); - let facing = self.camera.facing(centre, n); - // A closed solid hides the faces turned away. One seen almost exactly edge-on is kept and, - // unless it is a piece of a twisted wall, drawn first, under everything: on a curved - // surface rounding alone turns such faces either way, and dropping every other one would - // leave gaps. - let (n, behind) = if facing > 0.0 { - (n, false) - } else if self.closed && facing > -0.02 { - (n, true) - } else if self.closed { - return Ok(()); - } else { - (V3(-n.0, -n.1, -n.2), false) - }; - #[cfg(test)] - let view = contours.clone(); - let contours: Vec> = contours.iter().map(|c| c.iter().map(|v| self.camera.project(*v)).collect()).collect(); - if contours.iter().flatten().any(|p| !p.x.is_finite() || !p.y.is_finite()) { - return Err("The 3D projection is not finite".into()); - } - if contours.iter().map(|c| signed_area(c).abs()).sum::() < 1e-10 { - return Ok(()); - } - let depth = if behind && order != Order::Twisted { - centre.2 - self.camera.far - } else if order == Order::Last && facing > 0.0 { - f64::MAX - } else { - centre.2 - }; - self.facets.push(Facet { - contours, - depth, - brightness: self.camera.brightness(n), - illumination: crate::illumination(&self.camera.view, n, self.camera.light), - rgb_gradient: None, - #[cfg(test)] - view, - }); - if self.facets.len() > MAX_FACES { - return Err("The 3D surface is too complex; simplify the shape".into()); - } - Ok(()) - } - - /// A planar surface triangle with shared, smooth outward normals. Culling and sorting - /// still use its geometry, rather than the normals used only for lighting. - fn smooth_triangle(&mut self, vertices: [V3; 3], normals: [V3; 3], back: bool) -> Result<(), String> { - let contours = vec![vertices.to_vec()]; - let n = newell(&contours); - let n = if back { V3(-n.0, -n.1, -n.2) } else { n }; - let centre = vertices.iter().fold(V3::default(), |s, v| V3(s.0 + v.0 / 3.0, s.1 + v.1 / 3.0, s.2 + v.2 / 3.0)); - let flip = !self.closed && self.camera.facing(centre, n.unit()) <= 0.0; - let before = self.facets.len(); - self.face(contours, n, Order::Depth)?; - if self.facets.len() == before || !self.camera.view.shade { - return Ok(()); - } - let values: Vec<_> = normals - .into_iter() - .map(|n| crate::unclipped_illumination(&self.camera.view, if flip { V3(-n.0, -n.1, -n.2) } else { n }, self.camera.light)) - .collect(); - if let Some(face) = self.facets.last_mut() { - let points = face.contours.first().map(Vec::as_slice).unwrap_or_default(); - face.rgb_gradient = crate::RgbGradient::triangle(points, &values); - } - Ok(()) - } - - /// The band between two partner rings (same length), as one quad per edge, or for a twisted - /// wall with `cuts` (see [`cuts`]) as `cuts[i]` quads along edge `i`. `outward` is the - /// direction (before rotation) the band faces away from the solid, for orienting normals. - fn band(&mut self, a: &[V3], b: &[V3], cuts: &[usize], outward: impl Fn(usize) -> V3) -> Result<(), String> { - let n = a.len().min(b.len()); - // At 0 and 1, the corners themselves. - let lerp = |p: V3, q: V3, s: f64| { - if s <= 0.0 { - p - } else if s >= 1.0 { - q - } else { - V3(p.0 + (q.0 - p.0) * s, p.1 + (q.1 - p.1) * s, p.2 + (q.2 - p.2) * s) - } - }; - let order = if cuts.is_empty() { Order::Depth } else { Order::Twisted }; - for i in 0..n { - let j = (i + 1) % n; - let pieces = cuts.get(i).copied().unwrap_or(1).max(1); - for k in 0..pieces { - let (s0, s1) = (k as f64 / pieces as f64, (k + 1) as f64 / pieces as f64); - let quad = vec![lerp(a[i], a[j], s0), lerp(a[i], a[j], s1), lerp(b[i], b[j], s1), lerp(b[i], b[j], s0)]; - // A twist bends the strip between two slices into a saddle: its normal is the whole - // quad's average (Newell's), not one corner's triangle, which leans far to one side. - let mut normal = newell(std::slice::from_ref(&quad)); - if normal.dot(normal) < 1e-18 { - continue; - } - if normal.dot(self.camera.rotate(outward(i))) < 0.0 { - normal = V3(-normal.0, -normal.1, -normal.2); - } - self.face(vec![quad], normal, order)?; - } - } - Ok(()) - } - - fn finish(mut self) -> Vec { - self.facets.sort_by(|a, b| a.depth.total_cmp(&b.depth)); - self.facets - } -} - -/// How many pieces each edge of a twisted wall is drawn in, from two neighboring slices of the -/// wall, `a` and `b` (partner rings): as many as make each piece about as long along the edge as -/// the band between the slices is wide there, 1 to [`MAX_PIECES`]. Twisted walls pass in front of -/// one another, and a strip the length of its edge has no one depth where another crosses it. -fn cuts(a: &[V3], b: &[V3]) -> Vec { - let n = a.len().min(b.len()); - let length = |p: V3, q: V3| q.sub(p).dot(q.sub(p)).sqrt(); - (0..n) - .map(|i| { - let j = (i + 1) % n; - let (along, across) = (length(a[i], a[j]).max(length(b[i], b[j])), length(a[i], b[i]).max(length(a[j], b[j]))); - let count = (along / across.max(1e-9)).ceil(); - if count.is_finite() { (count as usize).clamp(1, MAX_PIECES) } else { MAX_PIECES } - }) - .collect() -} - -/// Lowers the [`cuts`] of every ring's edges (one list per ring) so that `bands` bands of each -/// make no more than `room` quads in all, or one per edge and band when that is more: the pieces -/// past each edge's first share the room in proportion. -fn fit(cuts: &mut [Vec], bands: usize, room: usize) { - let edges = cuts.iter().map(Vec::len).sum::().saturating_mul(bands); - let total = cuts.iter().flatten().sum::().saturating_mul(bands); - let room = room.max(edges); - if total > room { - let share = (room - edges) as f64 / (total - edges) as f64; - for c in cuts.iter_mut().flatten() { - *c = 1 + ((c.saturating_sub(1)) as f64 * share).floor() as usize; - } - } -} - -/// The area-weighted normal of (possibly non-planar) contours, by Newell's method: for contours -/// winding with their filled side on the left as seen from in front, it points to the front. -fn newell(contours: &[Vec]) -> V3 { - let mut n = V3::default(); - for c in contours { - for (a, b) in c.iter().zip(c.iter().cycle().skip(1)).take(c.len()) { - let x = a.cross(*b); - n = V3(n.0 + x.0, n.1 + x.1, n.2 + x.2); - } - } - n -} - -fn check_bounds(b: Rect) -> Result<(), String> { - if [b.x0, b.y0, b.x1, b.y1].iter().all(|v| v.is_finite() && v.abs() <= 1e8) { - Ok(()) - } else { - Err("3D effects require finite coordinates within 100 million points".into()) - } -} - -/// The art's outward direction at edge `i` of `ring` (before rotation), tilted towards the face -/// it borders by `z` (+1 front, -1 back, 0 a side wall). -fn outward(ring: &[Point], i: usize, z: f64) -> V3 { - let n = ring.len(); - let (lx, ly) = left(ring[i], ring[(i + 1) % n]); - V3(-lx, -ly, z) -} - -/// The cross-section of a side: each ring of the shape at depth fraction `t` (0 front, 1 back), -/// turned by the twist and scaled by the taper about the centre of `bounds`. -fn section(ring: &[Point], bounds: Rect, o: &Extrude, t: f64) -> Vec { - let turn = Turn::new(bounds, o.twist, o.taper, t); - ring.iter().map(|p| turn.apply(*p)).collect() -} - -/// The twist and taper at depth fraction `t` (0 front, 1 back): a turn by `twist · t` and a -/// scale from 1 to `taper`, about the centre of the bounds. -#[derive(Clone, Copy)] -struct Turn { - centre: Point, - sin: f64, - cos: f64, - scale: f64, -} - -impl Turn { - fn new(bounds: Rect, twist: f64, taper: f64, t: f64) -> Self { - let (sin, cos) = (twist * t).to_radians().sin_cos(); - Self { centre: bounds.center(), sin, cos, scale: 1.0 + (taper - 1.0) * t } - } - fn apply(&self, p: Point) -> Point { - let (x, y) = ((p.x - self.centre.x) * self.scale, (p.y - self.centre.y) * self.scale); - Point::new(self.centre.x + x * self.cos - y * self.sin, self.centre.y + x * self.sin + y * self.cos) - } -} - -/// One step of an extrusion's outline: how far it is inset (negative: outset) and its depth -/// fraction (0 front, 1 back), with the face it borders (+1 front bevel, -1 back bevel, 0 wall). -#[derive(Clone, Copy)] -struct Step { - inset: f64, - t: f64, - face: f64, -} - -/// The outline's steps from the front cap's edge to the back cap's edge. -fn steps(o: &Extrude) -> Vec { - let d = o.depth.max(1e-6); - let mut out = vec![]; - let bevel = o.bevel.filter(|b| b.width > 0.0 || b.height > 0.0); - // A bevel can't be deeper than the extrusion (half of it on each side with both sides). - let sides = if bevel.is_some_and(|b| b.both_sides) { 2.0 } else { 1.0 }; - let height = bevel.map_or(0.0, |b| b.height.min(o.depth / sides)); - // The bevel's cross-section across the width: (share across, share down). - let profile = |b: &Bevel| -> Vec<(f64, f64)> { - let mut pts = vec![(0.0, 0.0)]; - let n = b.repeat.max(1) as f64; - let run = (1.0 - b.space) / n; - let gap = if b.repeat > 1 { b.space / (n - 1.0) } else { 0.0 }; - for k in 0..b.repeat.max(1) { - let start = k as f64 * (run + gap); - for s in b.shape.samples().into_iter().skip(1) { - let (x, y) = b.shape.at(s); - pts.push((start + x * run, (k as f64 + y) / n)); - } - if b.repeat > 1 && k + 1 < b.repeat { - pts.push((start + run + gap, (k as f64 + 1.0) / n)); - } - } - if let Some(last) = pts.last_mut() { - last.0 = 1.0; - } - pts - }; - // With the bevel inside, the cap is inset by its width and the wall is the outline; outside, - // the cap is the outline and the wall is outset by the width. - let wall = match bevel { - Some(b) if !b.inside => -b.width, - _ => 0.0, - }; - // A twist turns a bevel's steps too: each step between two samples is cut into steps that - // turn no further than 6°, as the walls' slices do at most. - let refine = |pts: Vec<(f64, f64)>| -> Vec<(f64, f64)> { - if o.twist == 0.0 { - return pts; - } - let mut out = Vec::with_capacity(pts.len()); - for (k, p) in pts.iter().enumerate() { - if let Some(q) = k.checked_sub(1).and_then(|k| pts.get(k)) { - let turn = (p.1 - q.1).abs() * height / d * o.twist.abs(); - let n = ((turn / 6.0).ceil() as usize).clamp(1, 64); - for m in 1..n { - let s = m as f64 / n as f64; - out.push((q.0 + (p.0 - q.0) * s, q.1 + (p.1 - q.1) * s)); - } - } - out.push(*p); - } - out - }; - if let Some(b) = bevel { - for (x, y) in refine(profile(&b)) { - let inset = if b.inside { b.width * (1.0 - x) } else { -b.width * x }; - out.push(Step { inset, t: y * height / d, face: 1.0 }); - } - } else { - out.push(Step { inset: 0.0, t: 0.0, face: 1.0 }); - } - // Twist and taper bend the walls: cut them into slices. - let slices = if o.twist != 0.0 || o.taper != 1.0 { ((o.twist.abs() / 6.0).ceil() as usize).clamp(8, 64) } else { 1 }; - let back_height = if bevel.is_some_and(|b| b.both_sides) { height } else { 0.0 }; - let (t0, t1) = (height / d, 1.0 - back_height / d); - for k in 0..=slices { - let t = t0 + (t1 - t0) * k as f64 / slices as f64; - out.push(Step { inset: wall, t, face: 0.0 }); - } - if let Some(b) = bevel.filter(|b| b.both_sides) { - for (x, y) in refine(profile(&b)).into_iter().rev() { - let inset = if b.inside { b.width * (1.0 - x) } else { -b.width * x }; - out.push(Step { inset, t: 1.0 - y * height / d, face: -1.0 }); - } - } - out.dedup_by(|a, b| (a.inset - b.inset).abs() < 1e-9 && (a.t - b.t).abs() < 1e-12); - out -} - -/// Extrude & Bevel: the shape's outline extruded backwards, with optional twist, taper, bevels -/// and caps. -pub fn extrude(path: &PathData, rule: FillRule, bounds: Rect, options: Extrude, view: View) -> Result, String> { - check_bounds(bounds)?; - let o = options.sanitized(); - // No depth: the flat shape, seen from either side. - if o.depth <= 0.0 { - return plane(path, rule, bounds, view); - } - let rings = rings(path, rule, bounds.width().max(bounds.height()), MAX_SOLID_POINTS)?; - let camera = Camera::new(view, bounds, o.depth); - let mut b = Builder { camera: &camera, closed: o.cap, facets: vec![] }; - let steps = steps(&o); - let z = |t: f64| o.depth / 2.0 - o.depth * t; - let ring_at = - |ring: &[Point], s: &Step| -> Vec { section(&inset(ring, s.inset), bounds, &o, s.t).iter().map(|p| camera.place(*p, z(s.t))).collect() }; - // Every band of a twisted outline, bevels included, is cut as its walls' first slice is (with - // bevels and no walls, as its shallowest band is), within the facet limit. - let depth_of = |s: &[Step]| s.get(1).zip(s.first()).map_or(0.0, |(b, a)| (b.t - a.t).abs()); - let wall = steps - .windows(2) - .find(|s| s[0].face == 0.0 && s[1].face == 0.0) - .or_else(|| steps.windows(2).filter(|s| depth_of(s) > 1e-12).min_by(|a, c| depth_of(a).total_cmp(&depth_of(c)))); - let mut pieces: Vec> = rings - .iter() - .map(|ring| match wall { - Some([s0, s1]) if o.twist != 0.0 => cuts(&ring_at(ring, s0), &ring_at(ring, s1)), - _ => vec![], - }) - .collect(); - fit(&mut pieces, steps.len().saturating_sub(1), MAX_FACES); - for (ring, pieces) in rings.iter().zip(&pieces) { - let placed: Vec> = steps.iter().map(|s| ring_at(ring, s)).collect(); - for (k, pair) in placed.windows(2).enumerate() { - let ([a, c], Some(s0), Some(s1)) = (pair, steps.get(k), steps.get(k + 1)) else { continue }; - // A band between two front bevel steps borders the front, one reaching a back bevel - // step the back; the rest are walls. - let face = if s1.face < 0.0 { - -1.0 - } else if s0.face > 0.0 && s1.face > 0.0 { - 1.0 - } else { - 0.0 - }; - b.band(a, c, pieces, |i| outward(ring, i, face))?; - } - } - if o.cap { - // A cap turned toward the viewer covers the walls it overlaps, even those nearer than its - // center: twisted walls, and the walls of a fill inside the solid of its stroke. - for (s, n) in [(steps.first(), 1.0), (steps.last(), -1.0)] { - let Some(s) = s else { continue }; - b.face(rings.iter().map(|r| ring_at(r, s)).collect(), camera.rotate(V3(0.0, 0.0, n)), Order::Last)?; - } - } - Ok(b.finish()) -} - -/// Rotate: the flat art turned in 3D about the centre of `bounds` (which may be larger than the -/// shape's own, so that every part of one object turns about the same point). -pub fn plane(path: &PathData, rule: FillRule, bounds: Rect, view: View) -> Result, String> { - check_bounds(bounds)?; - let rings = rings(path, rule, bounds.width().max(bounds.height()), MAX_PLANE_POINTS)?; - let camera = Camera::new(view, bounds, 0.0); - let mut b = Builder { camera: &camera, closed: false, facets: vec![] }; - let contours: Vec> = rings.iter().map(|r| r.iter().map(|p| camera.place(*p, 0.0)).collect()).collect(); - // A plane is one flat face, seen from either side. Its normal is turned toward the viewer, so - // it is drawn last, as a cap is, and planes drawn in one order keep the order they are given in. - let up = camera.rotate(V3(0.0, 0.0, 1.0)); - let n = if camera.facing(V3::default(), up) < 0.0 { V3(-up.0, -up.1, -up.2) } else { up }; - b.face(contours, n, Order::Last)?; - Ok(b.finish()) -} - -/// The distance from `p` to the segment `a`–`b`. -fn segment_distance(p: Point, a: Point, b: Point) -> f64 { - let (dx, dy) = (b.x - a.x, b.y - a.y); - let l2 = dx * dx + dy * dy; - let t = if l2 < 1e-18 { 0.0 } else { (((p.x - a.x) * dx + (p.y - a.y) * dy) / l2).clamp(0.0, 1.0) }; - p.distance(Point::new(a.x + dx * t, a.y + dy * t)) -} - -/// The distance from `p` to the nearest edge of `rings`, edge by edge. -#[cfg(test)] -fn edge_distance(rings: &[Vec], p: Point) -> f64 { - rings - .iter() - .flat_map(|r| r.iter().zip(r.iter().cycle().skip(1)).take(r.len())) - .map(|(a, b)| segment_distance(p, *a, *b)) - .fold(f64::INFINITY, f64::min) -} - -/// The edges of `rings` filed in a grid of square buckets, so the nearest edge to a point is -/// found among the buckets around it rather than among all the edges. -struct Edges { - edges: Vec<(Point, Point)>, - normals: Vec<(Vec2, Vec2)>, - origin: Point, - size: f64, - n: usize, - buckets: Vec>, -} - -impl Edges { - fn new(rings: &[Vec]) -> Self { - const N: usize = 64; - let edges: Vec<(Point, Point)> = - rings.iter().flat_map(|r| r.iter().zip(r.iter().cycle().skip(1)).take(r.len()).map(|(a, b)| (*a, *b))).collect(); - let normals = rings - .iter() - .flat_map(|r| { - let normal = |i: usize| { - let (x, y) = left(r[i], r[(i + 1) % r.len()]); - Vec2::new(x, y) - }; - (0..r.len()).map(move |i| { - let n = normal(i); - let at = |other: Vec2| if n.dot(other) > 0.7 { (n + other).normalize() } else { n }; - (at(normal((i + r.len() - 1) % r.len())), at(normal((i + 1) % r.len()))) - }) - }) - .collect(); - let bounds = edges.iter().fold(None, |b: Option, (a, c)| Some(b.map_or(Rect::from_points(*a, *c), |b| b.union_pt(*a).union_pt(*c)))); - let bounds = bounds.unwrap_or(Rect::ZERO); - let size = (bounds.width().max(bounds.height()) / N as f64).max(1e-9); - let mut buckets = vec![vec![]; N * N]; - let index = |v: f64, o: f64| (((v - o) / size).floor().max(0.0) as usize).min(N - 1); - for (k, (a, b)) in edges.iter().enumerate() { - let (i0, i1) = (index(a.x.min(b.x), bounds.x0), index(a.x.max(b.x), bounds.x0)); - let (j0, j1) = (index(a.y.min(b.y), bounds.y0), index(a.y.max(b.y), bounds.y0)); - for j in j0..=j1 { - for i in i0..=i1 { - buckets[j * N + i].push(k as u32); - } - } - } - Self { edges, normals, origin: bounds.origin(), size, n: N, buckets } - } - - /// The bucket row or column of coordinate `v`, from the grid's origin `o` on that axis. - fn cell(&self, v: f64, o: f64) -> isize { - (((v - o) / self.size).floor() as isize).clamp(0, self.n as isize - 1) - } - - /// Is `p` within `tolerance` of an edge? Only the buckets `p` could share an edge with are - /// looked at. - fn near(&self, p: Point, tolerance: f64) -> bool { - self.touch(Rect::new(p.x - tolerance, p.y - tolerance, p.x + tolerance, p.y + tolerance)) && { - let n = self.n as isize; - let (i0, i1) = (self.cell(p.x - tolerance, self.origin.x), self.cell(p.x + tolerance, self.origin.x)); - let (j0, j1) = (self.cell(p.y - tolerance, self.origin.y), self.cell(p.y + tolerance, self.origin.y)); - (j0..=j1).any(|j| { - (i0..=i1).any(|i| { - self.buckets[(j * n + i) as usize] - .iter() - .filter_map(|k| self.edges.get(*k as usize)) - .any(|(a, b)| segment_distance(p, *a, *b) <= tolerance) - }) - }) - } - } - - /// Call `f` with each edge whose bounding box overlaps `r` (an edge in several of the - /// buckets `r` covers comes more than once). - fn each(&self, r: Rect, mut f: impl FnMut(Point, Point)) { - let n = self.n as isize; - for j in self.cell(r.y0, self.origin.y)..=self.cell(r.y1, self.origin.y) { - for i in self.cell(r.x0, self.origin.x)..=self.cell(r.x1, self.origin.x) { - for (a, b) in self.buckets[(j * n + i) as usize].iter().filter_map(|k| self.edges.get(*k as usize)) { - if a.x.min(b.x) <= r.x1 && a.x.max(b.x) >= r.x0 && a.y.min(b.y) <= r.y1 && a.y.max(b.y) >= r.y0 { - f(*a, *b); - } - } - } - } - } - - /// Does any edge pass through `r` (its bounding box, at least, overlapping it)? - fn touch(&self, r: Rect) -> bool { - let n = self.n as isize; - let (i0, i1) = (self.cell(r.x0, self.origin.x), self.cell(r.x1, self.origin.x)); - let (j0, j1) = (self.cell(r.y0, self.origin.y), self.cell(r.y1, self.origin.y)); - (j0..=j1).any(|j| { - (i0..=i1).any(|i| { - self.buckets[(j * n + i) as usize] - .iter() - .filter_map(|k| self.edges.get(*k as usize)) - .any(|(a, b)| a.x.min(b.x) <= r.x1 && a.x.max(b.x) >= r.x0 && a.y.min(b.y) <= r.y1 && a.y.max(b.y) >= r.y0) - }) - }) - } - - /// The distance from `p` to the nearest edge: the buckets are searched in growing squares - /// around `p`'s, until no unsearched bucket could hold a nearer edge. - fn distance(&self, p: Point) -> f64 { - self.closest(p).2 - } - - /// Closest outline point, its inward normal and distance. Curved outline samples share - /// corner normals, while a sharp source corner keeps each edge's own direction. - fn closest(&self, p: Point) -> (Point, Vec2, f64) { - let (fx, fy) = ((p.x - self.origin.x) / self.size, (p.y - self.origin.y) / self.size); - if !fx.is_finite() || !fy.is_finite() { - return (p, Vec2::ZERO, f64::INFINITY); - } - let n = self.n as isize; - let (ci, cj) = ((fx.floor() as isize).clamp(0, n - 1), (fy.floor() as isize).clamp(0, n - 1)); - // How far `p` is outside the grid: the nearest point of the grid to it lies in `p`'s bucket. - let outside = (fx.clamp(0.0, n as f64) - fx).hypot(fy.clamp(0.0, n as f64) - fy) * self.size; - let mut best = f64::INFINITY; - let mut nearest = (p, Vec2::ZERO); - for r in 0..n { - for j in (cj - r)..=(cj + r) { - for i in (ci - r)..=(ci + r) { - if (j - cj).abs() != r && (i - ci).abs() != r || i < 0 || j < 0 || i >= n || j >= n { - continue; - } - for k in &self.buckets[(j * n + i) as usize] { - if let Some((a, b)) = self.edges.get(*k as usize) { - let d = *b - *a; - let t = if d.hypot2() < 1e-18 { 0.0 } else { ((p - *a).dot(d) / d.hypot2()).clamp(0.0, 1.0) }; - let q = Point::new(a.x + d.x * t, a.y + d.y * t); - let distance = p.distance(q); - if distance < best { - best = distance; - let n = self.normals.get(*k as usize).map_or(Vec2::ZERO, |(a, b)| *a * (1.0 - t) + *b * t); - nearest = (q, if n.hypot2() > 0.0 { n.normalize() } else { n }); - } - } - } - } - } - // Buckets r + 1 away are at least r buckets from any point of `p`'s own, and that - // far across the grid from a `p` outside it (at right angles to its way in, at worst). - if best <= outside.hypot(r as f64 * self.size) { - break; - } - } - (nearest.0, nearest.1, best) - } -} - -/// Is `p` filled by `rings` (nonzero winding, as they are oriented)? -fn filled(rings: &[Vec], p: Point) -> bool { - rings.iter().filter(|r| contains(r, p)).count() % 2 == 1 -} - -/// `ring` cut to the rectangle `cell` (Sutherland–Hodgman: the rectangle is convex, the ring -/// needn't be). -fn clip(ring: &[Point], cell: Rect) -> Vec { - let mut pts = ring.to_vec(); - let sides: [fn(Point, Rect) -> f64; 4] = [|p, c| p.x - c.x0, |p, c| c.x1 - p.x, |p, c| p.y - c.y0, |p, c| c.y1 - p.y]; - for side in sides { - if pts.is_empty() { - break; - } - let mut out = vec![]; - for (a, b) in pts.iter().zip(pts.iter().cycle().skip(1)).take(pts.len()) { - let (da, db) = (side(*a, cell), side(*b, cell)); - if da >= 0.0 { - out.push(*a); - } - if (da >= 0.0) != (db >= 0.0) { - out.push(a.lerp(*b, da / (da - db))); - } - } - pts = out; - } - pts.dedup_by(|a, b| a.distance(*b) < 1e-9); - if pts.len() >= 3 && signed_area(&pts).abs() > 1e-12 { pts } else { vec![] } -} - -/// The shape's surface as a membrane under pressure, held at zero along its outline (and its -/// holes): the solution of ∇²φ = −1 on a grid, with each node's distance to the outline along the -/// grid lines taken into account where it is less than a grid step (Shortley–Weller), so the -/// outline needn't fall on the grid. √φ rises vertically from the outline and rounds smoothly -/// inside, without the creases a distance-to-the-edge height has along a shape's middle. -struct Membrane { - origin: Point, - step: f64, - nx: usize, - ny: usize, - /// φ at each node; outside the shape, extrapolated (negative) from the inside, so - /// interpolation crosses zero near the outline. - phi: Vec, - peak: f64, -} - -impl Membrane { - /// `edges` is the index of `rings`' edges (`Edges::new(rings)`). - fn solve(rings: &[Vec], edges: &Edges, origin: Point, step: f64, nx: usize, ny: usize) -> Self { - let node = |i: usize, j: usize| Point::new(origin.x + i as f64 * step, origin.y + j as f64 * step); - let idx = |i: usize, j: usize| j * nx + i; - // A node on the outline itself is held at zero with the outside. - let inside: Vec = (0..nx * ny) - .map(|k| { - let p = node(k % nx, k / nx); - filled(rings, p) && !edges.near(p, 1e-6 * step) - }) - .collect(); - // The nearest crossing of the outline from `p` along (dx, dy), within one step. - let crossing = |p: Point, dx: f64, dy: f64| -> f64 { - let mut best = step; - let reach = Rect::from_points(p, Point::new(p.x + dx * step, p.y + dy * step)); - edges.each(reach, |a, b| { - let t = if dy == 0.0 { - if (a.y > p.y) == (b.y > p.y) { - return; - } - (a.x + (p.y - a.y) / (b.y - a.y) * (b.x - a.x) - p.x) * dx - } else { - if (a.x > p.x) == (b.x > p.x) { - return; - } - (a.y + (p.x - a.x) / (b.x - a.x) * (b.y - a.y) - p.y) * dy - }; - if t > 1e-9 * step && t < best { - best = t; - } - }); - best - }; - // Each inside node's four neighbours: (weight, neighbour inside or None at the outline). - let mut stencil: Vec<[(f64, Option); 4]> = vec![[(0.0, None); 4]; nx * ny]; - for j in 0..ny { - for i in 0..nx { - let k = idx(i, j); - if !inside[k] { - continue; - } - let p = node(i, j); - let dirs = [(1.0, 0.0), (-1.0, 0.0), (0.0, 1.0), (0.0, -1.0)]; - let h: Vec = dirs.iter().map(|(dx, dy)| crossing(p, *dx, *dy)).collect(); - for (d, (dx, dy)) in dirs.iter().enumerate() { - let opp = h[d ^ 1]; - let w = 2.0 / (h[d] * (h[d] + opp)); - let (ni, nj) = (i as isize + *dx as isize, j as isize + *dy as isize); - let neighbour = (h[d] >= step && ni >= 0 && nj >= 0 && (ni as usize) < nx && (nj as usize) < ny) - .then(|| idx(ni as usize, nj as usize)) - .filter(|n| inside[*n]); - stencil[k][d] = (w, neighbour); - } - } - } - // Successive over-relaxation. - let mut phi = vec![0.0; nx * ny]; - let tolerance = 1e-9 * step * step * (nx * ny) as f64; - for _ in 0..2000 { - let mut change: f64 = 0.0; - for k in 0..nx * ny { - if !inside[k] { - continue; - } - let (mut sum, mut weight) = (1.0, 0.0); - for (w, n) in stencil[k] { - weight += w; - sum += w * n.map_or(0.0, |n| phi[n]); - } - let next = phi[k] + 1.8 * (sum / weight - phi[k]); - change = change.max((next - phi[k]).abs()); - phi[k] = next; - } - if change < tolerance { - break; - } - } - let peak = phi.iter().copied().fold(0.0, f64::max); - // Outside nodes next to the shape: φ continued past the outline at the slope of their - // inside neighbours, so values interpolated between them reach zero at the outline. - let mut ext = phi.clone(); - for j in 0..ny { - for i in 0..nx { - let k = idx(i, j); - if inside[k] { - continue; - } - let d_out = edges.distance(node(i, j)); - let mut slopes = vec![]; - for (di, dj) in [(1isize, 0isize), (-1, 0), (0, 1), (0, -1)] { - let (ni, nj) = (i as isize + di, j as isize + dj); - if ni < 0 || nj < 0 || ni as usize >= nx || nj as usize >= ny || !inside[idx(ni as usize, nj as usize)] { - continue; - } - let d_in = edges.distance(node(ni as usize, nj as usize)).max(1e-9); - slopes.push(phi[idx(ni as usize, nj as usize)] / d_in); - } - if !slopes.is_empty() { - ext[k] = -d_out * slopes.iter().sum::() / slopes.len() as f64; - } - } - } - Self { origin, step, nx, ny, phi: ext, peak } - } - - /// φ at `p`, interpolated between the four nodes around it. - fn at(&self, p: Point) -> f64 { - self.interpolate(p, |i, j| self.phi.get(j * self.nx + i).copied().unwrap_or(0.0)) - } - - /// Interpolate central node derivatives, rather than differentiating each bilinear cell. - /// The latter jumps at every grid line and shows the solver's grid in the lighting. - fn gradient(&self, p: Point) -> (f64, f64) { - let value = |i, j| self.phi.get(j * self.nx + i).copied().unwrap_or(0.0); - let dx = self.interpolate(p, |i, j| { - let (a, b) = (i.saturating_sub(1), (i + 1).min(self.nx - 1)); - (value(b, j) - value(a, j)) / ((b - a).max(1) as f64 * self.step) - }); - let dy = self.interpolate(p, |i, j| { - let (a, b) = (j.saturating_sub(1), (j + 1).min(self.ny - 1)); - (value(i, b) - value(i, a)) / ((b - a).max(1) as f64 * self.step) - }); - (dx, dy) - } - - fn interpolate(&self, p: Point, v: impl Fn(usize, usize) -> f64) -> f64 { - let fx = ((p.x - self.origin.x) / self.step).clamp(0.0, (self.nx - 1) as f64); - let fy = ((p.y - self.origin.y) / self.step).clamp(0.0, (self.ny - 1) as f64); - let (i, j) = ((fx.floor() as usize).min(self.nx.saturating_sub(2)), (fy.floor() as usize).min(self.ny.saturating_sub(2))); - let (tx, ty) = (fx - i as f64, fy - j as f64); - if self.nx < 2 || self.ny < 2 { - return v(0, 0); - } - let top = v(i, j) * (1.0 - tx) + v(i + 1, j) * tx; - let bottom = v(i, j + 1) * (1.0 - tx) + v(i + 1, j + 1) * tx; - top * (1.0 - ty) + bottom * ty - } -} - -/// A bounded ear clipping pass on a simple cell contour. Keep its collinear border nodes: -/// those shared samples can have different heights once lifted onto the membrane. A tiny -/// cell containing a complete hole, or a complex/degenerate contour, keeps its original fill. -fn cell_triangles(pieces: &[Vec]) -> Option> { - let mut triangles = vec![]; - for ring in pieces { - if ring.len() > 128 || signed_area(ring) <= 0.0 { - return None; - } - let mut ring = ring.clone(); - ring.dedup_by(|a, b| a.distance(*b) < 1e-9); - if ring.first().zip(ring.last()).is_some_and(|(a, b)| a.distance(*b) < 1e-9) { - ring.pop(); - } - while ring.len() > 3 { - let ear = (0..ring.len()).find(|&i| { - let n = ring.len(); - let (a, b, c) = (ring[(i + n - 1) % n], ring[i], ring[(i + 1) % n]); - if (b - a).cross(c - a) <= 1e-12 { - return false; - } - !ring.iter().enumerate().any(|(j, p)| { - j != i - && j != (i + n - 1) % n - && j != (i + 1) % n - && (b - a).cross(*p - a) >= -1e-12 - && (c - b).cross(*p - b) >= -1e-12 - && (a - c).cross(*p - c) >= -1e-12 - }) - })?; - let n = ring.len(); - triangles.push([ring[(ear + n - 1) % n], ring[ear], ring[(ear + 1) % n]]); - ring.remove(ear); - } - if let [a, b, c] = ring.as_slice() { - triangles.push([*a, *b, *c]); - } - } - Some(triangles) -} - -/// Most outline segments the booleans of one inflated object's mixed cells take together. -const MAX_MATERIAL_SEGMENTS: usize = 4_000_000; - -/// The part of `outline` (its rings, each with its bounds) inside `cell`, as rings (holes running -/// the other way). Only the rings whose bounds meet the cell's go into the boolean. -fn material(outline: &[(Rect, vectorcraft_geom::SubPath)], cell: Rect) -> Option>> { - let r = PathData::single(vectorcraft_geom::SubPath::polyline( - &[Point::new(cell.x0, cell.y0), Point::new(cell.x1, cell.y0), Point::new(cell.x1, cell.y1), Point::new(cell.x0, cell.y1)], - true, - )); - let near = PathData::new( - outline.iter().filter(|(b, _)| b.x0 <= cell.x1 && cell.x0 <= b.x1 && b.y0 <= cell.y1 && cell.y0 <= b.y1).map(|(_, s)| s.clone()).collect(), - ); - let rule = FillRule::NonZero; - let part = - vectorcraft_pathops::try_boolean(&r, rule, &near, rule, vectorcraft_pathops::BoolOp::Intersect, vectorcraft_pathops::DEFAULT_PRECISION) - .ok()?; - rings(&part, rule, cell.width(), MAX_PLANE_POINTS).ok() -} - -/// A shape to inflate: its rings and the membrane under them on the grid of the bounds. -struct Puff { - rings: Vec>, - /// The index of the rings' edges. - edges: Edges, - cell: f64, - nx: usize, - ny: usize, - m: Membrane, - /// How high the front rises at its peak. - height: f64, -} - -/// The rings of `path` and its membrane on the grid of `bounds`; `None` for a shape with no area. -fn puff(path: &PathData, rule: FillRule, bounds: Rect, o: &Inflate) -> Result, String> { - let rings = rings(path, rule, bounds.width().max(bounds.height()), MAX_SOLID_POINTS)?; - if rings.is_empty() { - return Ok(None); - } - let size = bounds.width().max(bounds.height()); - if size <= 0.0 { - return Ok(None); - } - const CELLS: usize = 48; - let cell = size / CELLS as f64; - let (nx, ny) = (((bounds.width() / cell).ceil() as usize).max(1), ((bounds.height() / cell).ceil() as usize).max(1)); - let edges = Edges::new(&rings); - let m = Membrane::solve(&rings, &edges, bounds.origin(), cell, nx + 1, ny + 1); - // A circle of radius r puffs into a half sphere of radius r at 50 % (its membrane peaks at r²/4). - let height = 2.0 * m.peak.sqrt() * o.volume / 50.0; - Ok(Some(Puff { rings, edges, cell, nx, ny, m, height })) -} - -/// Inflate: the shape puffed out like a cushion, its edge kept where the art's outline is. -pub fn inflate(path: &PathData, rule: FillRule, bounds: Rect, options: Inflate, view: View) -> Result, String> { - inflate_together(&[(path, rule)], bounds, options, view).pop().unwrap_or(Ok(vec![])) -} - -/// Inflate for several shapes of one object (a fill and the outlines of its strokes): each puffs -/// on its own over the grid of `bounds`, and one eye sees them all, as far away as the highest -/// puff needs with Perspective. One result for each shape, in order. -pub fn inflate_together(shapes: &[(&PathData, FillRule)], bounds: Rect, options: Inflate, view: View) -> Vec, String>> { - let o = options.sanitized(); - let puffs: Vec, String>> = - shapes.iter().map(|(path, rule)| check_bounds(bounds).and_then(|()| puff(path, *rule, bounds, &o))).collect(); - let reach = puffs.iter().flatten().flatten().map(|p| p.height).fold(0.0, f64::max); - puffs - .into_iter() - .map(|p| match p { - Ok(Some(p)) => inflated(p, bounds, o, view, reach), - Ok(None) => Ok(vec![]), - Err(e) => Err(e), - }) - .collect() -} - -/// The facets of one inflated shape, for a camera that allows for a puff `reach` high (or the -/// shape's own, if higher): `reach` is the highest puff of the shapes drawn together. -fn inflated(p: Puff, bounds: Rect, o: Inflate, view: View, reach: f64) -> Result, String> { - let Puff { rings, edges, cell, nx, ny, m, height } = p; - let outline: Vec<(Rect, vectorcraft_geom::SubPath)> = rings - .iter() - .filter_map(|r| { - let b = r.iter().skip(1).fold(Rect::from_points(*r.first()?, *r.first()?), |b, p| b.union_pt(*p)); - Some((b, vectorcraft_geom::SubPath::polyline(r, true))) - }) - .collect(); - // When no shape puffs (no volume, or shapes too narrow for the grid), every front (and a back - // that puffs too) lies flat in one plane: its cells are drawn last, as caps are. If one shape - // puffs, a flat front can be behind its puff, so the cells keep their depths. - let flat = reach <= 0.0; - let lift = |p: Point| { - if m.peak <= 0.0 || edges.near(p, 1e-9) { - return 0.0; - } - height * (m.at(p).max(0.0) / m.peak).sqrt() - }; - let normal = |p: Point| { - let (dx, dy) = m.gradient(p); - let (_, inward, distance) = edges.closest(p); - let t = (distance / (2.0 * cell)).clamp(0.0, 1.0); - let weight = (1.0 - t).powi(2) * (1.0 + 2.0 * t); - let gradient = Vec2::new(dx, dy); - let slope = gradient.dot(inward).max(0.0); - // The zero boundary condition also fixes the tangential derivative at the outline. - // Extrapolated outside solver nodes cannot supply it reliably. Blend to that known - // direction and its limiting phi near the rim, without showing the grid in highlights. - let gradient = gradient * (1.0 - weight) + inward * (slope * weight); - let phi = m.at(p).max(0.0) * (1.0 - weight) + slope * distance * weight; - // Multiplying the height-function normal by 2√(peak φ) avoids division by - // zero at the rim, where the limiting normal is horizontal. - let n = V3(-height * gradient.x, -height * gradient.y, 2.0 * (m.peak * phi).sqrt()).unit(); - if n.dot(n) > 0.5 { n } else { V3(0.0, 0.0, 1.0) } - }; - let smooth = view.shade && view.has_colored_lights(); - let camera = Camera::new(view, bounds, o.depth + 2.0 * height.max(reach)); - // Without its back cap the object is a shell open behind, its inside showing. - let mut b = Builder { camera: &camera, closed: o.cap, facets: vec![] }; - let half = o.depth / 2.0; - // The back: turned and scaled by the twist and taper, its puff too. - let back_turn = Turn::new(bounds, o.twist, o.taper, 1.0); - // The grid's cells, those the surface rises steeply across (near the outline) cut into - // halves or quarters each way, so the steep rim is drawn in facets as fine as its curve. - // Every corner sits on one fine lattice (a quarter cell), computed the same way wherever it - // is used, and a cell's border takes its finer neighbour's lattice points too: where cells - // of different fineness meet, their edges follow the surface through the same points, with - // no cracks. - const FINE: usize = 4; - let unit = cell / FINE as f64; - let at = |qx: usize, qy: usize| Point::new(bounds.x0 + qx as f64 * unit, bounds.y0 + qy as f64 * unit); - let base = |i: usize, j: usize| Rect::from_points(at(i * FINE, j * FINE), at((i + 1) * FINE, (j + 1) * FINE)); - // Segments the booleans of mixed cells may take in all ([`MAX_MATERIAL_SEGMENTS`]): each one - // builds the topology of every outline ring near its cell, so a self-crossing outline cost its - // crossings again in every cell it spans. Past it, cells take their clipped pieces. - let material_left = std::cell::Cell::new(MAX_MATERIAL_SEGMENTS); - // The shape's part of a rectangle, as rings. - let cut = |r: Rect| -> Vec> { - // A cell no edge crosses is wholly inside the shape or wholly outside it. - if !edges.touch(r) { - let corners = vec![Point::new(r.x0, r.y0), Point::new(r.x1, r.y0), Point::new(r.x1, r.y1), Point::new(r.x0, r.y1)]; - return if filled(&rings, r.center()) { vec![corners] } else { vec![] }; - } - let pieces: Vec> = rings.iter().map(|ring| clip(ring, r)).filter(|c| !c.is_empty()).collect(); - if pieces.iter().map(|p| signed_area(p)).sum::() <= 1e-9 * r.area() { - return vec![]; - } - // Where a hole's piece cancels part of an outline's, that only holds flat: lifted onto - // the surface the two no longer meet. Such a cell takes the shape's own part of it. - let mixed = pieces.iter().any(|p| signed_area(p) < 0.0) && pieces.iter().any(|p| signed_area(p) > 0.0); - if !mixed { - return pieces; - } - let near: usize = - outline.iter().filter(|(b, _)| b.x0 <= r.x1 && r.x0 <= b.x1 && b.y0 <= r.y1 && r.y0 <= b.y1).map(|(_, sp)| sp.anchors.len()).sum(); - match material_left.get().checked_sub(near) { - Some(left) => { - material_left.set(left); - material(&outline, r).unwrap_or(pieces) - } - None => pieces, - } - }; - // How many times each cell is halved each way (0, 1 or 2), from how steeply it rises. - let levels: Vec = (0..nx * ny) - .map(|k| { - let pieces = cut(base(k % nx, k / nx)); - let (lo, hi) = pieces.iter().flatten().map(|q| lift(*q)).fold((f64::INFINITY, f64::NEG_INFINITY), |(lo, hi), h| (lo.min(h), hi.max(h))); - let rise = if pieces.is_empty() { 0.0 } else { (hi - lo) / cell }; - if rise > 6.0 { - 2 - } else if rise > 2.0 { - 1 - } else { - 0 - } - }) - .collect(); - let level = - |i: isize, j: isize| -> u32 { if i < 0 || j < 0 || i as usize >= nx || j as usize >= ny { 0 } else { levels[j as usize * nx + i as usize] } }; - for j in 0..ny { - for i in 0..nx { - let l = level(i as isize, j as isize); - let step = FINE >> l; - // The lattice spacing along each side of this cell: its own or its neighbour's. - let (ii, jj) = (i as isize, j as isize); - let side_step = |dl: u32| FINE >> l.max(dl); - let (left, right, top, bottom) = - (side_step(level(ii - 1, jj)), side_step(level(ii + 1, jj)), side_step(level(ii, jj - 1)), side_step(level(ii, jj + 1))); - let (bx0, by0, bx1, by1) = (i * FINE, j * FINE, (i + 1) * FINE, (j + 1) * FINE); - for sy in (by0..by1).step_by(step) { - for sx in (bx0..bx1).step_by(step) { - let r = Rect::from_points(at(sx, sy), at(sx + step, sy + step)); - let pieces = cut(r); - if pieces.is_empty() { - continue; - } - // Lattice points along the cell's own border, between each edge's ends. - let stitch = |a: Point, b: Point| -> Vec { - let tol = 1e-7 * unit; - let along = |fixed_x: Option, fixed_y: Option, every: usize| -> Vec { - let mut out = vec![]; - match (fixed_x, fixed_y) { - (Some(qx), None) => { - let (lo, hi) = (a.y.min(b.y), a.y.max(b.y)); - for qy in (by0..=by1).step_by(every) { - let p = at(qx, qy); - if p.y > lo + tol && p.y < hi - tol { - out.push(p); - } - } - } - (None, Some(qy)) => { - let (lo, hi) = (a.x.min(b.x), a.x.max(b.x)); - for qx in (bx0..=bx1).step_by(every) { - let p = at(qx, qy); - if p.x > lo + tol && p.x < hi - tol { - out.push(p); - } - } - } - _ => {} - } - if (a.x > b.x) || (a.y > b.y) { - out.reverse(); - } - out - }; - let on = |v: f64, w: f64| (v - w).abs() <= tol; - let (x0, x1, y0, y1) = (at(bx0, by0).x, at(bx1, by1).x, at(bx0, by0).y, at(bx1, by1).y); - if on(a.x, x0) && on(b.x, x0) { - along(Some(bx0), None, left) - } else if on(a.x, x1) && on(b.x, x1) { - along(Some(bx1), None, right) - } else if on(a.y, y0) && on(b.y, y0) { - along(None, Some(by0), top) - } else if on(a.y, y1) && on(b.y, y1) { - along(None, Some(by1), bottom) - } else { - vec![] - } - }; - let pieces: Vec> = pieces - .iter() - .map(|p| { - let mut out = Vec::with_capacity(p.len()); - for (a, b) in p.iter().zip(p.iter().cycle().skip(1)).take(p.len()) { - out.push(*a); - out.extend(stitch(*a, *b)); - } - out - }) - .collect(); - if height > 0.0 - && smooth - && let Some(triangles) = cell_triangles(&pieces) - { - for points in triangles { - let front = points.map(|q| camera.place(q, half + lift(q))); - b.smooth_triangle(front, points.map(|q| camera.rotate(normal(q))), false)?; - if o.both_sides && o.cap { - let back = points.map(|q| camera.place(back_turn.apply(q), -half - lift(q) * back_turn.scale)); - let normals = points.map(|q| { - let n = normal(q); - // A uniform taper scales all three coordinates equally; - // twist turns the back normal about the shape's z-axis. - camera.rotate(V3(back_turn.cos * n.0 - back_turn.sin * n.1, back_turn.sin * n.0 + back_turn.cos * n.1, -n.2)) - }); - b.smooth_triangle(back, normals, true)?; - } - } - continue; - } - let lifts: Vec> = pieces.iter().map(|p| p.iter().map(|q| lift(*q)).collect()).collect(); - // Each piece lifted onto the surface. Its normal is the lifted piece's own average - // (Newell's), not the slope at one point: near the edge the surface turns - // vertical, and a point sample there misjudges which way a whole piece faces. - let front: Vec> = - pieces.iter().zip(&lifts).map(|(p, h)| p.iter().zip(h).map(|(q, h)| camera.place(*q, half + h)).collect()).collect(); - let n = newell(&front); - let order = if flat { Order::Last } else { Order::Depth }; - b.face(front, n, order)?; - if o.both_sides && o.cap { - // The back winds the same way round as the front, so its outward normal is - // reversed. - let back: Vec> = pieces - .iter() - .zip(&lifts) - .map(|(p, h)| p.iter().zip(h).map(|(q, h)| camera.place(back_turn.apply(*q), -half - h * back_turn.scale)).collect()) - .collect(); - let n = newell(&back); - b.face(back, V3(-n.0, -n.1, -n.2), order)?; - } - } - } - } - } - if !o.both_sides && o.cap { - let back: Vec> = rings.iter().map(|r| r.iter().map(|p| camera.place(back_turn.apply(*p), -half)).collect()).collect(); - // The back turned toward the viewer covers the walls it overlaps, even those nearer than - // its center: twisted walls, and the walls of a fill inside the solid of its stroke. - b.face(back, camera.rotate(V3(0.0, 0.0, -1.0)), Order::Last)?; - } - if o.depth > 0.0 { - // The sides, cut into slices that follow a twist or taper. - let slices = if o.twist != 0.0 || o.taper != 1.0 { ((o.twist.abs() / 6.0).ceil() as usize).clamp(8, 64) } else { 1 }; - let at = |ring: &[Point], k: usize| -> Vec { - let t = k as f64 / slices as f64; - let turn = Turn::new(bounds, o.twist, o.taper, t); - ring.iter().map(|p| camera.place(turn.apply(*p), half - o.depth * t)).collect() - }; - let mut pieces: Vec> = rings.iter().map(|r| if o.twist != 0.0 { cuts(&at(r, 0), &at(r, 1)) } else { vec![] }).collect(); - fit(&mut pieces, slices, MAX_FACES.saturating_sub(b.facets.len())); - for (ring, pieces) in rings.iter().zip(&pieces) { - for k in 0..slices { - b.band(&at(ring, k), &at(ring, k + 1), pieces, |i| outward(ring, i, 0.0))?; - } - } - } - Ok(b.finish()) -} - #[cfg(test)] mod tests; diff --git a/crates/three-d/src/solid/build.rs b/crates/three-d/src/solid/build.rs new file mode 100644 index 000000000..c9527edaf --- /dev/null +++ b/crates/three-d/src/solid/build.rs @@ -0,0 +1,460 @@ +use super::rings::{MAX_PIECES, MAX_PLANE_POINTS, MAX_SOLID_POINTS, inset, left, rings, signed_area}; +use super::*; + +/// Turns 3D points into document points and shades faces for one view. +pub(super) struct Camera { + pub(super) view: Revolve, + origin: Point, + /// The eye's distance in front of the centre; `None` without perspective. + eye: Option, + /// Deeper than any part of the object: faces drawn first go this far behind it. + far: f64, + pub(super) light: V3, +} + +impl Camera { + pub(super) fn new(view: View, bounds: Rect, depth: f64) -> Self { + let view = view.revolve(); + let [lx, ly, lz] = view.light_direction(); + let size = (bounds.width() / 2.0).max(bounds.height() / 2.0).max(depth / 2.0).max(1.0); + Self { view, origin: bounds.center(), eye: eye_distance(view.perspective, size), far: 16.0 * size, light: V3(lx, ly, lz) } + } + pub(super) fn rotate(&self, v: V3) -> V3 { + self.view.rotate(v) + } + /// A point of the art, `z` points in front of its plane, in view space. + pub(super) fn place(&self, p: Point, z: f64) -> V3 { + self.rotate(V3(p.x - self.origin.x, p.y - self.origin.y, z)) + } + pub(super) fn project(&self, v: V3) -> Point { + let f = perspective_scale(self.eye, v.2); + Point::new(self.origin.x + v.0 * f, self.origin.y + v.1 * f) + } + /// The cosine between a face's outward unit normal `n` and the direction from `at` (both in + /// view space) to the eye: positive when the face is turned towards the viewer. + pub(super) fn facing(&self, at: V3, n: V3) -> f64 { + let to_eye = self.eye.map_or(V3(0.0, 0.0, 1.0), |e| V3(0.0, 0.0, e).sub(at).unit()); + n.dot(to_eye) + } + pub(super) fn brightness(&self, n: V3) -> f64 { + crate::brightness(&self.view, n, self.light) + } +} + +/// Facets being collected for one solid. +pub(super) struct Builder<'a> { + pub(super) camera: &'a Camera, + /// A closed solid hides the faces turned away from the viewer; an open one shows their + /// insides. + pub(super) closed: bool, + pub(super) facets: Vec, +} + +/// Where a face goes in the order the facets are drawn, back to front. +#[derive(Clone, Copy, PartialEq, Eq)] +pub(super) enum Order { + /// The face is sorted by the depth of its center, except that a face seen almost edge-on goes + /// first, under everything. + Depth, + /// The face is sorted by the depth of its center even when it is seen almost edge-on, so that + /// a piece of a twisted wall at the solid's outline covers the walls behind it. + Twisted, + /// The face is drawn after every other face when it is turned toward the viewer. A flat end of + /// the solid has the whole solid behind its plane, so it covers whatever it overlaps. Its depth + /// is `f64::MAX`, so the faces drawn last by several solids sorted together keep the order they + /// were added in. + Last, +} + +impl Builder<'_> { + /// Add a face of view-space contours with outward normal `n` (view space), drawn in `order`. + pub(super) fn face(&mut self, contours: Vec>, n: V3, order: Order) -> Result<(), String> { + let count: usize = contours.iter().map(Vec::len).sum(); + if count < 3 { + return Ok(()); + } + let centre = contours.iter().flatten().fold(V3::default(), |s, v| V3(s.0 + v.0, s.1 + v.1, s.2 + v.2)); + let centre = V3(centre.0 / count as f64, centre.1 / count as f64, centre.2 / count as f64); + let n = n.unit(); + let facing = self.camera.facing(centre, n); + // A closed solid hides the faces turned away. One seen almost exactly edge-on is kept and, + // unless it is a piece of a twisted wall, drawn first, under everything: on a curved + // surface rounding alone turns such faces either way, and dropping every other one would + // leave gaps. + let (n, behind) = if facing > 0.0 { + (n, false) + } else if self.closed && facing > -0.02 { + (n, true) + } else if self.closed { + return Ok(()); + } else { + (V3(-n.0, -n.1, -n.2), false) + }; + #[cfg(test)] + let view = contours.clone(); + let contours: Vec> = contours.iter().map(|c| c.iter().map(|v| self.camera.project(*v)).collect()).collect(); + if contours.iter().flatten().any(|p| !p.x.is_finite() || !p.y.is_finite()) { + return Err("The 3D projection is not finite".into()); + } + if contours.iter().map(|c| signed_area(c).abs()).sum::() < 1e-10 { + return Ok(()); + } + let depth = if behind && order != Order::Twisted { + centre.2 - self.camera.far + } else if order == Order::Last && facing > 0.0 { + f64::MAX + } else { + centre.2 + }; + self.facets.push(Facet { + contours, + depth, + brightness: self.camera.brightness(n), + illumination: crate::illumination(&self.camera.view, n, self.camera.light), + rgb_gradient: None, + #[cfg(test)] + view, + }); + if self.facets.len() > MAX_FACES { + return Err("The 3D surface is too complex; simplify the shape".into()); + } + Ok(()) + } + + /// A planar surface triangle with shared, smooth outward normals. Culling and sorting + /// still use its geometry, rather than the normals used only for lighting. + pub(super) fn smooth_triangle(&mut self, vertices: [V3; 3], normals: [V3; 3], back: bool) -> Result<(), String> { + let contours = vec![vertices.to_vec()]; + let n = newell(&contours); + let n = if back { V3(-n.0, -n.1, -n.2) } else { n }; + let centre = vertices.iter().fold(V3::default(), |s, v| V3(s.0 + v.0 / 3.0, s.1 + v.1 / 3.0, s.2 + v.2 / 3.0)); + let flip = !self.closed && self.camera.facing(centre, n.unit()) <= 0.0; + let before = self.facets.len(); + self.face(contours, n, Order::Depth)?; + if self.facets.len() == before || !self.camera.view.shade { + return Ok(()); + } + let values: Vec<_> = normals + .into_iter() + .map(|n| crate::unclipped_illumination(&self.camera.view, if flip { V3(-n.0, -n.1, -n.2) } else { n }, self.camera.light)) + .collect(); + if let Some(face) = self.facets.last_mut() { + let points = face.contours.first().map(Vec::as_slice).unwrap_or_default(); + face.rgb_gradient = crate::RgbGradient::triangle(points, &values); + } + Ok(()) + } + + /// The band between two partner rings (same length), as one quad per edge, or for a twisted + /// wall with `cuts` (see [`cuts`]) as `cuts[i]` quads along edge `i`. `outward` is the + /// direction (before rotation) the band faces away from the solid, for orienting normals. + pub(super) fn band(&mut self, a: &[V3], b: &[V3], cuts: &[usize], outward: impl Fn(usize) -> V3) -> Result<(), String> { + let n = a.len().min(b.len()); + // At 0 and 1, the corners themselves. + let lerp = |p: V3, q: V3, s: f64| { + if s <= 0.0 { + p + } else if s >= 1.0 { + q + } else { + V3(p.0 + (q.0 - p.0) * s, p.1 + (q.1 - p.1) * s, p.2 + (q.2 - p.2) * s) + } + }; + let order = if cuts.is_empty() { Order::Depth } else { Order::Twisted }; + for i in 0..n { + let j = (i + 1) % n; + let pieces = cuts.get(i).copied().unwrap_or(1).max(1); + for k in 0..pieces { + let (s0, s1) = (k as f64 / pieces as f64, (k + 1) as f64 / pieces as f64); + let quad = vec![lerp(a[i], a[j], s0), lerp(a[i], a[j], s1), lerp(b[i], b[j], s1), lerp(b[i], b[j], s0)]; + // A twist bends the strip between two slices into a saddle: its normal is the whole + // quad's average (Newell's), not one corner's triangle, which leans far to one side. + let mut normal = newell(std::slice::from_ref(&quad)); + if normal.dot(normal) < 1e-18 { + continue; + } + if normal.dot(self.camera.rotate(outward(i))) < 0.0 { + normal = V3(-normal.0, -normal.1, -normal.2); + } + self.face(vec![quad], normal, order)?; + } + } + Ok(()) + } + + pub(super) fn finish(mut self) -> Vec { + self.facets.sort_by(|a, b| a.depth.total_cmp(&b.depth)); + self.facets + } +} + +/// How many pieces each edge of a twisted wall is drawn in, from two neighboring slices of the +/// wall, `a` and `b` (partner rings): as many as make each piece about as long along the edge as +/// the band between the slices is wide there, 1 to [`MAX_PIECES`]. Twisted walls pass in front of +/// one another, and a strip the length of its edge has no one depth where another crosses it. +pub(super) fn cuts(a: &[V3], b: &[V3]) -> Vec { + let n = a.len().min(b.len()); + let length = |p: V3, q: V3| q.sub(p).dot(q.sub(p)).sqrt(); + (0..n) + .map(|i| { + let j = (i + 1) % n; + let (along, across) = (length(a[i], a[j]).max(length(b[i], b[j])), length(a[i], b[i]).max(length(a[j], b[j]))); + let count = (along / across.max(1e-9)).ceil(); + if count.is_finite() { (count as usize).clamp(1, MAX_PIECES) } else { MAX_PIECES } + }) + .collect() +} + +/// Lowers the [`cuts`] of every ring's edges (one list per ring) so that `bands` bands of each +/// make no more than `room` quads in all, or one per edge and band when that is more: the pieces +/// past each edge's first share the room in proportion. +pub(super) fn fit(cuts: &mut [Vec], bands: usize, room: usize) { + let edges = cuts.iter().map(Vec::len).sum::().saturating_mul(bands); + let total = cuts.iter().flatten().sum::().saturating_mul(bands); + let room = room.max(edges); + if total > room { + let share = (room - edges) as f64 / (total - edges) as f64; + for c in cuts.iter_mut().flatten() { + *c = 1 + ((c.saturating_sub(1)) as f64 * share).floor() as usize; + } + } +} + +/// The area-weighted normal of (possibly non-planar) contours, by Newell's method: for contours +/// winding with their filled side on the left as seen from in front, it points to the front. +pub(super) fn newell(contours: &[Vec]) -> V3 { + let mut n = V3::default(); + for c in contours { + for (a, b) in c.iter().zip(c.iter().cycle().skip(1)).take(c.len()) { + let x = a.cross(*b); + n = V3(n.0 + x.0, n.1 + x.1, n.2 + x.2); + } + } + n +} + +pub(super) fn check_bounds(b: Rect) -> Result<(), String> { + if [b.x0, b.y0, b.x1, b.y1].iter().all(|v| v.is_finite() && v.abs() <= 1e8) { + Ok(()) + } else { + Err("3D effects require finite coordinates within 100 million points".into()) + } +} + +/// The art's outward direction at edge `i` of `ring` (before rotation), tilted towards the face +/// it borders by `z` (+1 front, -1 back, 0 a side wall). +pub(super) fn outward(ring: &[Point], i: usize, z: f64) -> V3 { + let n = ring.len(); + let (lx, ly) = left(ring[i], ring[(i + 1) % n]); + V3(-lx, -ly, z) +} + +/// The cross-section of a side: each ring of the shape at depth fraction `t` (0 front, 1 back), +/// turned by the twist and scaled by the taper about the centre of `bounds`. +pub(super) fn section(ring: &[Point], bounds: Rect, o: &Extrude, t: f64) -> Vec { + let turn = Turn::new(bounds, o.twist, o.taper, t); + ring.iter().map(|p| turn.apply(*p)).collect() +} + +/// The twist and taper at depth fraction `t` (0 front, 1 back): a turn by `twist · t` and a +/// scale from 1 to `taper`, about the centre of the bounds. +#[derive(Clone, Copy)] +pub(super) struct Turn { + centre: Point, + pub(super) sin: f64, + pub(super) cos: f64, + pub(super) scale: f64, +} + +impl Turn { + pub(super) fn new(bounds: Rect, twist: f64, taper: f64, t: f64) -> Self { + let (sin, cos) = (twist * t).to_radians().sin_cos(); + Self { centre: bounds.center(), sin, cos, scale: 1.0 + (taper - 1.0) * t } + } + pub(super) fn apply(&self, p: Point) -> Point { + let (x, y) = ((p.x - self.centre.x) * self.scale, (p.y - self.centre.y) * self.scale); + Point::new(self.centre.x + x * self.cos - y * self.sin, self.centre.y + x * self.sin + y * self.cos) + } +} + +/// One step of an extrusion's outline: how far it is inset (negative: outset) and its depth +/// fraction (0 front, 1 back), with the face it borders (+1 front bevel, -1 back bevel, 0 wall). +#[derive(Clone, Copy)] +struct Step { + inset: f64, + t: f64, + face: f64, +} + +/// The outline's steps from the front cap's edge to the back cap's edge. +fn steps(o: &Extrude) -> Vec { + let d = o.depth.max(1e-6); + let mut out = vec![]; + let bevel = o.bevel.filter(|b| b.width > 0.0 || b.height > 0.0); + // A bevel can't be deeper than the extrusion (half of it on each side with both sides). + let sides = if bevel.is_some_and(|b| b.both_sides) { 2.0 } else { 1.0 }; + let height = bevel.map_or(0.0, |b| b.height.min(o.depth / sides)); + // The bevel's cross-section across the width: (share across, share down). + let profile = |b: &Bevel| -> Vec<(f64, f64)> { + let mut pts = vec![(0.0, 0.0)]; + let n = b.repeat.max(1) as f64; + let run = (1.0 - b.space) / n; + let gap = if b.repeat > 1 { b.space / (n - 1.0) } else { 0.0 }; + for k in 0..b.repeat.max(1) { + let start = k as f64 * (run + gap); + for s in b.shape.samples().into_iter().skip(1) { + let (x, y) = b.shape.at(s); + pts.push((start + x * run, (k as f64 + y) / n)); + } + if b.repeat > 1 && k + 1 < b.repeat { + pts.push((start + run + gap, (k as f64 + 1.0) / n)); + } + } + if let Some(last) = pts.last_mut() { + last.0 = 1.0; + } + pts + }; + // With the bevel inside, the cap is inset by its width and the wall is the outline; outside, + // the cap is the outline and the wall is outset by the width. + let wall = match bevel { + Some(b) if !b.inside => -b.width, + _ => 0.0, + }; + // A twist turns a bevel's steps too: each step between two samples is cut into steps that + // turn no further than 6°, as the walls' slices do at most. + let refine = |pts: Vec<(f64, f64)>| -> Vec<(f64, f64)> { + if o.twist == 0.0 { + return pts; + } + let mut out = Vec::with_capacity(pts.len()); + for (k, p) in pts.iter().enumerate() { + if let Some(q) = k.checked_sub(1).and_then(|k| pts.get(k)) { + let turn = (p.1 - q.1).abs() * height / d * o.twist.abs(); + let n = ((turn / 6.0).ceil() as usize).clamp(1, 64); + for m in 1..n { + let s = m as f64 / n as f64; + out.push((q.0 + (p.0 - q.0) * s, q.1 + (p.1 - q.1) * s)); + } + } + out.push(*p); + } + out + }; + if let Some(b) = bevel { + for (x, y) in refine(profile(&b)) { + let inset = if b.inside { b.width * (1.0 - x) } else { -b.width * x }; + out.push(Step { inset, t: y * height / d, face: 1.0 }); + } + } else { + out.push(Step { inset: 0.0, t: 0.0, face: 1.0 }); + } + // Twist and taper bend the walls: cut them into slices. + let slices = if o.twist != 0.0 || o.taper != 1.0 { ((o.twist.abs() / 6.0).ceil() as usize).clamp(8, 64) } else { 1 }; + let back_height = if bevel.is_some_and(|b| b.both_sides) { height } else { 0.0 }; + let (t0, t1) = (height / d, 1.0 - back_height / d); + for k in 0..=slices { + let t = t0 + (t1 - t0) * k as f64 / slices as f64; + out.push(Step { inset: wall, t, face: 0.0 }); + } + if let Some(b) = bevel.filter(|b| b.both_sides) { + for (x, y) in refine(profile(&b)).into_iter().rev() { + let inset = if b.inside { b.width * (1.0 - x) } else { -b.width * x }; + out.push(Step { inset, t: 1.0 - y * height / d, face: -1.0 }); + } + } + out.dedup_by(|a, b| (a.inset - b.inset).abs() < 1e-9 && (a.t - b.t).abs() < 1e-12); + out +} + +/// Extrude & Bevel: the shape's outline extruded backwards, with optional twist, taper, bevels +/// and caps. +pub fn extrude(path: &PathData, rule: FillRule, bounds: Rect, options: Extrude, view: View) -> Result, String> { + check_bounds(bounds)?; + let o = options.sanitized(); + // No depth: the flat shape, seen from either side. + if o.depth <= 0.0 { + return plane(path, rule, bounds, view); + } + let rings = rings(path, rule, bounds.width().max(bounds.height()), MAX_SOLID_POINTS)?; + let camera = Camera::new(view, bounds, o.depth); + let mut b = Builder { camera: &camera, closed: o.cap, facets: vec![] }; + let steps = steps(&o); + let z = |t: f64| o.depth / 2.0 - o.depth * t; + let ring_at = + |ring: &[Point], s: &Step| -> Vec { section(&inset(ring, s.inset), bounds, &o, s.t).iter().map(|p| camera.place(*p, z(s.t))).collect() }; + // Every band of a twisted outline, bevels included, is cut as its walls' first slice is (with + // bevels and no walls, as its shallowest band is), within the facet limit. + let depth_of = |s: &[Step]| s.get(1).zip(s.first()).map_or(0.0, |(b, a)| (b.t - a.t).abs()); + let wall = steps + .windows(2) + .find(|s| s[0].face == 0.0 && s[1].face == 0.0) + .or_else(|| steps.windows(2).filter(|s| depth_of(s) > 1e-12).min_by(|a, c| depth_of(a).total_cmp(&depth_of(c)))); + let mut pieces: Vec> = rings + .iter() + .map(|ring| match wall { + Some([s0, s1]) if o.twist != 0.0 => cuts(&ring_at(ring, s0), &ring_at(ring, s1)), + _ => vec![], + }) + .collect(); + fit(&mut pieces, steps.len().saturating_sub(1), MAX_FACES); + for (ring, pieces) in rings.iter().zip(&pieces) { + let placed: Vec> = steps.iter().map(|s| ring_at(ring, s)).collect(); + for (k, pair) in placed.windows(2).enumerate() { + let ([a, c], Some(s0), Some(s1)) = (pair, steps.get(k), steps.get(k + 1)) else { continue }; + // A band between two front bevel steps borders the front, one reaching a back bevel + // step the back; the rest are walls. + let face = if s1.face < 0.0 { + -1.0 + } else if s0.face > 0.0 && s1.face > 0.0 { + 1.0 + } else { + 0.0 + }; + b.band(a, c, pieces, |i| outward(ring, i, face))?; + } + } + if o.cap { + // A cap turned toward the viewer covers the walls it overlaps, even those nearer than its + // center: twisted walls, and the walls of a fill inside the solid of its stroke. + for (s, n) in [(steps.first(), 1.0), (steps.last(), -1.0)] { + let Some(s) = s else { continue }; + b.face(rings.iter().map(|r| ring_at(r, s)).collect(), camera.rotate(V3(0.0, 0.0, n)), Order::Last)?; + } + } + Ok(b.finish()) +} + +/// Rotate: the flat art turned in 3D about the centre of `bounds` (which may be larger than the +/// shape's own, so that every part of one object turns about the same point). +pub fn plane(path: &PathData, rule: FillRule, bounds: Rect, view: View) -> Result, String> { + check_bounds(bounds)?; + let rings = rings(path, rule, bounds.width().max(bounds.height()), MAX_PLANE_POINTS)?; + let camera = Camera::new(view, bounds, 0.0); + let mut b = Builder { camera: &camera, closed: false, facets: vec![] }; + let contours: Vec> = rings.iter().map(|r| r.iter().map(|p| camera.place(*p, 0.0)).collect()).collect(); + // A plane is one flat face, seen from either side. Its normal is turned toward the viewer, so + // it is drawn last, as a cap is, and planes drawn in one order keep the order they are given in. + let up = camera.rotate(V3(0.0, 0.0, 1.0)); + let n = if camera.facing(V3::default(), up) < 0.0 { V3(-up.0, -up.1, -up.2) } else { up }; + b.face(contours, n, Order::Last)?; + Ok(b.finish()) +} + +/// The distance from `p` to the segment `a`–`b`. +pub(super) fn segment_distance(p: Point, a: Point, b: Point) -> f64 { + let (dx, dy) = (b.x - a.x, b.y - a.y); + let l2 = dx * dx + dy * dy; + let t = if l2 < 1e-18 { 0.0 } else { (((p.x - a.x) * dx + (p.y - a.y) * dy) / l2).clamp(0.0, 1.0) }; + p.distance(Point::new(a.x + dx * t, a.y + dy * t)) +} + +/// The distance from `p` to the nearest edge of `rings`, edge by edge. +#[cfg(test)] +pub(super) fn edge_distance(rings: &[Vec], p: Point) -> f64 { + rings + .iter() + .flat_map(|r| r.iter().zip(r.iter().cycle().skip(1)).take(r.len())) + .map(|(a, b)| segment_distance(p, *a, *b)) + .fold(f64::INFINITY, f64::min) +} diff --git a/crates/three-d/src/solid/inflate.rs b/crates/three-d/src/solid/inflate.rs new file mode 100644 index 000000000..be5aec370 --- /dev/null +++ b/crates/three-d/src/solid/inflate.rs @@ -0,0 +1,305 @@ +use super::build::{Builder, Camera, Order, Turn, check_bounds, cuts, fit, newell, outward}; +use super::membrane::{Edges, MAX_MATERIAL_SEGMENTS, Membrane, cell_triangles, clip, filled, material}; +use super::rings::{MAX_SOLID_POINTS, rings, signed_area}; +use super::*; + +/// A shape to inflate: its rings and the membrane under them on the grid of the bounds. +struct Puff { + rings: Vec>, + /// The index of the rings' edges. + edges: Edges, + cell: f64, + nx: usize, + ny: usize, + m: Membrane, + /// How high the front rises at its peak. + height: f64, +} + +/// The rings of `path` and its membrane on the grid of `bounds`; `None` for a shape with no area. +fn puff(path: &PathData, rule: FillRule, bounds: Rect, o: &Inflate) -> Result, String> { + let rings = rings(path, rule, bounds.width().max(bounds.height()), MAX_SOLID_POINTS)?; + if rings.is_empty() { + return Ok(None); + } + let size = bounds.width().max(bounds.height()); + if size <= 0.0 { + return Ok(None); + } + const CELLS: usize = 48; + let cell = size / CELLS as f64; + let (nx, ny) = (((bounds.width() / cell).ceil() as usize).max(1), ((bounds.height() / cell).ceil() as usize).max(1)); + let edges = Edges::new(&rings); + let m = Membrane::solve(&rings, &edges, bounds.origin(), cell, nx + 1, ny + 1); + // A circle of radius r puffs into a half sphere of radius r at 50 % (its membrane peaks at r²/4). + let height = 2.0 * m.peak.sqrt() * o.volume / 50.0; + Ok(Some(Puff { rings, edges, cell, nx, ny, m, height })) +} + +/// Inflate: the shape puffed out like a cushion, its edge kept where the art's outline is. +pub fn inflate(path: &PathData, rule: FillRule, bounds: Rect, options: Inflate, view: View) -> Result, String> { + inflate_together(&[(path, rule)], bounds, options, view).pop().unwrap_or(Ok(vec![])) +} + +/// Inflate for several shapes of one object (a fill and the outlines of its strokes): each puffs +/// on its own over the grid of `bounds`, and one eye sees them all, as far away as the highest +/// puff needs with Perspective. One result for each shape, in order. +pub fn inflate_together(shapes: &[(&PathData, FillRule)], bounds: Rect, options: Inflate, view: View) -> Vec, String>> { + let o = options.sanitized(); + let puffs: Vec, String>> = + shapes.iter().map(|(path, rule)| check_bounds(bounds).and_then(|()| puff(path, *rule, bounds, &o))).collect(); + let reach = puffs.iter().flatten().flatten().map(|p| p.height).fold(0.0, f64::max); + puffs + .into_iter() + .map(|p| match p { + Ok(Some(p)) => inflated(p, bounds, o, view, reach), + Ok(None) => Ok(vec![]), + Err(e) => Err(e), + }) + .collect() +} + +/// The facets of one inflated shape, for a camera that allows for a puff `reach` high (or the +/// shape's own, if higher): `reach` is the highest puff of the shapes drawn together. +fn inflated(p: Puff, bounds: Rect, o: Inflate, view: View, reach: f64) -> Result, String> { + let Puff { rings, edges, cell, nx, ny, m, height } = p; + let outline: Vec<(Rect, vectorcraft_geom::SubPath)> = rings + .iter() + .filter_map(|r| { + let b = r.iter().skip(1).fold(Rect::from_points(*r.first()?, *r.first()?), |b, p| b.union_pt(*p)); + Some((b, vectorcraft_geom::SubPath::polyline(r, true))) + }) + .collect(); + // When no shape puffs (no volume, or shapes too narrow for the grid), every front (and a back + // that puffs too) lies flat in one plane: its cells are drawn last, as caps are. If one shape + // puffs, a flat front can be behind its puff, so the cells keep their depths. + let flat = reach <= 0.0; + let lift = |p: Point| { + if m.peak <= 0.0 || edges.near(p, 1e-9) { + return 0.0; + } + height * (m.at(p).max(0.0) / m.peak).sqrt() + }; + let normal = |p: Point| { + let (dx, dy) = m.gradient(p); + let (_, inward, distance) = edges.closest(p); + let t = (distance / (2.0 * cell)).clamp(0.0, 1.0); + let weight = (1.0 - t).powi(2) * (1.0 + 2.0 * t); + let gradient = Vec2::new(dx, dy); + let slope = gradient.dot(inward).max(0.0); + // The zero boundary condition also fixes the tangential derivative at the outline. + // Extrapolated outside solver nodes cannot supply it reliably. Blend to that known + // direction and its limiting phi near the rim, without showing the grid in highlights. + let gradient = gradient * (1.0 - weight) + inward * (slope * weight); + let phi = m.at(p).max(0.0) * (1.0 - weight) + slope * distance * weight; + // Multiplying the height-function normal by 2√(peak φ) avoids division by + // zero at the rim, where the limiting normal is horizontal. + let n = V3(-height * gradient.x, -height * gradient.y, 2.0 * (m.peak * phi).sqrt()).unit(); + if n.dot(n) > 0.5 { n } else { V3(0.0, 0.0, 1.0) } + }; + let smooth = view.shade && view.has_colored_lights(); + let camera = Camera::new(view, bounds, o.depth + 2.0 * height.max(reach)); + // Without its back cap the object is a shell open behind, its inside showing. + let mut b = Builder { camera: &camera, closed: o.cap, facets: vec![] }; + let half = o.depth / 2.0; + // The back: turned and scaled by the twist and taper, its puff too. + let back_turn = Turn::new(bounds, o.twist, o.taper, 1.0); + // The grid's cells, those the surface rises steeply across (near the outline) cut into + // halves or quarters each way, so the steep rim is drawn in facets as fine as its curve. + // Every corner sits on one fine lattice (a quarter cell), computed the same way wherever it + // is used, and a cell's border takes its finer neighbour's lattice points too: where cells + // of different fineness meet, their edges follow the surface through the same points, with + // no cracks. + const FINE: usize = 4; + let unit = cell / FINE as f64; + let at = |qx: usize, qy: usize| Point::new(bounds.x0 + qx as f64 * unit, bounds.y0 + qy as f64 * unit); + let base = |i: usize, j: usize| Rect::from_points(at(i * FINE, j * FINE), at((i + 1) * FINE, (j + 1) * FINE)); + // Segments the booleans of mixed cells may take in all ([`MAX_MATERIAL_SEGMENTS`]): each one + // builds the topology of every outline ring near its cell, so a self-crossing outline cost its + // crossings again in every cell it spans. Past it, cells take their clipped pieces. + let material_left = std::cell::Cell::new(MAX_MATERIAL_SEGMENTS); + // The shape's part of a rectangle, as rings. + let cut = |r: Rect| -> Vec> { + // A cell no edge crosses is wholly inside the shape or wholly outside it. + if !edges.touch(r) { + let corners = vec![Point::new(r.x0, r.y0), Point::new(r.x1, r.y0), Point::new(r.x1, r.y1), Point::new(r.x0, r.y1)]; + return if filled(&rings, r.center()) { vec![corners] } else { vec![] }; + } + let pieces: Vec> = rings.iter().map(|ring| clip(ring, r)).filter(|c| !c.is_empty()).collect(); + if pieces.iter().map(|p| signed_area(p)).sum::() <= 1e-9 * r.area() { + return vec![]; + } + // Where a hole's piece cancels part of an outline's, that only holds flat: lifted onto + // the surface the two no longer meet. Such a cell takes the shape's own part of it. + let mixed = pieces.iter().any(|p| signed_area(p) < 0.0) && pieces.iter().any(|p| signed_area(p) > 0.0); + if !mixed { + return pieces; + } + let near: usize = + outline.iter().filter(|(b, _)| b.x0 <= r.x1 && r.x0 <= b.x1 && b.y0 <= r.y1 && r.y0 <= b.y1).map(|(_, sp)| sp.anchors.len()).sum(); + match material_left.get().checked_sub(near) { + Some(left) => { + material_left.set(left); + material(&outline, r).unwrap_or(pieces) + } + None => pieces, + } + }; + // How many times each cell is halved each way (0, 1 or 2), from how steeply it rises. + let levels: Vec = (0..nx * ny) + .map(|k| { + let pieces = cut(base(k % nx, k / nx)); + let (lo, hi) = pieces.iter().flatten().map(|q| lift(*q)).fold((f64::INFINITY, f64::NEG_INFINITY), |(lo, hi), h| (lo.min(h), hi.max(h))); + let rise = if pieces.is_empty() { 0.0 } else { (hi - lo) / cell }; + if rise > 6.0 { + 2 + } else if rise > 2.0 { + 1 + } else { + 0 + } + }) + .collect(); + let level = + |i: isize, j: isize| -> u32 { if i < 0 || j < 0 || i as usize >= nx || j as usize >= ny { 0 } else { levels[j as usize * nx + i as usize] } }; + for j in 0..ny { + for i in 0..nx { + let l = level(i as isize, j as isize); + let step = FINE >> l; + // The lattice spacing along each side of this cell: its own or its neighbour's. + let (ii, jj) = (i as isize, j as isize); + let side_step = |dl: u32| FINE >> l.max(dl); + let (left, right, top, bottom) = + (side_step(level(ii - 1, jj)), side_step(level(ii + 1, jj)), side_step(level(ii, jj - 1)), side_step(level(ii, jj + 1))); + let (bx0, by0, bx1, by1) = (i * FINE, j * FINE, (i + 1) * FINE, (j + 1) * FINE); + for sy in (by0..by1).step_by(step) { + for sx in (bx0..bx1).step_by(step) { + let r = Rect::from_points(at(sx, sy), at(sx + step, sy + step)); + let pieces = cut(r); + if pieces.is_empty() { + continue; + } + // Lattice points along the cell's own border, between each edge's ends. + let stitch = |a: Point, b: Point| -> Vec { + let tol = 1e-7 * unit; + let along = |fixed_x: Option, fixed_y: Option, every: usize| -> Vec { + let mut out = vec![]; + match (fixed_x, fixed_y) { + (Some(qx), None) => { + let (lo, hi) = (a.y.min(b.y), a.y.max(b.y)); + for qy in (by0..=by1).step_by(every) { + let p = at(qx, qy); + if p.y > lo + tol && p.y < hi - tol { + out.push(p); + } + } + } + (None, Some(qy)) => { + let (lo, hi) = (a.x.min(b.x), a.x.max(b.x)); + for qx in (bx0..=bx1).step_by(every) { + let p = at(qx, qy); + if p.x > lo + tol && p.x < hi - tol { + out.push(p); + } + } + } + _ => {} + } + if (a.x > b.x) || (a.y > b.y) { + out.reverse(); + } + out + }; + let on = |v: f64, w: f64| (v - w).abs() <= tol; + let (x0, x1, y0, y1) = (at(bx0, by0).x, at(bx1, by1).x, at(bx0, by0).y, at(bx1, by1).y); + if on(a.x, x0) && on(b.x, x0) { + along(Some(bx0), None, left) + } else if on(a.x, x1) && on(b.x, x1) { + along(Some(bx1), None, right) + } else if on(a.y, y0) && on(b.y, y0) { + along(None, Some(by0), top) + } else if on(a.y, y1) && on(b.y, y1) { + along(None, Some(by1), bottom) + } else { + vec![] + } + }; + let pieces: Vec> = pieces + .iter() + .map(|p| { + let mut out = Vec::with_capacity(p.len()); + for (a, b) in p.iter().zip(p.iter().cycle().skip(1)).take(p.len()) { + out.push(*a); + out.extend(stitch(*a, *b)); + } + out + }) + .collect(); + if height > 0.0 + && smooth + && let Some(triangles) = cell_triangles(&pieces) + { + for points in triangles { + let front = points.map(|q| camera.place(q, half + lift(q))); + b.smooth_triangle(front, points.map(|q| camera.rotate(normal(q))), false)?; + if o.both_sides && o.cap { + let back = points.map(|q| camera.place(back_turn.apply(q), -half - lift(q) * back_turn.scale)); + let normals = points.map(|q| { + let n = normal(q); + // A uniform taper scales all three coordinates equally; + // twist turns the back normal about the shape's z-axis. + camera.rotate(V3(back_turn.cos * n.0 - back_turn.sin * n.1, back_turn.sin * n.0 + back_turn.cos * n.1, -n.2)) + }); + b.smooth_triangle(back, normals, true)?; + } + } + continue; + } + let lifts: Vec> = pieces.iter().map(|p| p.iter().map(|q| lift(*q)).collect()).collect(); + // Each piece lifted onto the surface. Its normal is the lifted piece's own average + // (Newell's), not the slope at one point: near the edge the surface turns + // vertical, and a point sample there misjudges which way a whole piece faces. + let front: Vec> = + pieces.iter().zip(&lifts).map(|(p, h)| p.iter().zip(h).map(|(q, h)| camera.place(*q, half + h)).collect()).collect(); + let n = newell(&front); + let order = if flat { Order::Last } else { Order::Depth }; + b.face(front, n, order)?; + if o.both_sides && o.cap { + // The back winds the same way round as the front, so its outward normal is + // reversed. + let back: Vec> = pieces + .iter() + .zip(&lifts) + .map(|(p, h)| p.iter().zip(h).map(|(q, h)| camera.place(back_turn.apply(*q), -half - h * back_turn.scale)).collect()) + .collect(); + let n = newell(&back); + b.face(back, V3(-n.0, -n.1, -n.2), order)?; + } + } + } + } + } + if !o.both_sides && o.cap { + let back: Vec> = rings.iter().map(|r| r.iter().map(|p| camera.place(back_turn.apply(*p), -half)).collect()).collect(); + // The back turned toward the viewer covers the walls it overlaps, even those nearer than + // its center: twisted walls, and the walls of a fill inside the solid of its stroke. + b.face(back, camera.rotate(V3(0.0, 0.0, -1.0)), Order::Last)?; + } + if o.depth > 0.0 { + // The sides, cut into slices that follow a twist or taper. + let slices = if o.twist != 0.0 || o.taper != 1.0 { ((o.twist.abs() / 6.0).ceil() as usize).clamp(8, 64) } else { 1 }; + let at = |ring: &[Point], k: usize| -> Vec { + let t = k as f64 / slices as f64; + let turn = Turn::new(bounds, o.twist, o.taper, t); + ring.iter().map(|p| camera.place(turn.apply(*p), half - o.depth * t)).collect() + }; + let mut pieces: Vec> = rings.iter().map(|r| if o.twist != 0.0 { cuts(&at(r, 0), &at(r, 1)) } else { vec![] }).collect(); + fit(&mut pieces, slices, MAX_FACES.saturating_sub(b.facets.len())); + for (ring, pieces) in rings.iter().zip(&pieces) { + for k in 0..slices { + b.band(&at(ring, k), &at(ring, k + 1), pieces, |i| outward(ring, i, 0.0))?; + } + } + } + Ok(b.finish()) +} diff --git a/crates/three-d/src/solid/membrane.rs b/crates/three-d/src/solid/membrane.rs new file mode 100644 index 000000000..296b659df --- /dev/null +++ b/crates/three-d/src/solid/membrane.rs @@ -0,0 +1,400 @@ +use super::build::segment_distance; +use super::rings::{MAX_PLANE_POINTS, contains, left, rings, signed_area}; +use super::*; + +/// The edges of `rings` filed in a grid of square buckets, so the nearest edge to a point is +/// found among the buckets around it rather than among all the edges. +pub(super) struct Edges { + edges: Vec<(Point, Point)>, + normals: Vec<(Vec2, Vec2)>, + origin: Point, + size: f64, + n: usize, + buckets: Vec>, +} + +impl Edges { + pub(super) fn new(rings: &[Vec]) -> Self { + const N: usize = 64; + let edges: Vec<(Point, Point)> = + rings.iter().flat_map(|r| r.iter().zip(r.iter().cycle().skip(1)).take(r.len()).map(|(a, b)| (*a, *b))).collect(); + let normals = rings + .iter() + .flat_map(|r| { + let normal = |i: usize| { + let (x, y) = left(r[i], r[(i + 1) % r.len()]); + Vec2::new(x, y) + }; + (0..r.len()).map(move |i| { + let n = normal(i); + let at = |other: Vec2| if n.dot(other) > 0.7 { (n + other).normalize() } else { n }; + (at(normal((i + r.len() - 1) % r.len())), at(normal((i + 1) % r.len()))) + }) + }) + .collect(); + let bounds = edges.iter().fold(None, |b: Option, (a, c)| Some(b.map_or(Rect::from_points(*a, *c), |b| b.union_pt(*a).union_pt(*c)))); + let bounds = bounds.unwrap_or(Rect::ZERO); + let size = (bounds.width().max(bounds.height()) / N as f64).max(1e-9); + let mut buckets = vec![vec![]; N * N]; + let index = |v: f64, o: f64| (((v - o) / size).floor().max(0.0) as usize).min(N - 1); + for (k, (a, b)) in edges.iter().enumerate() { + let (i0, i1) = (index(a.x.min(b.x), bounds.x0), index(a.x.max(b.x), bounds.x0)); + let (j0, j1) = (index(a.y.min(b.y), bounds.y0), index(a.y.max(b.y), bounds.y0)); + for j in j0..=j1 { + for i in i0..=i1 { + buckets[j * N + i].push(k as u32); + } + } + } + Self { edges, normals, origin: bounds.origin(), size, n: N, buckets } + } + + /// The bucket row or column of coordinate `v`, from the grid's origin `o` on that axis. + fn cell(&self, v: f64, o: f64) -> isize { + (((v - o) / self.size).floor() as isize).clamp(0, self.n as isize - 1) + } + + /// Is `p` within `tolerance` of an edge? Only the buckets `p` could share an edge with are + /// looked at. + pub(super) fn near(&self, p: Point, tolerance: f64) -> bool { + self.touch(Rect::new(p.x - tolerance, p.y - tolerance, p.x + tolerance, p.y + tolerance)) && { + let n = self.n as isize; + let (i0, i1) = (self.cell(p.x - tolerance, self.origin.x), self.cell(p.x + tolerance, self.origin.x)); + let (j0, j1) = (self.cell(p.y - tolerance, self.origin.y), self.cell(p.y + tolerance, self.origin.y)); + (j0..=j1).any(|j| { + (i0..=i1).any(|i| { + self.buckets[(j * n + i) as usize] + .iter() + .filter_map(|k| self.edges.get(*k as usize)) + .any(|(a, b)| segment_distance(p, *a, *b) <= tolerance) + }) + }) + } + } + + /// Call `f` with each edge whose bounding box overlaps `r` (an edge in several of the + /// buckets `r` covers comes more than once). + fn each(&self, r: Rect, mut f: impl FnMut(Point, Point)) { + let n = self.n as isize; + for j in self.cell(r.y0, self.origin.y)..=self.cell(r.y1, self.origin.y) { + for i in self.cell(r.x0, self.origin.x)..=self.cell(r.x1, self.origin.x) { + for (a, b) in self.buckets[(j * n + i) as usize].iter().filter_map(|k| self.edges.get(*k as usize)) { + if a.x.min(b.x) <= r.x1 && a.x.max(b.x) >= r.x0 && a.y.min(b.y) <= r.y1 && a.y.max(b.y) >= r.y0 { + f(*a, *b); + } + } + } + } + } + + /// Does any edge pass through `r` (its bounding box, at least, overlapping it)? + pub(super) fn touch(&self, r: Rect) -> bool { + let n = self.n as isize; + let (i0, i1) = (self.cell(r.x0, self.origin.x), self.cell(r.x1, self.origin.x)); + let (j0, j1) = (self.cell(r.y0, self.origin.y), self.cell(r.y1, self.origin.y)); + (j0..=j1).any(|j| { + (i0..=i1).any(|i| { + self.buckets[(j * n + i) as usize] + .iter() + .filter_map(|k| self.edges.get(*k as usize)) + .any(|(a, b)| a.x.min(b.x) <= r.x1 && a.x.max(b.x) >= r.x0 && a.y.min(b.y) <= r.y1 && a.y.max(b.y) >= r.y0) + }) + }) + } + + /// The distance from `p` to the nearest edge: the buckets are searched in growing squares + /// around `p`'s, until no unsearched bucket could hold a nearer edge. + pub(super) fn distance(&self, p: Point) -> f64 { + self.closest(p).2 + } + + /// Closest outline point, its inward normal and distance. Curved outline samples share + /// corner normals, while a sharp source corner keeps each edge's own direction. + pub(super) fn closest(&self, p: Point) -> (Point, Vec2, f64) { + let (fx, fy) = ((p.x - self.origin.x) / self.size, (p.y - self.origin.y) / self.size); + if !fx.is_finite() || !fy.is_finite() { + return (p, Vec2::ZERO, f64::INFINITY); + } + let n = self.n as isize; + let (ci, cj) = ((fx.floor() as isize).clamp(0, n - 1), (fy.floor() as isize).clamp(0, n - 1)); + // How far `p` is outside the grid: the nearest point of the grid to it lies in `p`'s bucket. + let outside = (fx.clamp(0.0, n as f64) - fx).hypot(fy.clamp(0.0, n as f64) - fy) * self.size; + let mut best = f64::INFINITY; + let mut nearest = (p, Vec2::ZERO); + for r in 0..n { + for j in (cj - r)..=(cj + r) { + for i in (ci - r)..=(ci + r) { + if (j - cj).abs() != r && (i - ci).abs() != r || i < 0 || j < 0 || i >= n || j >= n { + continue; + } + for k in &self.buckets[(j * n + i) as usize] { + if let Some((a, b)) = self.edges.get(*k as usize) { + let d = *b - *a; + let t = if d.hypot2() < 1e-18 { 0.0 } else { ((p - *a).dot(d) / d.hypot2()).clamp(0.0, 1.0) }; + let q = Point::new(a.x + d.x * t, a.y + d.y * t); + let distance = p.distance(q); + if distance < best { + best = distance; + let n = self.normals.get(*k as usize).map_or(Vec2::ZERO, |(a, b)| *a * (1.0 - t) + *b * t); + nearest = (q, if n.hypot2() > 0.0 { n.normalize() } else { n }); + } + } + } + } + } + // Buckets r + 1 away are at least r buckets from any point of `p`'s own, and that + // far across the grid from a `p` outside it (at right angles to its way in, at worst). + if best <= outside.hypot(r as f64 * self.size) { + break; + } + } + (nearest.0, nearest.1, best) + } +} + +/// Is `p` filled by `rings` (nonzero winding, as they are oriented)? +pub(super) fn filled(rings: &[Vec], p: Point) -> bool { + rings.iter().filter(|r| contains(r, p)).count() % 2 == 1 +} + +/// `ring` cut to the rectangle `cell` (Sutherland–Hodgman: the rectangle is convex, the ring +/// needn't be). +pub(super) fn clip(ring: &[Point], cell: Rect) -> Vec { + let mut pts = ring.to_vec(); + let sides: [fn(Point, Rect) -> f64; 4] = [|p, c| p.x - c.x0, |p, c| c.x1 - p.x, |p, c| p.y - c.y0, |p, c| c.y1 - p.y]; + for side in sides { + if pts.is_empty() { + break; + } + let mut out = vec![]; + for (a, b) in pts.iter().zip(pts.iter().cycle().skip(1)).take(pts.len()) { + let (da, db) = (side(*a, cell), side(*b, cell)); + if da >= 0.0 { + out.push(*a); + } + if (da >= 0.0) != (db >= 0.0) { + out.push(a.lerp(*b, da / (da - db))); + } + } + pts = out; + } + pts.dedup_by(|a, b| a.distance(*b) < 1e-9); + if pts.len() >= 3 && signed_area(&pts).abs() > 1e-12 { pts } else { vec![] } +} + +/// The shape's surface as a membrane under pressure, held at zero along its outline (and its +/// holes): the solution of ∇²φ = −1 on a grid, with each node's distance to the outline along the +/// grid lines taken into account where it is less than a grid step (Shortley–Weller), so the +/// outline needn't fall on the grid. √φ rises vertically from the outline and rounds smoothly +/// inside, without the creases a distance-to-the-edge height has along a shape's middle. +pub(super) struct Membrane { + origin: Point, + step: f64, + nx: usize, + ny: usize, + /// φ at each node; outside the shape, extrapolated (negative) from the inside, so + /// interpolation crosses zero near the outline. + phi: Vec, + pub(super) peak: f64, +} + +impl Membrane { + /// `edges` is the index of `rings`' edges (`Edges::new(rings)`). + pub(super) fn solve(rings: &[Vec], edges: &Edges, origin: Point, step: f64, nx: usize, ny: usize) -> Self { + let node = |i: usize, j: usize| Point::new(origin.x + i as f64 * step, origin.y + j as f64 * step); + let idx = |i: usize, j: usize| j * nx + i; + // A node on the outline itself is held at zero with the outside. + let inside: Vec = (0..nx * ny) + .map(|k| { + let p = node(k % nx, k / nx); + filled(rings, p) && !edges.near(p, 1e-6 * step) + }) + .collect(); + // The nearest crossing of the outline from `p` along (dx, dy), within one step. + let crossing = |p: Point, dx: f64, dy: f64| -> f64 { + let mut best = step; + let reach = Rect::from_points(p, Point::new(p.x + dx * step, p.y + dy * step)); + edges.each(reach, |a, b| { + let t = if dy == 0.0 { + if (a.y > p.y) == (b.y > p.y) { + return; + } + (a.x + (p.y - a.y) / (b.y - a.y) * (b.x - a.x) - p.x) * dx + } else { + if (a.x > p.x) == (b.x > p.x) { + return; + } + (a.y + (p.x - a.x) / (b.x - a.x) * (b.y - a.y) - p.y) * dy + }; + if t > 1e-9 * step && t < best { + best = t; + } + }); + best + }; + // Each inside node's four neighbours: (weight, neighbour inside or None at the outline). + let mut stencil: Vec<[(f64, Option); 4]> = vec![[(0.0, None); 4]; nx * ny]; + for j in 0..ny { + for i in 0..nx { + let k = idx(i, j); + if !inside[k] { + continue; + } + let p = node(i, j); + let dirs = [(1.0, 0.0), (-1.0, 0.0), (0.0, 1.0), (0.0, -1.0)]; + let h: Vec = dirs.iter().map(|(dx, dy)| crossing(p, *dx, *dy)).collect(); + for (d, (dx, dy)) in dirs.iter().enumerate() { + let opp = h[d ^ 1]; + let w = 2.0 / (h[d] * (h[d] + opp)); + let (ni, nj) = (i as isize + *dx as isize, j as isize + *dy as isize); + let neighbour = (h[d] >= step && ni >= 0 && nj >= 0 && (ni as usize) < nx && (nj as usize) < ny) + .then(|| idx(ni as usize, nj as usize)) + .filter(|n| inside[*n]); + stencil[k][d] = (w, neighbour); + } + } + } + // Successive over-relaxation. + let mut phi = vec![0.0; nx * ny]; + let tolerance = 1e-9 * step * step * (nx * ny) as f64; + for _ in 0..2000 { + let mut change: f64 = 0.0; + for k in 0..nx * ny { + if !inside[k] { + continue; + } + let (mut sum, mut weight) = (1.0, 0.0); + for (w, n) in stencil[k] { + weight += w; + sum += w * n.map_or(0.0, |n| phi[n]); + } + let next = phi[k] + 1.8 * (sum / weight - phi[k]); + change = change.max((next - phi[k]).abs()); + phi[k] = next; + } + if change < tolerance { + break; + } + } + let peak = phi.iter().copied().fold(0.0, f64::max); + // Outside nodes next to the shape: φ continued past the outline at the slope of their + // inside neighbours, so values interpolated between them reach zero at the outline. + let mut ext = phi.clone(); + for j in 0..ny { + for i in 0..nx { + let k = idx(i, j); + if inside[k] { + continue; + } + let d_out = edges.distance(node(i, j)); + let mut slopes = vec![]; + for (di, dj) in [(1isize, 0isize), (-1, 0), (0, 1), (0, -1)] { + let (ni, nj) = (i as isize + di, j as isize + dj); + if ni < 0 || nj < 0 || ni as usize >= nx || nj as usize >= ny || !inside[idx(ni as usize, nj as usize)] { + continue; + } + let d_in = edges.distance(node(ni as usize, nj as usize)).max(1e-9); + slopes.push(phi[idx(ni as usize, nj as usize)] / d_in); + } + if !slopes.is_empty() { + ext[k] = -d_out * slopes.iter().sum::() / slopes.len() as f64; + } + } + } + Self { origin, step, nx, ny, phi: ext, peak } + } + + /// φ at `p`, interpolated between the four nodes around it. + pub(super) fn at(&self, p: Point) -> f64 { + self.interpolate(p, |i, j| self.phi.get(j * self.nx + i).copied().unwrap_or(0.0)) + } + + /// Interpolate central node derivatives, rather than differentiating each bilinear cell. + /// The latter jumps at every grid line and shows the solver's grid in the lighting. + pub(super) fn gradient(&self, p: Point) -> (f64, f64) { + let value = |i, j| self.phi.get(j * self.nx + i).copied().unwrap_or(0.0); + let dx = self.interpolate(p, |i, j| { + let (a, b) = (i.saturating_sub(1), (i + 1).min(self.nx - 1)); + (value(b, j) - value(a, j)) / ((b - a).max(1) as f64 * self.step) + }); + let dy = self.interpolate(p, |i, j| { + let (a, b) = (j.saturating_sub(1), (j + 1).min(self.ny - 1)); + (value(i, b) - value(i, a)) / ((b - a).max(1) as f64 * self.step) + }); + (dx, dy) + } + + fn interpolate(&self, p: Point, v: impl Fn(usize, usize) -> f64) -> f64 { + let fx = ((p.x - self.origin.x) / self.step).clamp(0.0, (self.nx - 1) as f64); + let fy = ((p.y - self.origin.y) / self.step).clamp(0.0, (self.ny - 1) as f64); + let (i, j) = ((fx.floor() as usize).min(self.nx.saturating_sub(2)), (fy.floor() as usize).min(self.ny.saturating_sub(2))); + let (tx, ty) = (fx - i as f64, fy - j as f64); + if self.nx < 2 || self.ny < 2 { + return v(0, 0); + } + let top = v(i, j) * (1.0 - tx) + v(i + 1, j) * tx; + let bottom = v(i, j + 1) * (1.0 - tx) + v(i + 1, j + 1) * tx; + top * (1.0 - ty) + bottom * ty + } +} + +/// A bounded ear clipping pass on a simple cell contour. Keep its collinear border nodes: +/// those shared samples can have different heights once lifted onto the membrane. A tiny +/// cell containing a complete hole, or a complex/degenerate contour, keeps its original fill. +pub(super) fn cell_triangles(pieces: &[Vec]) -> Option> { + let mut triangles = vec![]; + for ring in pieces { + if ring.len() > 128 || signed_area(ring) <= 0.0 { + return None; + } + let mut ring = ring.clone(); + ring.dedup_by(|a, b| a.distance(*b) < 1e-9); + if ring.first().zip(ring.last()).is_some_and(|(a, b)| a.distance(*b) < 1e-9) { + ring.pop(); + } + while ring.len() > 3 { + let ear = (0..ring.len()).find(|&i| { + let n = ring.len(); + let (a, b, c) = (ring[(i + n - 1) % n], ring[i], ring[(i + 1) % n]); + if (b - a).cross(c - a) <= 1e-12 { + return false; + } + !ring.iter().enumerate().any(|(j, p)| { + j != i + && j != (i + n - 1) % n + && j != (i + 1) % n + && (b - a).cross(*p - a) >= -1e-12 + && (c - b).cross(*p - b) >= -1e-12 + && (a - c).cross(*p - c) >= -1e-12 + }) + })?; + let n = ring.len(); + triangles.push([ring[(ear + n - 1) % n], ring[ear], ring[(ear + 1) % n]]); + ring.remove(ear); + } + if let [a, b, c] = ring.as_slice() { + triangles.push([*a, *b, *c]); + } + } + Some(triangles) +} + +/// Most outline segments the booleans of one inflated object's mixed cells take together. +pub(super) const MAX_MATERIAL_SEGMENTS: usize = 4_000_000; + +/// The part of `outline` (its rings, each with its bounds) inside `cell`, as rings (holes running +/// the other way). Only the rings whose bounds meet the cell's go into the boolean. +pub(super) fn material(outline: &[(Rect, vectorcraft_geom::SubPath)], cell: Rect) -> Option>> { + let r = PathData::single(vectorcraft_geom::SubPath::polyline( + &[Point::new(cell.x0, cell.y0), Point::new(cell.x1, cell.y0), Point::new(cell.x1, cell.y1), Point::new(cell.x0, cell.y1)], + true, + )); + let near = PathData::new( + outline.iter().filter(|(b, _)| b.x0 <= cell.x1 && cell.x0 <= b.x1 && b.y0 <= cell.y1 && cell.y0 <= b.y1).map(|(_, s)| s.clone()).collect(), + ); + let rule = FillRule::NonZero; + let part = + vectorcraft_pathops::try_boolean(&r, rule, &near, rule, vectorcraft_pathops::BoolOp::Intersect, vectorcraft_pathops::DEFAULT_PRECISION) + .ok()?; + rings(&part, rule, cell.width(), MAX_PLANE_POINTS).ok() +} diff --git a/crates/three-d/src/solid/rings.rs b/crates/three-d/src/solid/rings.rs new file mode 100644 index 000000000..57335e4c6 --- /dev/null +++ b/crates/three-d/src/solid/rings.rs @@ -0,0 +1,163 @@ +use super::*; + +/// Outline points Extrude & Bevel and Inflate sample at most: every point costs walls, bevel +/// bands and distance tests. +pub(super) const MAX_SOLID_POINTS: usize = 2048; + +/// Outline points Rotate samples at most: it only turns the outline, so an outlined paragraph fits. +pub(super) const MAX_PLANE_POINTS: usize = 65536; + +/// Pieces a twisted wall's strip is cut into along one edge at most. +pub(super) const MAX_PIECES: usize = 64; + +/// How many straight pieces follow `c` within `tolerance` (Wang's bound on a cubic's distance from +/// its chords), 1 to 32. +fn pieces(c: vectorcraft_geom::CubicBez, tolerance: f64) -> u32 { + let d = |a: Point, b: Point, c: Point| (a.to_vec2() - 2.0 * b.to_vec2() + c.to_vec2()).hypot(); + let bend = d(c.p0, c.p1, c.p2).max(d(c.p1, c.p2, c.p3)); + let n = (0.75 * bend / tolerance).sqrt().ceil(); + if n.is_finite() { (n as u32).clamp(1, 32) } else { 32 } +} + +/// The rings of a shape: closed polygons with the filled side on their left, (-dy, dx): outer +/// outlines run one way and holes the other, so nonzero winding fills them as the art is. Curves +/// are followed within a tolerance of 1/2000 of `size` (at least 0.05 pt), so small curves take +/// few points and large ones more; more than `limit` points in all is refused. +/// +/// Under the nonzero `rule` an outline with filled art on both sides is no edge of the shape and +/// is left out: a contour inside another that runs the same way doesn't cut a hole, as it would +/// if nesting alone decided (#1724). +pub(super) fn rings(path: &PathData, rule: FillRule, size: f64, limit: usize) -> Result>, String> { + let tolerance = (size / 2000.0).max(0.05); + let mut count = 0; + let mut out = vec![]; + for sp in &path.subpaths { + let segments = if sp.closed { sp.anchors.len() } else { sp.anchors.len().saturating_sub(1) }; + if segments > limit { + return Err("The shape is too complex for a 3D effect".into()); + } + let mut ring: Vec = vec![]; + if let Some(first) = sp.anchors.first() { + ring.push(first.p); + } + for i in 0..segments { + let curve = sp.segment(i); + let samples = if sp.segment_is_line(i) { 1 } else { pieces(curve, tolerance) }; + count += samples as usize; + if count > limit { + return Err("The shape is too complex for a 3D effect".into()); + } + for k in 1..=samples { + ring.push(curve.eval(f64::from(k) / f64::from(samples))); + } + } + if ring.iter().any(|p| !p.x.is_finite() || !p.y.is_finite() || p.x.abs() > 1e8 || p.y.abs() > 1e8) { + return Err("The shape contains invalid coordinates".into()); + } + ring.dedup_by(|a, b| a.distance(*b) < 1e-9); + while ring.len() > 1 && ring.first().zip(ring.last()).is_some_and(|(a, b)| a.distance(*b) < 1e-9) { + ring.pop(); + } + if ring.len() >= 3 && signed_area(&ring).abs() > 1e-9 { + out.push(ring); + } + } + // Which rings hold which: each ring against the others, its edges walked only where its + // bounds hold the point, and within [`MAX_CONTAINMENT_EDGES`] edges in all. + let too_complex = || String::from("The shape is too complex for a 3D effect"); + let mut work = 0usize; + let mut within = |out: &[Vec], bounds: &[Rect], i: usize| -> Result, String> { + let p = out.get(i).and_then(|r| r.first()).copied().unwrap_or_default(); + let mut held_by = vec![]; + for (j, (o, b)) in out.iter().zip(bounds).enumerate() { + if j == i || !(b.x0 <= p.x && p.x <= b.x1 && b.y0 <= p.y && p.y <= b.y1) { + continue; + } + work = work.saturating_add(o.len()); + if work > MAX_CONTAINMENT_EDGES { + return Err(too_complex()); + } + if contains(o, p) { + held_by.push(j); + } + } + Ok(held_by) + }; + let bounds_of = |out: &[Vec]| -> Vec { + out.iter().map(|r| r.iter().map(|q| Rect::from_points(*q, *q)).reduce(|a, b| a.union(b)).unwrap_or_default()).collect() + }; + if rule == FillRule::NonZero { + let turn = |r: &[Point]| if signed_area(r) > 0.0 { 1 } else { -1 }; + let bounds = bounds_of(&out); + let mut edge = Vec::with_capacity(out.len()); + for (i, r) in out.iter().enumerate() { + // The winding number just outside the ring, and just inside it. + let outside: i32 = within(&out, &bounds, i)?.into_iter().filter_map(|j| out.get(j)).map(|o| turn(o)).sum(); + edge.push((outside != 0) != (outside + turn(r) != 0)); + } + let mut edge = edge.into_iter(); + out.retain(|_| edge.next().unwrap_or(true)); + } + // A ring inside an odd number of others is a hole. + let bounds = bounds_of(&out); + let holes: Vec = (0..out.len()).map(|i| within(&out, &bounds, i).map(|held| held.len() % 2 == 1)).collect::>()?; + for (r, hole) in out.iter_mut().zip(holes) { + if (signed_area(r) > 0.0) == hole { + r.reverse(); + } + } + Ok(out) +} + +/// Most ring edges [`rings`] walks to tell which rings hold which: every ring is tested against +/// the others, so many rings of many points each would otherwise take billions of steps. +const MAX_CONTAINMENT_EDGES: usize = 50_000_000; + +pub(super) fn signed_area(r: &[Point]) -> f64 { + r.iter().zip(r.iter().cycle().skip(1)).take(r.len()).map(|(a, b)| a.x * b.y - b.x * a.y).sum::() / 2.0 +} + +pub(super) fn contains(r: &[Point], p: Point) -> bool { + let mut inside = false; + for (a, b) in r.iter().zip(r.iter().cycle().skip(1)).take(r.len()) { + if (a.y > p.y) != (b.y > p.y) && p.x < a.x + (p.y - a.y) / (b.y - a.y) * (b.x - a.x) { + inside = !inside; + } + } + inside +} + +/// The unit normal on the left of `a`→`b`: towards the filled side. +pub(super) fn left(a: Point, b: Point) -> (f64, f64) { + let (dx, dy) = (b.x - a.x, b.y - a.y); + let l = dx.hypot(dy); + if l < 1e-12 { (0.0, 0.0) } else { (-dy / l, dx / l) } +} + +/// `ring` moved `w` points towards its filled side (negative: away from it), each corner along +/// its bisector, its mitre limited to four times `w` so sharp corners don't shoot off. Each point +/// stays the partner of its source point, so rings at different insets join face to face. +pub(super) fn inset(ring: &[Point], w: f64) -> Vec { + if w == 0.0 { + return ring.to_vec(); + } + let n = ring.len(); + (0..n) + .map(|i| { + let p = ring[i]; + let prev = ring[(i + n - 1) % n]; + let next = ring[(i + 1) % n]; + let (ax, ay) = left(prev, p); + let (bx, by) = left(p, next); + let (mut mx, mut my) = (ax + bx, ay + by); + let ml = mx.hypot(my); + if ml < 1e-9 { + (mx, my) = (ax, ay); + } else { + (mx, my) = (mx / ml, my / ml); + } + let cos = (mx * ax + my * ay).max(0.25); + Point::new(p.x + mx * w / cos, p.y + my * w / cos) + }) + .collect() +}