effectcraft/crates/path/src/interp.rs

596 lines
22 KiB
Rust
Raw Blame History

This file contains ambiguous Unicode characters

This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.

//! Smart mask interpolation (the Mask Interpolation panel): give two mask shapes the same vertex
//! count with a sensible vertex correspondence, and compute in-between shapes.
//!
//! - **Vertex insertion** splits Bezier segments with de Casteljau subdivision, so both shapes
//! keep their exact outline while gaining vertices; vertices are distributed by **arc length**
//! so the i-th vertex of each shape sits at a comparable fraction of the outline.
//! - **Correspondence** picks the cyclic offset (and, for closed shapes, the winding) of the
//! second shape that minimises a matching cost: normalised position differences plus, when
//! Quality is above 0, the χ² distance of log-polar **shape context** histograms (Belongie,
//! Malik & Puzicha, "Shape Matching and Object Recognition Using Shape Contexts", PAMI 2002).
//! *First Vertices Match* pins the offset to 0; *Use 1:1 Vertex Matches* keeps the original
//! vertices when the counts agree.
//! - **In-betweens** move vertices along straight lines (*Use Linear Vertex Paths*) or, by
//! default, with the shapes' best-fit rotation and scale (Procrustes) interpolated so turning
//! shapes turn instead of collapsing; *Bending Resistance* weights that rigid component
//! against the straight-line one.
use effectcraft_keyframe::ShapePath;
/// Matching Method.
#[derive(Clone, Copy, Debug, Default, PartialEq, Eq)]
pub enum Matching {
#[default]
Auto,
Curve,
Polyline,
}
/// Add Mask Shape Vertices.
#[derive(Clone, Copy, Debug, PartialEq)]
pub enum AddVertices {
/// Pixels between vertices.
Pixels(f64),
/// Total vertices.
Total(usize),
/// One vertex every this percentage of the outline.
Percent(f64),
}
/// Mask Interpolation options.
#[derive(Clone, Copy, Debug, PartialEq)]
pub struct InterpOpts {
pub linear: bool,
/// Bending Resistance, 0–100.
pub bending_resistance: f64,
/// Quality, 0–100.
pub quality: f64,
pub add_vertices: Option<AddVertices>,
pub matching: Matching,
pub one_to_one: bool,
pub first_vertices_match: bool,
}
impl Default for InterpOpts {
fn default() -> Self {
InterpOpts {
linear: false,
bending_resistance: 50.0,
quality: 50.0,
add_vertices: Some(AddVertices::Pixels(15.0)),
matching: Matching::Auto,
one_to_one: false,
first_vertices_match: true,
}
}
}
type P = [f64; 2];
fn add(a: P, b: P) -> P {
[a[0] + b[0], a[1] + b[1]]
}
fn sub(a: P, b: P) -> P {
[a[0] - b[0], a[1] - b[1]]
}
fn mul(a: P, k: f64) -> P {
[a[0] * k, a[1] * k]
}
fn lerp(a: P, b: P, t: f64) -> P {
[a[0] + (b[0] - a[0]) * t, a[1] + (b[1] - a[1]) * t]
}
/// Segment `i` (vertex i → i+1, wrapping for closed paths) as cubic control points.
fn segment(p: &ShapePath, i: usize) -> [P; 4] {
let n = p.len();
let j = (i + 1) % n;
let (a, b) = (p.vertices[i], p.vertices[j]);
let o = p.out_tangents.get(i).copied().unwrap_or([0.0; 2]);
let it = p.in_tangents.get(j).copied().unwrap_or([0.0; 2]);
[a, add(a, o), add(b, it), b]
}
fn seg_count(p: &ShapePath) -> usize {
let n = p.len();
if n < 2 {
0
} else if p.closed {
n
} else {
n - 1
}
}
fn cubic(c: &[P; 4], t: f64) -> P {
let u = 1.0 - t;
let (a, b, cc, d) = (u * u * u, 3.0 * u * u * t, 3.0 * u * t * t, t * t * t);
[a * c[0][0] + b * c[1][0] + cc * c[2][0] + d * c[3][0], a * c[0][1] + b * c[1][1] + cc * c[2][1] + d * c[3][1]]
}
fn seg_len(c: &[P; 4]) -> f64 {
let mut l = 0.0;
let mut prev = c[0];
for k in 1..=24 {
let q = cubic(c, k as f64 / 24.0);
l += (q[0] - prev[0]).hypot(q[1] - prev[1]);
prev = q;
}
l
}
/// Outline length.
pub fn length(p: &ShapePath) -> f64 {
(0..seg_count(p)).map(|i| seg_len(&segment(p, i))).sum()
}
/// Split cubic `c` at parameter values `ts` (increasing, in (0, 1)) into `ts.len() + 1` cubics.
fn split_many(c: &[P; 4], ts: &[f64]) -> Vec<[P; 4]> {
let mut out = vec![];
let mut rest = *c;
let mut t0 = 0.0;
for &t in ts {
let u = ((t - t0) / (1.0 - t0)).clamp(0.0, 1.0);
let (l, r) = split(&rest, u);
out.push(l);
rest = r;
t0 = t;
}
out.push(rest);
out
}
fn split(c: &[P; 4], t: f64) -> ([P; 4], [P; 4]) {
let p01 = lerp(c[0], c[1], t);
let p12 = lerp(c[1], c[2], t);
let p23 = lerp(c[2], c[3], t);
let a = lerp(p01, p12, t);
let b = lerp(p12, p23, t);
let m = lerp(a, b, t);
([c[0], p01, a, m], [m, b, p23, c[3]])
}
/// Parameter values splitting `c` into `k` pieces of (approximately) equal arc length.
fn arc_params(c: &[P; 4], k: usize) -> Vec<f64> {
if k <= 1 {
return vec![];
}
const S: usize = 64;
let mut cum = vec![0.0; S + 1];
let mut prev = c[0];
for i in 1..=S {
let q = cubic(c, i as f64 / S as f64);
cum[i] = cum[i - 1] + (q[0] - prev[0]).hypot(q[1] - prev[1]);
prev = q;
}
let total = cum[S];
if total < 1e-12 {
return (1..k).map(|j| j as f64 / k as f64).collect();
}
(1..k)
.map(|j| {
let target = total * j as f64 / k as f64;
let i = cum.partition_point(|v| *v < target).clamp(1, S);
let f = (target - cum[i - 1]) / (cum[i] - cum[i - 1]).max(1e-12);
((i - 1) as f64 + f) / S as f64
})
.collect()
}
/// Rebuild a ShapePath from consecutive cubics (closed: the last ends at the first's start).
fn from_cubics(cs: &[[P; 4]], closed: bool) -> ShapePath {
let n = if closed { cs.len() } else { cs.len() + 1 };
let mut v = vec![[0.0; 2]; n];
let mut ins = vec![[0.0; 2]; n];
let mut outs = vec![[0.0; 2]; n];
for (i, c) in cs.iter().enumerate() {
let j = if closed { (i + 1) % n } else { i + 1 };
v[i] = c[0];
outs[i] = sub(c[1], c[0]);
ins[j] = sub(c[2], c[3]);
if !closed && j == n - 1 {
v[j] = c[3];
}
}
ShapePath { vertices: v, in_tangents: ins, out_tangents: outs, closed, feather: Vec::new() }
}
/// Exactly the same outline with `total` vertices (≥ the current count), the extra vertices
/// distributed over segments in proportion to their length.
pub fn with_vertex_count(p: &ShapePath, total: usize) -> ShapePath {
let m = seg_count(p);
let base = if p.closed { m } else { m + 1 };
if m == 0 || total <= base {
return p.clone();
}
let segs: Vec<[P; 4]> = (0..m).map(|i| segment(p, i)).collect();
let lens: Vec<f64> = segs.iter().map(seg_len).collect();
let sum: f64 = lens.iter().sum::<f64>().max(1e-12);
let pieces = total - base + m; // total segments wanted
// Largest remainder allocation, at least one piece per segment.
let extra = pieces - m;
let want: Vec<f64> = lens.iter().map(|l| l / sum * extra as f64).collect();
let mut k: Vec<usize> = want.iter().map(|w| 1 + w.floor() as usize).collect();
let mut left = pieces - k.iter().sum::<usize>();
let mut order: Vec<usize> = (0..m).collect();
order.sort_by(|a, b| (want[*b] - want[*b].floor()).total_cmp(&(want[*a] - want[*a].floor())).then(a.cmp(b)));
for i in order.iter().cycle() {
if left == 0 {
break;
}
k[*i] += 1;
left -= 1;
}
let mut cubics = vec![];
for (c, ki) in segs.iter().zip(&k) {
cubics.extend(split_many(c, &arc_params(c, *ki)));
}
from_cubics(&cubics, p.closed)
}
/// The same outline as a polyline with `total` vertices at equal arc-length spacing, starting at
/// vertex 0.
pub fn polyline(p: &ShapePath, total: usize) -> ShapePath {
let m = seg_count(p);
if m == 0 || total < 2 {
return ShapePath::polygon(&p.vertices, p.closed);
}
// Dense samples with cumulative length.
let mut pts = vec![p.vertices[0]];
for i in 0..m {
let c = segment(p, i);
for k in 1..=32 {
pts.push(cubic(&c, k as f64 / 32.0));
}
}
let mut cum = vec![0.0];
for w in pts.windows(2) {
let l = cum.last().copied().unwrap_or(0.0) + (w[1][0] - w[0][0]).hypot(w[1][1] - w[0][1]);
cum.push(l);
}
let total_len = *cum.last().unwrap_or(&0.0);
let n = total;
let div = if p.closed { n as f64 } else { (n - 1) as f64 };
let out: Vec<P> = (0..n)
.map(|j| {
let target = total_len * j as f64 / div;
let i = cum.partition_point(|v| *v < target).clamp(1, cum.len() - 1);
let f = (target - cum[i - 1]) / (cum[i] - cum[i - 1]).max(1e-12);
lerp(pts[i - 1], pts[i], f)
})
.collect();
ShapePath::polygon(&out, p.closed)
}
fn has_curves(p: &ShapePath) -> bool {
p.in_tangents.iter().chain(&p.out_tangents).any(|t| t[0].abs() > 1e-9 || t[1].abs() > 1e-9)
}
/// Rotate a closed path's vertex order to start at `k`.
fn rotate(p: &ShapePath, k: usize) -> ShapePath {
let n = p.len();
if n == 0 || k.is_multiple_of(n) {
return p.clone();
}
let r = |v: &Vec<P>| (0..n).map(|i| v[(i + k) % n]).collect();
ShapePath {
vertices: r(&p.vertices),
in_tangents: r(&p.in_tangents),
out_tangents: r(&p.out_tangents),
closed: p.closed,
feather: p.feather.iter().map(|f| effectcraft_keyframe::FeatherPoint { segment: (f.segment + n - k % n) % n, ..*f }).collect(),
}
}
/// The same path traversed the other way (closed paths keep vertex 0 first).
pub fn reversed(p: &ShapePath) -> ShapePath {
let n = p.len();
if n == 0 {
return p.clone();
}
let idx: Vec<usize> = if p.closed { (0..n).map(|i| (n - i) % n).collect() } else { (0..n).rev().collect() };
ShapePath {
vertices: idx.iter().map(|i| p.vertices[*i]).collect(),
in_tangents: idx.iter().map(|i| p.out_tangents[*i]).collect(),
out_tangents: idx.iter().map(|i| p.in_tangents[*i]).collect(),
closed: p.closed,
// Feather points sit on segments that reversal renumbers: they are not carried.
feather: Vec::new(),
}
}
fn signed_area(v: &[P]) -> f64 {
let n = v.len();
(0..n).map(|i| v[i][0] * v[(i + 1) % n][1] - v[(i + 1) % n][0] * v[i][1]).sum::<f64>() * 0.5
}
/// Centroid and RMS radius of points.
fn frame(v: &[P]) -> (P, f64) {
let n = v.len().max(1) as f64;
let c = [v.iter().map(|p| p[0]).sum::<f64>() / n, v.iter().map(|p| p[1]).sum::<f64>() / n];
let r = (v.iter().map(|p| (p[0] - c[0]).powi(2) + (p[1] - c[1]).powi(2)).sum::<f64>() / n).sqrt().max(1e-9);
(c, r)
}
const RBINS: usize = 5;
const ABINS: usize = 12;
/// Log-polar shape context histograms of every point (normalised by the RMS radius).
fn shape_contexts(v: &[P]) -> Vec<[f64; RBINS * ABINS]> {
let (_, r) = frame(v);
let n = v.len();
(0..n)
.map(|i| {
let mut h = [0.0; RBINS * ABINS];
for j in 0..n {
if i == j {
continue;
}
let d = sub(v[j], v[i]);
let dist = d[0].hypot(d[1]) / r;
// log radius in [1/8, 2] → 5 bins.
let lr = (dist.max(1e-6).log2() + 3.0) / 4.0 * RBINS as f64;
if !(0.0..RBINS as f64).contains(&lr) {
continue;
}
let a = (d[1].atan2(d[0]) + std::f64::consts::PI) / std::f64::consts::TAU * ABINS as f64;
let ai = (a as usize).min(ABINS - 1);
h[lr as usize * ABINS + ai] += 1.0;
}
let s: f64 = h.iter().sum::<f64>().max(1.0);
h.iter_mut().for_each(|x| *x /= s);
h
})
.collect()
}
fn chi2(a: &[f64], b: &[f64]) -> f64 {
a.iter().zip(b).map(|(x, y)| if x + y > 1e-12 { (x - y) * (x - y) / (x + y) } else { 0.0 }).sum::<f64>() * 0.5
}
/// Matching cost of `a[i] ↔ b[(i + k) % n]`.
fn cost(a: &[P], b: &[P], sa: Option<&[[f64; RBINS * ABINS]]>, sb: Option<&[[f64; RBINS * ABINS]]>, k: usize) -> f64 {
let n = a.len();
let (ca, ra) = frame(a);
let (cb, rb) = frame(b);
let mut c = 0.0;
for i in 0..n {
let j = (i + k) % n;
let pa = mul(sub(a[i], ca), 1.0 / ra);
let pb = mul(sub(b[j], cb), 1.0 / rb);
c += (pa[0] - pb[0]).powi(2) + (pa[1] - pb[1]).powi(2);
if let (Some(sa), Some(sb)) = (sa, sb) {
c += chi2(&sa[i], &sb[j]);
}
}
c / n as f64
}
/// Give `a` and `b` the same vertex count with matching vertices (the shapes that the first
/// and last keys get, and the input of [`interpolate`]).
pub fn correspond(a: &ShapePath, b: &ShapePath, o: &InterpOpts) -> (ShapePath, ShapePath) {
if a.is_empty() || b.is_empty() {
return (a.clone(), b.clone());
}
let closed = a.closed && b.closed;
let curves = match o.matching {
Matching::Curve => true,
Matching::Polyline => false,
Matching::Auto => has_curves(a) || has_curves(b),
};
let (mut ra, mut rb);
if o.one_to_one && a.len() == b.len() {
(ra, rb) = (a.clone(), b.clone());
} else {
let (la, lb) = (length(a), length(b));
let mut n = a.len().max(b.len());
if let Some(av) = o.add_vertices {
let want = match av {
AddVertices::Pixels(px) => (la.max(lb) / px.max(1.0)).ceil() as usize,
AddVertices::Total(t) => t,
AddVertices::Percent(pc) => (100.0 / pc.clamp(0.1, 100.0)).ceil() as usize,
};
n = n.max(want.min(2000));
}
if curves {
ra = with_vertex_count(a, n);
rb = with_vertex_count(b, n);
} else {
ra = polyline(a, n);
rb = polyline(b, n);
}
// Proportional splitting can leave counts unequal when one shape already had more
// vertices than its share: fall back to arc-length polylines then.
if ra.len() != rb.len() {
ra = polyline(a, n);
rb = polyline(b, n);
}
}
if closed && ra.len() == rb.len() && ra.len() >= 3 {
// Same winding.
if signed_area(&ra.vertices).signum() != signed_area(&rb.vertices).signum() {
rb = reversed(&rb);
}
if !o.first_vertices_match {
let n = ra.len();
let use_sc = o.quality > 0.0;
let (sa, sb) = if use_sc { (Some(shape_contexts(&ra.vertices)), Some(shape_contexts(&rb.vertices))) } else { (None, None) };
// Quality limits the offsets examined exhaustively (coarse step), then refines.
let step = if o.quality >= 50.0 { 1 } else { ((n as f64 / 24.0).ceil() as usize).max(1) };
let mut best = (f64::MAX, 0usize);
let mut k = 0;
while k < n {
let c = cost(&ra.vertices, &rb.vertices, sa.as_deref(), sb.as_deref(), k);
if c < best.0 {
best = (c, k);
}
k += step;
}
if step > 1 {
let k0 = best.1;
for d in 1..step {
for k in [(k0 + d) % n, (k0 + n - d) % n] {
let c = cost(&ra.vertices, &rb.vertices, sa.as_deref(), sb.as_deref(), k);
if c < best.0 {
best = (c, k);
}
}
}
}
rb = rotate(&rb, best.1);
}
}
(ra, rb)
}
/// The in-between shape at fraction `u` of two corresponded shapes (equal vertex counts).
pub fn interpolate(a: &ShapePath, b: &ShapePath, u: f64, o: &InterpOpts) -> ShapePath {
if a.len() != b.len() || a.is_empty() {
return if u < 1.0 { a.clone() } else { b.clone() };
}
let n = a.len();
let lin = |pa: &[P], pb: &[P]| -> Vec<P> { pa.iter().zip(pb).map(|(x, y)| lerp(*x, *y, u)).collect() };
if o.linear {
return ShapePath {
vertices: lin(&a.vertices, &b.vertices),
in_tangents: lin(&a.in_tangents, &b.in_tangents),
out_tangents: lin(&a.out_tangents, &b.out_tangents),
closed: a.closed,
feather: a.feather.clone(),
};
}
// Procrustes: b ≈ s R (a − ca) + cb.
let (ca, _) = frame(&a.vertices);
let (cb, _) = frame(&b.vertices);
let (mut ss, mut dot, mut cross) = (0.0, 0.0, 0.0);
for i in 0..n {
let (p, q) = (sub(a.vertices[i], ca), sub(b.vertices[i], cb));
ss += p[0] * p[0] + p[1] * p[1];
dot += p[0] * q[0] + p[1] * q[1];
cross += p[0] * q[1] - p[1] * q[0];
}
let theta = cross.atan2(dot);
let s = if ss > 1e-12 { dot.hypot(cross) / ss } else { 1.0 };
let rot = |v: P, ang: f64, k: f64| {
let (sn, cs) = ang.sin_cos();
[k * (cs * v[0] - sn * v[1]), k * (sn * v[0] + cs * v[1])]
};
let (th_u, s_u) = (theta * u, s.max(1e-9).powf(u));
let c_u = lerp(ca, cb, u);
let w = (o.bending_resistance / 100.0).clamp(0.0, 1.0);
let mut out = ShapePath { vertices: vec![], in_tangents: vec![], out_tangents: vec![], closed: a.closed, feather: a.feather.clone() };
for i in 0..n {
// Rigid path: a's vertex carried by the interpolated similarity, plus the residual
// (b − fully transformed a) blended in linearly.
let pa = sub(a.vertices[i], ca);
let full = add(rot(pa, theta, s), cb);
let resid = sub(b.vertices[i], full);
let rigid = add(add(rot(pa, th_u, s_u), c_u), mul(resid, u));
let linear = lerp(a.vertices[i], b.vertices[i], u);
out.vertices.push(lerp(linear, rigid, w));
for (src_a, src_b, dst) in [(&a.in_tangents, &b.in_tangents, &mut out.in_tangents), (&a.out_tangents, &b.out_tangents, &mut out.out_tangents)] {
let ta = src_a[i];
let tr = sub(src_b[i], rot(ta, theta, s));
let rigid_t = add(rot(ta, th_u, s_u), mul(tr, u));
dst.push(lerp(lerp(ta, src_b[i], u), rigid_t, w));
}
}
out
}
/// Do any two vertex paths (straight segments from `a[i]` to `b[i]`) cross? (A sanity measure of
/// the correspondence: good matches don't make vertices pass through each other.)
pub fn paths_cross(a: &ShapePath, b: &ShapePath) -> usize {
let n = a.len().min(b.len());
let mut c = 0;
let inter = |p1: P, p2: P, p3: P, p4: P| {
let d = |a: P, b: P, c: P| (b[0] - a[0]) * (c[1] - a[1]) - (b[1] - a[1]) * (c[0] - a[0]);
let (d1, d2, d3, d4) = (d(p3, p4, p1), d(p3, p4, p2), d(p1, p2, p3), d(p1, p2, p4));
d1 * d2 < 0.0 && d3 * d4 < 0.0
};
for i in 0..n {
for j in i + 1..n {
if inter(a.vertices[i], b.vertices[i], a.vertices[j], b.vertices[j]) {
c += 1;
}
}
}
c
}
#[cfg(test)]
mod tests {
use super::*;
fn star(c: P, r: f64, k: usize, rot: f64) -> ShapePath {
let pts: Vec<P> = (0..2 * k)
.map(|i| {
let a = rot + i as f64 * std::f64::consts::PI / k as f64;
let rr = if i % 2 == 0 { r } else { r * 0.5 };
[c[0] + rr * a.cos(), c[1] + rr * a.sin()]
})
.collect();
ShapePath::polygon(&pts, true)
}
#[test]
fn splitting_keeps_the_outline() {
let e = ShapePath::ellipse([50.0, 50.0], 80.0, 40.0);
let m = with_vertex_count(&e, 23);
assert_eq!(m.len(), 23);
assert!((length(&m) - length(&e)).abs() < 0.05, "{} vs {}", length(&m), length(&e));
// Points on the split path lie on the ellipse.
for i in 0..m.len() {
let v = m.vertices[i];
let q = ((v[0] - 50.0) / 40.0).powi(2) + ((v[1] - 50.0) / 20.0).powi(2);
assert!((q - 1.0).abs() < 0.01, "{v:?}");
}
let r = reversed(&ShapePath::rect([0.0, 0.0], 10.0, 10.0));
assert_eq!(r.vertices[0], [-5.0, -5.0]);
assert_eq!(r.vertices[1], [-5.0, 5.0]);
}
#[test]
fn correspondence_finds_the_rotation_and_paths_dont_cross() {
// A star rotated by 2 points' worth and with its first vertex elsewhere.
let a = star([100.0, 100.0], 60.0, 5, 0.0);
let b0 = star([140.0, 110.0], 50.0, 5, 0.3);
let b = rotate(&b0, 3);
let o = InterpOpts { first_vertices_match: false, add_vertices: None, ..Default::default() };
let (ra, rb) = correspond(&a, &b, &o);
assert_eq!(ra.len(), rb.len());
// The true correspondence is recovered (b0's vertex 0 matches a's vertex 0).
assert!((rb.vertices[0][0] - b0.vertices[0][0]).abs() < 1e-9 && (rb.vertices[0][1] - b0.vertices[0][1]).abs() < 1e-9);
assert!(paths_cross(&ra, &rb) <= 1);
// With First Vertices Match the given (bad) start is kept: paths cross.
let o2 = InterpOpts { first_vertices_match: true, add_vertices: None, ..Default::default() };
let (xa, xb) = correspond(&a, &b, &o2);
assert!(paths_cross(&xa, &xb) > 5);
// In-betweens keep the vertex count and stay star-like (area between the ends).
let mid = interpolate(&ra, &rb, 0.5, &o);
assert_eq!(mid.len(), ra.len());
let (aa, ab, am) = (signed_area(&ra.vertices).abs(), signed_area(&rb.vertices).abs(), signed_area(&mid.vertices).abs());
assert!(am > ab.min(aa) * 0.8 && am < aa.max(ab) * 1.2, "{aa} {ab} {am}");
// Different vertex counts with added vertices.
let c = ShapePath::ellipse([100.0, 100.0], 120.0, 120.0);
let o3 = InterpOpts { first_vertices_match: false, ..Default::default() };
let (pa, pb) = correspond(&a, &c, &o3);
assert_eq!(pa.len(), pb.len());
assert!(pa.len() >= 25);
assert!(paths_cross(&pa, &pb) <= pa.len() / 10, "{}", paths_cross(&pa, &pb));
}
#[test]
fn rigid_interpolation_turns_instead_of_collapsing() {
let a = ShapePath::polygon(&[[-50.0, -10.0], [50.0, -10.0], [50.0, 10.0], [-50.0, 10.0]], true);
// The same bar turned 90°.
let b = ShapePath::polygon(&[[10.0, -50.0], [10.0, 50.0], [-10.0, 50.0], [-10.0, -50.0]], true);
let o = InterpOpts { bending_resistance: 100.0, ..Default::default() };
let m = interpolate(&a, &b, 0.5, &o);
let len = (m.vertices[1][0] - m.vertices[0][0]).hypot(m.vertices[1][1] - m.vertices[0][1]);
assert!((len - 100.0).abs() < 1e-6, "{len}");
let lin = interpolate(&a, &b, 0.5, &InterpOpts { linear: true, ..Default::default() });
let l2 = (lin.vertices[1][0] - lin.vertices[0][0]).hypot(lin.vertices[1][1] - lin.vertices[0][1]);
assert!(l2 < 80.0);
}
}