diff --git a/_config.yml b/_config.yml
index d6b2c72..73f5e6d 100644
--- a/_config.yml
+++ b/_config.yml
@@ -41,3 +41,7 @@ plugins:
# - vendor/cache/
# - vendor/gems/
# - vendor/ruby/
+
+exclude:
+ - Cargo.toml
+ - Cargo.lock
diff --git a/_figures/.gitignore b/_figures/.gitignore
new file mode 100644
index 0000000..4f96631
--- /dev/null
+++ b/_figures/.gitignore
@@ -0,0 +1,2 @@
+/target
+/dist
diff --git a/_figures/beztoy/.gitignore b/_figures/beztoy/.gitignore
new file mode 100644
index 0000000..4f96631
--- /dev/null
+++ b/_figures/beztoy/.gitignore
@@ -0,0 +1,2 @@
+/target
+/dist
diff --git a/_figures/beztoy/Cargo.toml b/_figures/beztoy/Cargo.toml
new file mode 100644
index 0000000..93aacef
--- /dev/null
+++ b/_figures/beztoy/Cargo.toml
@@ -0,0 +1,12 @@
+[package]
+name = "beztoy"
+license = "Apache-2.0"
+version = "0.1.0"
+edition = "2021"
+
+[dependencies]
+console_error_panic_hook = "0.1"
+wasm-bindgen = "0.2.87"
+web-sys = "0.3.64"
+#xilem_svg = { git = "https://github.com/linebender/xilem", rev = "71d1db04dce2dce2264817177575067bec803932"}
+xilem_web = { git = "https://github.com/Philipp-M/xilem", rev = "01838144540210ea14a7b337584c2dd7ff7cf5a3" }
diff --git a/_figures/beztoy/README.md b/_figures/beztoy/README.md
new file mode 100644
index 0000000..f54a358
--- /dev/null
+++ b/_figures/beztoy/README.md
@@ -0,0 +1,5 @@
+# beztoy
+
+A toy for experimenting with Bézier curves, intended to become a testbed for Euler spiral based stroke expansion.
+
+Run with `trunk serve`, then navigate to the webpage.
diff --git a/_figures/beztoy/Trunk.toml b/_figures/beztoy/Trunk.toml
new file mode 100644
index 0000000..af4efa9
--- /dev/null
+++ b/_figures/beztoy/Trunk.toml
@@ -0,0 +1,4 @@
+[build]
+#public_url = "https://levien.com/tmp/beztoy/"
+
+
diff --git a/_figures/beztoy/beztoy.tar.gz b/_figures/beztoy/beztoy.tar.gz
new file mode 100644
index 0000000..7c1fa8d
Binary files /dev/null and b/_figures/beztoy/beztoy.tar.gz differ
diff --git a/_figures/beztoy/index.html b/_figures/beztoy/index.html
new file mode 100644
index 0000000..a73a2ff
--- /dev/null
+++ b/_figures/beztoy/index.html
@@ -0,0 +1,14 @@
+
+
+
+
+
+ See raphlinus#100 for more context.
+
diff --git a/_figures/beztoy/src/euler.rs b/_figures/beztoy/src/euler.rs
new file mode 100644
index 0000000..723bf9c
--- /dev/null
+++ b/_figures/beztoy/src/euler.rs
@@ -0,0 +1,339 @@
+// Copyright 2023 the raphlinus.github.io Authors
+// SPDX-License-Identifier: Apache-2.0
+
+//! Calculations and utilities for Euler spirals
+
+use xilem_web::svg::kurbo::{CubicBez, Point, Vec2};
+
+#[derive(Debug)]
+pub struct CubicParams {
+ pub th0: f64,
+ pub th1: f64,
+ pub d0: f64,
+ pub d1: f64,
+}
+
+#[derive(Debug)]
+pub struct EulerParams {
+ pub th0: f64,
+ pub th1: f64,
+ pub k0: f64,
+ pub k1: f64,
+ pub ch: f64,
+}
+
+#[derive(Debug)]
+pub struct EulerSeg {
+ pub p0: Point,
+ pub p1: Point,
+ pub params: EulerParams,
+}
+
+pub struct CubicToEulerIter {
+ c: CubicBez,
+ tolerance: f64,
+ // [t0 * dt .. (t0 + 1) * dt] is the range we're currently considering
+ t0: u64,
+ dt: f64,
+ last_p: Vec2,
+ last_q: Vec2,
+ last_t: f64,
+}
+
+impl CubicParams {
+ /// Compute parameters from endpoints and derivatives.
+ pub fn from_points_derivs(p0: Vec2, p1: Vec2, q0: Vec2, q1: Vec2, dt: f64) -> Self {
+ let chord = p1 - p0;
+ // Robustness note: we must protect this function from being called when the
+ // chord length is (near-)zero.
+ let scale = dt / chord.length_squared();
+ let h0 = Vec2::new(
+ q0.x * chord.x + q0.y * chord.y,
+ q0.y * chord.x - q0.x * chord.y,
+ );
+ let th0 = h0.atan2();
+ let d0 = h0.length() * scale;
+ let h1 = Vec2::new(
+ q1.x * chord.x + q1.y * chord.y,
+ q1.x * chord.y - q1.y * chord.x,
+ );
+ let th1 = h1.atan2();
+ let d1 = h1.length() * scale;
+ // Robustness note: we may want to clamp the magnitude of the angles to
+ // a bit less than pi. Perhaps here, perhaps downstream.
+ CubicParams { th0, th1, d0, d1 }
+ }
+
+ pub fn from_cubic(c: CubicBez) -> Self {
+ let chord = c.p3 - c.p0;
+ // TODO: if chord is 0, we have a problem
+ let d01 = c.p1 - c.p0;
+ let h0 = Vec2::new(
+ d01.x * chord.x + d01.y * chord.y,
+ d01.y * chord.x - d01.x * chord.y,
+ );
+ let th0 = h0.atan2();
+ let d0 = h0.hypot() / chord.hypot2();
+ let d23 = c.p3 - c.p2;
+ let h1 = Vec2::new(
+ d23.x * chord.x + d23.y * chord.y,
+ d23.x * chord.y - d23.y * chord.x,
+ );
+ let th1 = h1.atan2();
+ let d1 = h1.hypot() / chord.hypot2();
+ CubicParams { th0, th1, d0, d1 }
+ }
+
+ // Estimated error of GH to Euler spiral
+ //
+ // Return value is normalized to chord - to get actual error, multiply
+ // by chord.
+ pub fn est_euler_err(&self) -> f64 {
+ // Potential optimization: work with unit vector rather than angle
+ let cth0 = self.th0.cos();
+ let cth1 = self.th1.cos();
+ if cth0 * cth1 < 0.0 {
+ // Rationale: this happens when fitting a cusp or near-cusp with
+ // a near 180 degree u-turn. The actual ES is bounded in that case.
+ // Further subdivision won't reduce the angles if actually a cusp.
+ return 2.0;
+ }
+ let e0 = (2. / 3.) / (1.0 + cth0);
+ let e1 = (2. / 3.) / (1.0 + cth1);
+ let s0 = self.th0.sin();
+ let s1 = self.th1.sin();
+ // Note: some other versions take sin of s0 + s1 instead. Those are incorrect.
+ // Strangely, calibration is the same, but more work could be done.
+ let s01 = cth0 * s1 + cth1 * s0;
+ let amin = 0.15 * (2. * e0 * s0 + 2. * e1 * s1 - e0 * e1 * s01);
+ let a = 0.15 * (2. * self.d0 * s0 + 2. * self.d1 * s1 - self.d0 * self.d1 * s01);
+ let aerr = (a - amin).abs();
+ let symm = (self.th0 + self.th1).abs();
+ let asymm = (self.th0 - self.th1).abs();
+ let dist = (self.d0 - e0).hypot(self.d1 - e1);
+ let ctr = 3.7e-6 * symm.powi(5) + 6e-3 * asymm * symm.powi(2);
+ let halo_symm = 5e-3 * symm * dist;
+ let halo_asymm = 7e-2 * asymm * dist;
+ 1.25 * ctr + 1.55 * aerr + halo_symm + halo_asymm
+ }
+}
+
+impl EulerParams {
+ pub fn from_angles(th0: f64, th1: f64) -> EulerParams {
+ let k0 = th0 + th1;
+ let dth = th1 - th0;
+ let d2 = dth * dth;
+ let k2 = k0 * k0;
+ let mut a = 6.0;
+ a -= d2 * (1. / 70.);
+ a -= (d2 * d2) * (1. / 10780.);
+ a += (d2 * d2 * d2) * 2.769178184818219e-07;
+ let b = -0.1 + d2 * (1. / 4200.) + d2 * d2 * 1.6959677820260655e-05;
+ let c = -1. / 1400. + d2 * 6.84915970574303e-05 - k2 * 7.936475029053326e-06;
+ a += (b + c * k2) * k2;
+ let k1 = dth * a;
+
+ // calculation of chord
+ let mut ch = 1.0;
+ ch -= d2 * (1. / 40.);
+ ch += (d2 * d2) * 0.00034226190482569864;
+ ch -= (d2 * d2 * d2) * 1.9349474568904524e-06;
+ let b = -1. / 24. + d2 * 0.0024702380951963226 - d2 * d2 * 3.7297408997537985e-05;
+ let c = 1. / 1920. - d2 * 4.87350869747975e-05 - k2 * 3.1001936068463107e-06;
+ ch += (b + c * k2) * k2;
+ EulerParams {
+ th0,
+ th1,
+ k0,
+ k1,
+ ch,
+ }
+ }
+
+ pub fn eval_th(&self, t: f64) -> f64 {
+ (self.k0 + 0.5 * self.k1 * (t - 1.0)) * t - self.th0
+ }
+
+ /// Evaluate the curve at the given parameter.
+ ///
+ /// The parameter is in the range 0..1, and the result goes from (0, 0) to (1, 0).
+ fn eval(&self, t: f64) -> Point {
+ let thm = self.eval_th(t * 0.5);
+ let k0 = self.k0;
+ let k1 = self.k1;
+ let (u, v) = integ_euler_10((k0 + k1 * (0.5 * t - 0.5)) * t, k1 * t * t);
+ let s = t / self.ch * thm.sin();
+ let c = t / self.ch * thm.cos();
+ let x = u * c - v * s;
+ let y = -v * c - u * s;
+ Point::new(x, y)
+ }
+
+ fn eval_with_offset(&self, t: f64, offset: f64) -> Point {
+ let th = self.eval_th(t);
+ let v = Vec2::new(offset * th.sin(), offset * th.cos());
+ self.eval(t) + v
+ }
+
+ // Determine whether a render as a single cubic will be adequate
+ pub fn cubic_ok(&self) -> bool {
+ self.th0.abs() < 1.0 && self.th1.abs() < 1.0
+ }
+}
+
+impl EulerSeg {
+ pub fn from_params(p0: Point, p1: Point, params: EulerParams) -> Self {
+ EulerSeg { p0, p1, params }
+ }
+
+ /// Use two-parabola approximation.
+ pub fn to_cubic(&self) -> CubicBez {
+ let (s0, c0) = self.params.th0.sin_cos();
+ let (s1, c1) = self.params.th1.sin_cos();
+ let d0 = (2. / 3.) / (1.0 + c0);
+ let d1 = (2. / 3.) / (1.0 + c1);
+ let chord = self.p1 - self.p0;
+ let p1 = self.p0 + d0 * Vec2::new(chord.x * c0 - chord.y * s0, chord.y * c0 + chord.x * s0);
+ let p2 = self.p1 - d1 * Vec2::new(chord.x * c1 + chord.y * s1, chord.y * c1 - chord.x * s1);
+ CubicBez::new(self.p0, p1, p2, self.p1)
+ }
+
+ #[allow(unused)]
+ pub fn eval(&self, t: f64) -> Point {
+ let Point { x, y } = self.params.eval(t);
+ let chord = self.p1 - self.p0;
+ Point::new(
+ self.p0.x + chord.x * x - chord.y * y,
+ self.p0.y + chord.x * y + chord.y * x,
+ )
+ }
+
+ pub fn eval_with_offset(&self, t: f64, offset: f64) -> Point {
+ let chord = self.p1 - self.p0;
+ let scaled = offset / chord.hypot();
+ let Point { x, y } = self.params.eval_with_offset(t, scaled);
+ Point::new(
+ self.p0.x + chord.x * x - chord.y * y,
+ self.p0.y + chord.x * y + chord.y * x,
+ )
+ }
+}
+
+/// Evaluate both the point and derivative of a cubic bezier.
+fn eval_cubic_and_deriv(c: &CubicBez, t: f64) -> (Vec2, Vec2) {
+ let p0 = c.p0.to_vec2();
+ let p1 = c.p1.to_vec2();
+ let p2 = c.p2.to_vec2();
+ let p3 = c.p3.to_vec2();
+ let m = 1.0 - t;
+ let mm = m * m;
+ let mt = m * t;
+ let tt = t * t;
+ let p = p0 * (mm * m) + (p1 * (3.0 * mm) + p2 * (3.0 * mt) + p3 * tt) * t;
+ let q = (p1 - p0) * mm + (p2 - p1) * (2.0 * mt) + (p3 - p2) * tt;
+ (p, q)
+}
+
+impl Iterator for CubicToEulerIter {
+ type Item = EulerSeg;
+
+ fn next(&mut self) -> Option {
+ let t0 = (self.t0 as f64) * self.dt;
+ if t0 == 1.0 {
+ return None;
+ }
+ loop {
+ let mut t1 = t0 + self.dt;
+ let p0 = self.last_p;
+ let q0 = self.last_q;
+ let (mut p1, mut q1) = eval_cubic_and_deriv(&self.c, t1);
+ if q1.length_squared() < DERIV_THRESH.powi(2) {
+ let (new_p1, new_q1) = eval_cubic_and_deriv(&self.c, t1 - DERIV_EPS);
+ q1 = new_q1;
+ if t1 < 1. {
+ p1 = new_p1;
+ t1 -= DERIV_EPS;
+ }
+ }
+ // TODO: robustness
+ let actual_dt = t1 - self.last_t;
+ let cubic_params = CubicParams::from_points_derivs(p0, p1, q0, q1, actual_dt);
+ let est_err: f64 = cubic_params.est_euler_err();
+ let err = est_err * (p0 - p1).hypot();
+ if err <= self.tolerance {
+ self.t0 += 1;
+ let shift = self.t0.trailing_zeros();
+ self.t0 >>= shift;
+ self.dt *= (1 << shift) as f64;
+ let euler_params = EulerParams::from_angles(cubic_params.th0, cubic_params.th1);
+ let es = EulerSeg::from_params(p0.to_point(), p1.to_point(), euler_params);
+ self.last_p = p1;
+ self.last_q = q1;
+ self.last_t = t1;
+ return Some(es);
+ }
+ self.t0 *= 2;
+ self.dt *= 0.5;
+ }
+ }
+}
+
+/// Threshold below which a derivative is considered too small.
+const DERIV_THRESH: f64 = 1e-6;
+/// Amount to nudge t when derivative is near-zero.
+const DERIV_EPS: f64 = 1e-6;
+
+impl CubicToEulerIter {
+ pub fn new(c: CubicBez, tolerance: f64) -> Self {
+ let mut last_q = c.p1 - c.p0;
+ // TODO: tweak
+ if last_q.length_squared() < DERIV_THRESH.powi(2) {
+ last_q = eval_cubic_and_deriv(&c, DERIV_EPS).1;
+ }
+ CubicToEulerIter {
+ c,
+ tolerance,
+ t0: 0,
+ dt: 1.0,
+ last_p: c.p0.to_vec2(),
+ last_q,
+ last_t: 0.0,
+ }
+ }
+}
+
+/// Integrate Euler spiral.
+///
+/// TODO: investigate needed accuracy. We might be able to get away
+/// with 8th order.
+fn integ_euler_10(k0: f64, k1: f64) -> (f64, f64) {
+ let t1_1 = k0;
+ let t1_2 = 0.5 * k1;
+ let t2_2 = t1_1 * t1_1;
+ let t2_3 = 2. * (t1_1 * t1_2);
+ let t2_4 = t1_2 * t1_2;
+ let t3_4 = t2_2 * t1_2 + t2_3 * t1_1;
+ let t3_6 = t2_4 * t1_2;
+ let t4_4 = t2_2 * t2_2;
+ let t4_5 = 2. * (t2_2 * t2_3);
+ let t4_6 = 2. * (t2_2 * t2_4) + t2_3 * t2_3;
+ let t4_7 = 2. * (t2_3 * t2_4);
+ let t4_8 = t2_4 * t2_4;
+ let t5_6 = t4_4 * t1_2 + t4_5 * t1_1;
+ let t5_8 = t4_6 * t1_2 + t4_7 * t1_1;
+ let t6_6 = t4_4 * t2_2;
+ let t6_7 = t4_4 * t2_3 + t4_5 * t2_2;
+ let t6_8 = t4_4 * t2_4 + t4_5 * t2_3 + t4_6 * t2_2;
+ let t7_8 = t6_6 * t1_2 + t6_7 * t1_1;
+ let t8_8 = t6_6 * t2_2;
+ let mut u = 1.;
+ u -= (1. / 24.) * t2_2 + (1. / 160.) * t2_4;
+ u += (1. / 1920.) * t4_4 + (1. / 10752.) * t4_6 + (1. / 55296.) * t4_8;
+ u -= (1. / 322560.) * t6_6 + (1. / 1658880.) * t6_8;
+ u += (1. / 92897280.) * t8_8;
+ let mut v = (1. / 12.) * t1_2;
+ v -= (1. / 480.) * t3_4 + (1. / 2688.) * t3_6;
+ v += (1. / 53760.) * t5_6 + (1. / 276480.) * t5_8;
+ v -= (1. / 11612160.) * t7_8;
+ (u, v)
+}
diff --git a/_figures/beztoy/src/flatten.rs b/_figures/beztoy/src/flatten.rs
new file mode 100644
index 0000000..13278b6
--- /dev/null
+++ b/_figures/beztoy/src/flatten.rs
@@ -0,0 +1,91 @@
+// Copyright 2023 the raphlinus.github.io Authors
+// SPDX-License-Identifier: Apache-2.0
+
+//! Math for flattening of Euler spiral parallel curve
+
+use std::f64::consts::FRAC_PI_4;
+
+use xilem_web::svg::kurbo::BezPath;
+
+use crate::euler::EulerSeg;
+
+pub fn flatten_offset(iter: impl Iterator- , offset: f64) -> BezPath {
+ let mut result = BezPath::new();
+ let mut first = true;
+ for es in iter {
+ if core::mem::take(&mut first) {
+ result.move_to(es.eval_with_offset(0.0, offset));
+ }
+ let scale = es.p0.distance(es.p1);
+ let tol = 1.0;
+ let (k0, k1) = (es.params.k0 - 0.5 * es.params.k1, es.params.k1);
+ // compute forward integral to determine number of subdivisions
+ let dist_scaled = offset / scale;
+ let a = -2.0 * dist_scaled * k1;
+ let b = -1.0 - 2.0 * dist_scaled * k0;
+ let int0 = espc_int_approx(b);
+ let int1 = espc_int_approx(a + b);
+ let integral = int1 - int0;
+ let k_peak = k0 - k1 * b / a;
+ let integrand_peak = (k_peak * (k_peak * dist_scaled + 1.0)).abs().sqrt();
+ let scaled_int = integral * integrand_peak / a;
+ let n_frac = 0.5 * (scale / tol).sqrt() * scaled_int;
+ let n = n_frac.ceil();
+ for i in 0..n as usize {
+ let t = (i + 1) as f64 / n;
+ let inv = espc_int_inv_approx(integral * t + int0);
+ let s = (inv - b) / a;
+ result.line_to(es.eval_with_offset(s, offset));
+ }
+ }
+ result
+}
+
+const BREAK1: f64 = 0.8;
+const BREAK2: f64 = 1.25;
+const BREAK3: f64 = 2.1;
+const SIN_SCALE: f64 = 1.0976991822760038;
+const QUAD_A1: f64 = 0.6406;
+const QUAD_B1: f64 = -0.81;
+const QUAD_C1: f64 = 0.9148117935952064;
+const QUAD_A2: f64 = 0.5;
+const QUAD_B2: f64 = -0.156;
+const QUAD_C2: f64 = 0.16145779359520596;
+
+fn espc_int_approx(x: f64) -> f64 {
+ let y = x.abs();
+ let a = if y < BREAK1 {
+ (SIN_SCALE * y).sin() * (1.0 / SIN_SCALE)
+ } else if y < BREAK2 {
+ (8.0f64.sqrt() / 3.0) * (y - 1.0) * (y - 1.0).abs().sqrt() + FRAC_PI_4
+ } else {
+ let (a, b, c) = if y < BREAK3 {
+ (QUAD_A1, QUAD_B1, QUAD_C1)
+ } else {
+ (QUAD_A2, QUAD_B2, QUAD_C2)
+ };
+ a * y * y + b * y + c
+ };
+ a.copysign(x)
+}
+
+fn espc_int_inv_approx(x: f64) -> f64 {
+ let y = x.abs();
+ let a = if y < 0.7010707591262915 {
+ (x * SIN_SCALE).asin() * (1.0 / SIN_SCALE)
+ } else if y < 0.903249293595206 {
+ let b = y - FRAC_PI_4;
+ let u = b.abs().powf(2. / 3.).copysign(b);
+ u * (9.0f64 / 8.).cbrt() + 1.0
+ } else {
+ let (u, v, w) = if y < 2.038857793595206 {
+ const B: f64 = 0.5 * QUAD_B1 / QUAD_A1;
+ (B * B - QUAD_C1 / QUAD_A1, 1.0 / QUAD_A1, B)
+ } else {
+ const B: f64 = 0.5 * QUAD_B2 / QUAD_A2;
+ (B * B - QUAD_C2 / QUAD_A2, 1.0 / QUAD_A2, B)
+ };
+ (u + v * y).sqrt() - w
+ };
+ a.copysign(x)
+}
diff --git a/_figures/beztoy/src/main.rs b/_figures/beztoy/src/main.rs
new file mode 100644
index 0000000..d24ae65
--- /dev/null
+++ b/_figures/beztoy/src/main.rs
@@ -0,0 +1,158 @@
+// Copyright 2023 the raphlinus.github.io Authors
+// SPDX-License-Identifier: Apache-2.0
+
+//! An interactive toy for experimenting with rendering of Bézier paths,
+//! including Euler spiral based stroke expansion.
+
+mod euler;
+mod flatten;
+
+use xilem_web::{svg::{kurbo::{Point, BezPath, CubicBez, PathEl, Circle, Line, Shape}, peniko::Color}, PointerMsg, View, App, elements::svg::{g, svg}, document_body, interfaces::*};
+
+use crate::{
+ euler::{CubicParams, CubicToEulerIter},
+ flatten::flatten_offset,
+};
+
+#[derive(Default)]
+struct AppState {
+ p0: Point,
+ p1: Point,
+ p2: Point,
+ p3: Point,
+ grab: GrabState,
+}
+
+#[derive(Default)]
+struct GrabState {
+ is_down: bool,
+ id: i32,
+ dx: f64,
+ dy: f64,
+}
+
+impl GrabState {
+ fn handle(&mut self, pt: &mut Point, p: &PointerMsg) {
+ match p {
+ PointerMsg::Down(e) => {
+ if e.button == 0 {
+ self.dx = pt.x - e.x;
+ self.dy = pt.y - e.y;
+ self.id = e.id;
+ self.is_down = true;
+ }
+ }
+ PointerMsg::Move(e) => {
+ if self.is_down && self.id == e.id {
+ pt.x = self.dx + e.x;
+ pt.y = self.dy + e.y;
+ }
+ }
+ PointerMsg::Up(e) => {
+ if self.id == e.id {
+ self.is_down = false;
+ }
+ }
+ }
+ }
+}
+
+// https://iamkate.com/data/12-bit-rainbow/
+const RAINBOW_PALETTE: [Color; 12] = [
+ Color::rgb8(0x88, 0x11, 0x66),
+ Color::rgb8(0xaa, 0x33, 0x55),
+ Color::rgb8(0xcc, 0x66, 0x66),
+ Color::rgb8(0xee, 0x99, 0x44),
+ Color::rgb8(0xee, 0xdd, 0x00),
+ Color::rgb8(0x99, 0xdd, 0x55),
+ Color::rgb8(0x44, 0xdd, 0x88),
+ Color::rgb8(0x22, 0xcc, 0xbb),
+ Color::rgb8(0x00, 0xbb, 0xcc),
+ Color::rgb8(0x00, 0x99, 0xcc),
+ Color::rgb8(0x33, 0x66, 0xbb),
+ Color::rgb8(0x66, 0x33, 0x99),
+];
+
+fn app_logic(state: &mut AppState) -> impl View {
+ let mut path = BezPath::new();
+ path.move_to(state.p0);
+ path.curve_to(state.p1, state.p2, state.p3);
+ let stroke = xilem_web::svg::kurbo::Stroke::new(2.0);
+ let stroke_thick = xilem_web::svg::kurbo::Stroke::new(8.0);
+ let stroke_thin = xilem_web::svg::kurbo::Stroke::new(1.0);
+ const NONE: Color = Color::TRANSPARENT;
+ const HANDLE_RADIUS: f64 = 4.0;
+ let c = CubicBez::new(state.p0, state.p1, state.p2, state.p3);
+ let params = CubicParams::from_cubic(c);
+ let err = params.est_euler_err();
+ let mut spirals = vec![];
+ const TOL: f64 = 1.0;
+ for (i, es) in CubicToEulerIter::new(c, TOL).enumerate() {
+ let path = if es.params.cubic_ok() {
+ es.to_cubic().into_path(1.0)
+ } else {
+ // Janky rendering, we should be more sophisticated
+ // and subdivide into cubics with appropriate bounds
+ let mut path = BezPath::new();
+ const N: usize = 20;
+ path.move_to(es.p0);
+ for i in 1..N {
+ let t = i as f64 / N as f64;
+ path.line_to(es.eval(t));
+ }
+ path.line_to(es.p1);
+ path
+ };
+ let color = RAINBOW_PALETTE[(i * 7) % 12];
+ spirals.push(path.stroke(color, stroke_thick.clone()).fill(NONE));
+ }
+ let offset = 40.0;
+ let flat = flatten_offset(CubicToEulerIter::new(c, TOL), offset);
+ let flat2 = flatten_offset(CubicToEulerIter::new(c, TOL), -offset);
+ let mut flat_pts = vec![];
+ for seg in flat.elements().iter().chain(flat2.elements().iter()) {
+ match seg {
+ PathEl::MoveTo(p) | PathEl::LineTo(p) => {
+ let circle = Circle::new(*p, 2.0).fill(Color::BLACK);
+ flat_pts.push(circle);
+ }
+ _ => (),
+ }
+ }
+ svg(g((
+ g(spirals),
+ path.stroke(Color::BLACK, stroke_thin.clone()).fill(NONE),
+ flat.stroke(Color::BLUE, stroke_thin.clone()).fill(NONE),
+ flat2.stroke(Color::PURPLE, stroke_thin).fill(NONE),
+ g(flat_pts),
+ Line::new(state.p0, state.p1)
+ .stroke(Color::BLUE, stroke.clone()),
+ Line::new(state.p2, state.p3)
+ .stroke(Color::BLUE, stroke.clone()),
+ Line::new((790., 300.), (790., 300. - 1000. * err))
+ .stroke(Color::RED, stroke.clone()),
+ g((
+ Circle::new(state.p0, HANDLE_RADIUS)
+ .pointer(|s: &mut AppState, msg| s.grab.handle(&mut s.p0, &msg)),
+ Circle::new(state.p1, HANDLE_RADIUS)
+ .pointer(|s: &mut AppState, msg| s.grab.handle(&mut s.p1, &msg)),
+ Circle::new(state.p2, HANDLE_RADIUS)
+ .pointer(|s: &mut AppState, msg| s.grab.handle(&mut s.p2, &msg)),
+ Circle::new(state.p3, HANDLE_RADIUS)
+ .pointer(|s: &mut AppState, msg| s.grab.handle(&mut s.p3, &msg)),
+ )),
+ )))
+ .attr("width", 800)
+ .attr("height", 600)
+}
+
+pub fn main() {
+ console_error_panic_hook::set_once();
+ let mut state = AppState::default();
+ state.p0 = Point::new(100.0, 100.0);
+ state.p1 = Point::new(300.0, 150.0);
+ state.p2 = Point::new(500.0, 150.0);
+ state.p3 = Point::new(700.0, 150.0);
+ let app = App::new(state, app_logic);
+ app.run(&document_body());
+}
diff --git a/_figures/grid_robustness.drawio b/_figures/grid_robustness.drawio
new file mode 100644
index 0000000..9add686
--- /dev/null
+++ b/_figures/grid_robustness.drawio
@@ -0,0 +1 @@
+5VZdb5swFP01PK4CjEn62NJsU6d2raKozxa+AasGR45Tkv76mWADDkmzRam6dXmI7HM/bJ9zdIWHkmL9TZJFficocC/06dpDN14YXmJf/9fApgHiUdwAmWS0gYIOmLJXMKCpy1aMwtJJVEJwxRYumIqyhFQ5GJFSVG7aXHD31AXJYABMU8KH6BOjKm/QsX1WjX8HluX25MA3kYLYZAMsc0JF1YPQxEOJFEI1q2KdAK+5s7w0dV8PRNuLSSjV7xTMHuR88nr/YxYr/Pzz7v7xdnb7xXR5IXxlHmwuqzaWASlWJYW6ie+h6ypnCqYLktbRSkuusVwVXO8CvZwzzhPBhdzWIophTCONL5UUz9CL+NufjpgLgFSwPviyoOVL+wxEAUpudIotiA3FxmPh2OyrTrE2J++rZUFiXJK1vTsi9cJw+Qe8hsd5hZJe1QbVu1KU4PLokt7UAh1Y9ShDPQbwHgIsJoETxV7c9vtIMSc8CKYP7gQIXAEQ3uF1KVYyBVPV9+hOI+0LV8lwp5EiMgM1aLTVqH326bLFZ5VtYPokcUz/pqCHbfVhKuMdcS7xBT5N52B0tNU7K40+w+CLxn/d4Iv+q8EXnWvw4Q8efPh9B5/v16PvHxA0GEV2ELWj6URJ97QanUtUve0+HJv07usbTX4B
\ No newline at end of file
diff --git a/_figures/simplify_figs/.DS_Store b/_figures/simplify_figs/.DS_Store
new file mode 100644
index 0000000..81229e8
Binary files /dev/null and b/_figures/simplify_figs/.DS_Store differ
diff --git a/_figures/xilem_view.drawio b/_figures/xilem_view.drawio
new file mode 100644
index 0000000..e4dc1f3
--- /dev/null
+++ b/_figures/xilem_view.drawio
@@ -0,0 +1 @@
+zVZNc5swEP01HONByHz4mDiue/FMZtxpm6MKG9AEEBXC4P76Ckvi027STDu2L9a+Xa2W994ILLzOmi0nRbJjEaSWY0eNhR8tx0EIL+VfixwV4nm+AmJOI13UA3v6CzRoa7SiEZSjQsFYKmgxBkOW5xCKEUY4Z/W47IWl41MLEsMM2IcknaPfaCQShQau3eOfgcaJORnZOpMRU6yBMiERqwcQ3lh4zRkTapU1a0hb8gwvat+nC9luMA65eM+G7et2x35w3uyeANYlL3Jvc6e7HEha6Qe2HC+V/R5emGwrpxZHTYX3s2ImcVeehLqXBSgomj4pV3H7/3UvSPhqesmhVDuV1Hx0nR3OqjyCdk5bpuuECtgXJGyztbSVxBKRpTJC3e4DcAHNRSpQR7B0JrAMBD/KErNhpTUxptRh3SuMjGzJQF1TR7Sp4q5zz7tcaOr/QgbnogxlQfKPy/AFGjEQQTW7EREmGqDg2iLg/yTCQyUEy29VBie4NR2WZ3SYkAR5dN9e7jIKU1KWNBzzMiaxVUm/WVAgY8kMP35vkwvXhM+69hQ8NqPoqCM1BkSzN8aEazkqq3gIb9+7gvAYxFsXw1y7gTbuGWkMxiElgh7G457TS5/wxOjp0m8mVphaw7RQj6l3DV89k0YOnjTCk0aKh1mjk326x/64o/xbdRQ0VKhtvqvDZ31gu+53tcEVbYivaUMHowVeYWTj5TJwfHfl/tlM73Ul8v2FjVf9z5u43f9HJpVh/5mnyvuPZbz5DQ==
\ No newline at end of file