Skip to content

Commit 7bd7e66

Browse files
raphlinusraphlinus
authored
Merge pull request #105 from linebender/flatten
Add BezPath flattening
2 parents 6dbb8ef + 9d580ce commit 7bd7e66

File tree

3 files changed

+286
-34
lines changed

3 files changed

+286
-34
lines changed

‎src/bezpath.rs

+183-5
Original file line numberDiff line numberDiff line change
@@ -1,5 +1,6 @@
11
//! Bézier paths (up to cubic).
22
3+
use std::iter::FromIterator;
34
use std::ops::{Mul, Range};
45

56
use arrayvec::ArrayVec;
@@ -114,6 +115,79 @@ impl BezPath {
114115
BezPath::segments_of_slice(&self.0)
115116
}
116117

118+
/// Flatten the path, invoking the callback repeatedly.
119+
///
120+
/// Flattening is the action of approximating a curve with a succession of line segments.
121+
///
122+
/// <svg xmlns="http://www.w3.org/2000/svg" viewBox="0 0 120 30" height="30mm" width="120mm">
123+
/// <path d="M26.7 24.94l.82-11.15M44.46 5.1L33.8 7.34" fill="none" stroke="#55d400" stroke-width=".5"/>
124+
/// <path d="M26.7 24.94c.97-11.13 7.17-17.6 17.76-19.84M75.27 24.94l1.13-5.5 2.67-5.48 4-4.42L88 6.7l5.02-1.6" fill="none" stroke="#000"/>
125+
/// <path d="M77.57 19.37a1.1 1.1 0 0 1-1.08 1.08 1.1 1.1 0 0 1-1.1-1.08 1.1 1.1 0 0 1 1.08-1.1 1.1 1.1 0 0 1 1.1 1.1" color="#000" fill="none" stroke="#030303" stroke-linecap="round" stroke-opacity=".5"/>
126+
/// <path d="M77.57 19.37a1.1 1.1 0 0 1-1.08 1.08 1.1 1.1 0 0 1-1.1-1.08 1.1 1.1 0 0 1 1.08-1.1 1.1 1.1 0 0 1 1.1 1.1" color="#000" fill="#fff"/>
127+
/// <path d="M80.22 13.93a1.1 1.1 0 0 1-1.1 1.1 1.1 1.1 0 0 1-1.08-1.1 1.1 1.1 0 0 1 1.1-1.08 1.1 1.1 0 0 1 1.08 1.08" color="#000" fill="none" stroke="#030303" stroke-linecap="round" stroke-opacity=".5"/>
128+
/// <path d="M80.22 13.93a1.1 1.1 0 0 1-1.1 1.1 1.1 1.1 0 0 1-1.08-1.1 1.1 1.1 0 0 1 1.1-1.08 1.1 1.1 0 0 1 1.08 1.08" color="#000" fill="#fff"/>
129+
/// <path d="M84.08 9.55a1.1 1.1 0 0 1-1.08 1.1 1.1 1.1 0 0 1-1.1-1.1 1.1 1.1 0 0 1 1.1-1.1 1.1 1.1 0 0 1 1.08 1.1" color="#000" fill="none" stroke="#030303" stroke-linecap="round" stroke-opacity=".5"/>
130+
/// <path d="M84.08 9.55a1.1 1.1 0 0 1-1.08 1.1 1.1 1.1 0 0 1-1.1-1.1 1.1 1.1 0 0 1 1.1-1.1 1.1 1.1 0 0 1 1.08 1.1" color="#000" fill="#fff"/>
131+
/// <path d="M89.1 6.66a1.1 1.1 0 0 1-1.08 1.1 1.1 1.1 0 0 1-1.08-1.1 1.1 1.1 0 0 1 1.08-1.08 1.1 1.1 0 0 1 1.1 1.08" color="#000" fill="none" stroke="#030303" stroke-linecap="round" stroke-opacity=".5"/>
132+
/// <path d="M89.1 6.66a1.1 1.1 0 0 1-1.08 1.1 1.1 1.1 0 0 1-1.08-1.1 1.1 1.1 0 0 1 1.08-1.08 1.1 1.1 0 0 1 1.1 1.08" color="#000" fill="#fff"/>
133+
/// <path d="M94.4 5a1.1 1.1 0 0 1-1.1 1.1A1.1 1.1 0 0 1 92.23 5a1.1 1.1 0 0 1 1.08-1.08A1.1 1.1 0 0 1 94.4 5" color="#000" fill="none" stroke="#030303" stroke-linecap="round" stroke-opacity=".5"/>
134+
/// <path d="M94.4 5a1.1 1.1 0 0 1-1.1 1.1A1.1 1.1 0 0 1 92.23 5a1.1 1.1 0 0 1 1.08-1.08A1.1 1.1 0 0 1 94.4 5" color="#000" fill="#fff"/>
135+
/// <path d="M76.44 25.13a1.1 1.1 0 0 1-1.1 1.1 1.1 1.1 0 0 1-1.08-1.1 1.1 1.1 0 0 1 1.1-1.1 1.1 1.1 0 0 1 1.08 1.1" color="#000" fill="none" stroke="#030303" stroke-linecap="round" stroke-opacity=".5"/>
136+
/// <path d="M76.44 25.13a1.1 1.1 0 0 1-1.1 1.1 1.1 1.1 0 0 1-1.08-1.1 1.1 1.1 0 0 1 1.1-1.1 1.1 1.1 0 0 1 1.08 1.1" color="#000" fill="#fff"/>
137+
/// <path d="M27.78 24.9a1.1 1.1 0 0 1-1.08 1.08 1.1 1.1 0 0 1-1.1-1.08 1.1 1.1 0 0 1 1.1-1.1 1.1 1.1 0 0 1 1.08 1.1" color="#000" fill="none" stroke="#030303" stroke-linecap="round" stroke-opacity=".5"/>
138+
/// <path d="M27.78 24.9a1.1 1.1 0 0 1-1.08 1.08 1.1 1.1 0 0 1-1.1-1.08 1.1 1.1 0 0 1 1.1-1.1 1.1 1.1 0 0 1 1.08 1.1" color="#000" fill="#fff"/>
139+
/// <path d="M45.4 5.14a1.1 1.1 0 0 1-1.08 1.1 1.1 1.1 0 0 1-1.1-1.1 1.1 1.1 0 0 1 1.1-1.08 1.1 1.1 0 0 1 1.1 1.08" color="#000" fill="none" stroke="#030303" stroke-linecap="round" stroke-opacity=".5"/>
140+
/// <path d="M45.4 5.14a1.1 1.1 0 0 1-1.08 1.1 1.1 1.1 0 0 1-1.1-1.1 1.1 1.1 0 0 1 1.1-1.08 1.1 1.1 0 0 1 1.1 1.08" color="#000" fill="#fff"/>
141+
/// <path d="M28.67 13.8a1.1 1.1 0 0 1-1.1 1.08 1.1 1.1 0 0 1-1.08-1.08 1.1 1.1 0 0 1 1.08-1.1 1.1 1.1 0 0 1 1.1 1.1" color="#000" fill="none" stroke="#030303" stroke-linecap="round" stroke-opacity=".5"/>
142+
/// <path d="M28.67 13.8a1.1 1.1 0 0 1-1.1 1.08 1.1 1.1 0 0 1-1.08-1.08 1.1 1.1 0 0 1 1.08-1.1 1.1 1.1 0 0 1 1.1 1.1" color="#000" fill="#fff"/>
143+
/// <path d="M35 7.32a1.1 1.1 0 0 1-1.1 1.1 1.1 1.1 0 0 1-1.08-1.1 1.1 1.1 0 0 1 1.1-1.1A1.1 1.1 0 0 1 35 7.33" color="#000" fill="none" stroke="#030303" stroke-linecap="round" stroke-opacity=".5"/>
144+
/// <path d="M35 7.32a1.1 1.1 0 0 1-1.1 1.1 1.1 1.1 0 0 1-1.08-1.1 1.1 1.1 0 0 1 1.1-1.1A1.1 1.1 0 0 1 35 7.33" color="#000" fill="#fff"/>
145+
/// <text style="line-height:6.61458302px" x="35.74" y="284.49" font-size="5.29" font-family="Sans" letter-spacing="0" word-spacing="0" fill="#b3b3b3" stroke-width=".26" transform="translate(19.595 -267)">
146+
/// <tspan x="35.74" y="284.49" font-size="10.58">→</tspan>
147+
/// </text>
148+
/// </svg>
149+
///
150+
/// The tolerance value controls the maximum distance between the curved input
151+
/// segments and their polyline approximations. (In technical terms, this is the
152+
/// Hausdorff distance). The algorithm attempts to bound this distance between
153+
/// by `tolerance` but this is not absolutely guaranteed. The appropriate value
154+
/// depends on the use, but for antialiasted rendering, a value of 0.25 has been
155+
/// determined to give good results. The number of segments tends to scale as the
156+
/// inverse square root of tolerance.
157+
///
158+
/// <svg viewBox="0 0 47.5 13.2" height="100" width="350" xmlns="http://www.w3.org/2000/svg">
159+
/// <path d="M-2.44 9.53c16.27-8.5 39.68-7.93 52.13 1.9" fill="none" stroke="#dde9af" stroke-width="4.6"/>
160+
/// <path d="M-1.97 9.3C14.28 1.03 37.36 1.7 49.7 11.4" fill="none" stroke="#00d400" stroke-width=".57" stroke-linecap="round" stroke-dasharray="4.6, 2.291434"/>
161+
/// <path d="M-1.94 10.46L6.2 6.08l28.32-1.4 15.17 6.74" fill="none" stroke="#000" stroke-width=".6"/>
162+
/// <path d="M6.83 6.57a.9.9 0 0 1-1.25.15.9.9 0 0 1-.15-1.25.9.9 0 0 1 1.25-.15.9.9 0 0 1 .15 1.25" color="#000" stroke="#000" stroke-width=".57" stroke-linecap="round" stroke-opacity=".5"/>
163+
/// <path d="M35.35 5.3a.9.9 0 0 1-1.25.15.9.9 0 0 1-.15-1.25.9.9 0 0 1 1.25-.15.9.9 0 0 1 .15 1.24" color="#000" stroke="#000" stroke-width=".6" stroke-opacity=".5"/>
164+
/// <g fill="none" stroke="#ff7f2a" stroke-width=".26">
165+
/// <path d="M20.4 3.8l.1 1.83M19.9 4.28l.48-.56.57.52M21.02 5.18l-.5.56-.6-.53" stroke-width=".2978872"/>
166+
/// </g>
167+
/// </svg>
168+
///
169+
/// The callback will be called in order with each element of the generated
170+
/// path. Because the result is made of polylines, these will be straight-line
171+
/// path elements only, no curves.
172+
///
173+
/// This algorithm is based on the blog post [Flattening quadratic Béziers]
174+
/// but with some refinements. For one, there is a more careful approximation
175+
/// at cusps. For two, the algorithm is extended to work with cubic Béziers
176+
/// as well, by first subdividing into quadratics and then computing the
177+
/// subdivision of each quadratic. However, as a clever trick, these quadratics
178+
/// are subdivided fractionally, and their endpoints are not included.
179+
///
180+
/// TODO: write a paper explaining this in more detail.
181+
///
182+
/// Note: the [`flatten`](fn.flatten.html) function provides the same
183+
/// functionality but works with slices and other [`PathEl`] iterators.
184+
///
185+
/// [Flattening quadratic Béziers]: https://raphlinus.github.io/graphics/curves/2019/12/23/flatten-quadbez.html
186+
/// [`PathEl`]: enum.PathEl.html
187+
pub fn flatten(&self, tolerance: f64, callback: impl FnMut(PathEl)) {
188+
flatten(self, tolerance, callback);
189+
}
190+
117191
// TODO: expose as pub method? Maybe should be a trait so slice.segments() works?
118192
fn segments_of_slice<'a>(slice: &'a [PathEl]) -> BezPathSegs<'a> {
119193
let first = match slice.get(0) {
@@ -194,17 +268,17 @@ impl BezPath {
194268
}
195269
}
196270

197-
impl std::iter::FromIterator<PathEl> for BezPath {
271+
impl FromIterator<PathEl> for BezPath {
198272
fn from_iter<T: IntoIterator<Item = PathEl>>(iter: T) -> Self {
199273
let el_vec: Vec<_> = iter.into_iter().collect();
200274
BezPath::from_vec(el_vec)
201275
}
202276
}
203277

204-
// this has weird semantics; signature assumes taking ownership but impl'd on a reference
205-
// NOTE: after removing this, we should impl IntoIterator for BezPath (with no reference)
206-
// and that impl should just call `self.0.into_iter()`
207-
#[deprecated(since = "0.5.6", note = "use BezPath::iter instead")]
278+
/// Allow iteration over references to `BezPath`.
279+
///
280+
/// Note: the semantics are slightly different than simply iterating over the
281+
/// slice, as it returns `PathEl` items, rather than references.
208282
impl<'a> IntoIterator for &'a BezPath {
209283
type Item = PathEl;
210284
type IntoIter = std::iter::Cloned<std::slice::Iter<'a, PathEl>>;
@@ -214,6 +288,110 @@ impl<'a> IntoIterator for &'a BezPath {
214288
}
215289
}
216290

291+
impl IntoIterator for BezPath {
292+
type Item = PathEl;
293+
type IntoIter = std::vec::IntoIter<PathEl>;
294+
295+
fn into_iter(self) -> Self::IntoIter {
296+
self.0.into_iter()
297+
}
298+
}
299+
300+
/// Proportion of tolerance budget that goes to cubic to quadratic conversion.
301+
const TO_QUAD_TOL: f64 = 0.1;
302+
303+
/// Flatten the path, invoking the callback repeatedly.
304+
///
305+
/// See [`BezPath::flatten`](struct.BezPath.html#method.flatten) for more discussion.
306+
/// This signature is a bit more general, allowing flattening of `&[PathEl]` slices
307+
/// and other iterators yielding `PathEl`.
308+
pub fn flatten(
309+
path: impl IntoIterator<Item = PathEl>,
310+
tolerance: f64,
311+
mut callback: impl FnMut(PathEl),
312+
) {
313+
let sqrt_tol = tolerance.sqrt();
314+
let mut last_pt = None;
315+
let mut quad_buf = Vec::new();
316+
for el in path {
317+
match el {
318+
PathEl::MoveTo(p) => {
319+
last_pt = Some(p);
320+
callback(PathEl::MoveTo(p));
321+
}
322+
PathEl::LineTo(p) => {
323+
last_pt = Some(p);
324+
callback(PathEl::LineTo(p));
325+
}
326+
PathEl::QuadTo(p1, p2) => {
327+
if let Some(p0) = last_pt {
328+
let q = QuadBez::new(p0, p1, p2);
329+
let params = q.estimate_subdiv(sqrt_tol);
330+
let n = ((0.5 * params.val / sqrt_tol).ceil() as usize).max(1);
331+
let step = 1.0 / (n as f64);
332+
for i in 1..(n - 1) {
333+
let u = (i as f64) * step;
334+
let t = q.determine_subdiv_t(&params, u);
335+
let p = q.eval(t);
336+
callback(PathEl::LineTo(p));
337+
}
338+
callback(PathEl::LineTo(p2));
339+
}
340+
last_pt = Some(p2);
341+
}
342+
PathEl::CurveTo(p1, p2, p3) => {
343+
if let Some(p0) = last_pt {
344+
let c = CubicBez::new(p0, p1, p2, p3);
345+
346+
// Subdivide into quadratics, and estimate the number of
347+
// subdivisions required for each, summing to arrive at an
348+
// estimate for the number of subdivisions for the cubic.
349+
// Also retain these parameters for later.
350+
let iter = c.to_quads(tolerance * TO_QUAD_TOL);
351+
quad_buf.clear();
352+
quad_buf.reserve(iter.size_hint().0);
353+
let sqrt_remain_tol = sqrt_tol * (1.0 - TO_QUAD_TOL).sqrt();
354+
let mut sum = 0.0;
355+
for (_, _, q) in iter {
356+
let params = q.estimate_subdiv(sqrt_remain_tol);
357+
sum += params.val;
358+
quad_buf.push((q, params));
359+
}
360+
let n = ((0.5 * sum / sqrt_remain_tol).ceil() as usize).max(1);
361+
362+
// Iterate through the quadratics, outputting the points of
363+
// subdivisions that fall within that quadratic.
364+
let step = sum / (n as f64);
365+
let mut i = 1;
366+
let mut val_sum = 0.0;
367+
for (q, params) in &quad_buf {
368+
let mut target = (i as f64) * step;
369+
let recip_val = params.val.recip();
370+
while target < val_sum + params.val {
371+
let u = (target - val_sum) * recip_val;
372+
let t = q.determine_subdiv_t(&params, u);
373+
let p = q.eval(t);
374+
callback(PathEl::LineTo(p));
375+
i += 1;
376+
if i == n + 1 {
377+
break;
378+
}
379+
target = (i as f64) * step;
380+
}
381+
val_sum += params.val;
382+
}
383+
callback(PathEl::LineTo(p3));
384+
}
385+
last_pt = Some(p3);
386+
}
387+
PathEl::ClosePath => {
388+
last_pt = None;
389+
callback(PathEl::ClosePath);
390+
}
391+
}
392+
}
393+
}
394+
217395
impl Mul<PathEl> for Affine {
218396
type Output = PathEl;
219397

‎src/cubicbez.rs

+34-29
Original file line numberDiff line numberDiff line change
@@ -26,8 +26,8 @@ pub struct CubicBez {
2626
/// An iterator which produces quadratic Bézier segments.
2727
struct ToQuads {
2828
c: CubicBez,
29-
max_hypot2: f64,
30-
t: f64,
29+
i: usize,
30+
n: usize,
3131
}
3232

3333
impl CubicBez {
@@ -49,16 +49,29 @@ impl CubicBez {
4949
///
5050
/// Note that the resulting quadratic Béziers are not in general G1 continuous;
5151
/// they are optimized for minimizing distance error.
52+
///
53+
/// Also note that this iterator may produce zero quadratics when the control points
54+
/// are equally spaced and co-linear.
5255
#[inline]
5356
pub fn to_quads(&self, accuracy: f64) -> impl Iterator<Item = (f64, f64, QuadBez)> {
57+
// The maximum error, as a vector from the cubic to the best approximating
58+
// quadratic, is proportional to the third derivative, which is constant
59+
// across the segment. Thus, the error scales down as the third power of
60+
// the number of subdivisions. Our strategy then is to subdivide `t` evenly.
61+
//
62+
// This is an overestimate of the error because only the component
63+
// perpendicular to the first derivative is important. But the simplicity is
64+
// appealing.
65+
5466
// This magic number is the square of 36 / sqrt(3).
5567
// See: http://caffeineowl.com/graphics/2d/vectorial/cubic2quad01.html
5668
let max_hypot2 = 432.0 * accuracy * accuracy;
57-
ToQuads {
58-
c: *self,
59-
max_hypot2,
60-
t: 0.0,
61-
}
69+
let p1x2 = 3.0 * self.p1.to_vec2() - self.p0.to_vec2();
70+
let p2x2 = 3.0 * self.p2.to_vec2() - self.p3.to_vec2();
71+
let err = (p2x2 - p1x2).hypot2();
72+
let n = (err / max_hypot2).powf(1. / 6.0).ceil() as usize;
73+
74+
ToQuads { c: *self, n, i: 0 }
6275
}
6376
}
6477

@@ -245,30 +258,22 @@ impl Iterator for ToQuads {
245258
type Item = (f64, f64, QuadBez);
246259

247260
fn next(&mut self) -> Option<(f64, f64, QuadBez)> {
248-
let t0 = self.t;
249-
let mut t1 = 1.0;
250-
if t0 == t1 {
261+
if self.i == self.n {
251262
return None;
252263
}
253-
loop {
254-
let seg = self.c.subsegment(t0..t1);
255-
// Compute error for candidate quadratic.
256-
let p1x2 = 3.0 * seg.p1.to_vec2() - seg.p0.to_vec2();
257-
let p2x2 = 3.0 * seg.p2.to_vec2() - seg.p3.to_vec2();
258-
let err = (p2x2 - p1x2).hypot2();
259-
if err < self.max_hypot2 {
260-
let result = QuadBez::new(seg.p0, ((p1x2 + p2x2) / 4.0).to_point(), seg.p3);
261-
self.t = t1;
262-
return Some((t0, t1, result));
263-
} else {
264-
let shrink = if t1 == 1.0 && err < 64.0 * self.max_hypot2 {
265-
0.5
266-
} else {
267-
0.999_999 * (self.max_hypot2 / err).powf(1. / 6.0)
268-
};
269-
t1 = t0 + shrink * (t1 - t0);
270-
}
271-
}
264+
let t0 = self.i as f64 / self.n as f64;
265+
let t1 = (self.i + 1) as f64 / self.n as f64;
266+
let seg = self.c.subsegment(t0..t1);
267+
let p1x2 = 3.0 * seg.p1.to_vec2() - seg.p0.to_vec2();
268+
let p2x2 = 3.0 * seg.p2.to_vec2() - seg.p3.to_vec2();
269+
let result = QuadBez::new(seg.p0, ((p1x2 + p2x2) / 4.0).to_point(), seg.p3);
270+
self.i += 1;
271+
Some((t0, t1, result))
272+
}
273+
274+
fn size_hint(&self) -> (usize, Option<usize>) {
275+
let remaining = self.n - self.i;
276+
(remaining, Some(remaining))
272277
}
273278
}
274279

‎src/quadbez.rs

+69
Original file line numberDiff line numberDiff line change
@@ -44,6 +44,75 @@ impl QuadBez {
4444
self.p2,
4545
)
4646
}
47+
48+
/// Estimate the number of subdivisions for flattening.
49+
pub(crate) fn estimate_subdiv(&self, sqrt_tol: f64) -> FlattenParams {
50+
// Determine transformation to $y = x^2$ parabola.
51+
let d01 = self.p1 - self.p0;
52+
let d12 = self.p2 - self.p1;
53+
let dd = d01 - d12;
54+
let cross = (self.p2 - self.p0).cross(dd);
55+
let x0 = d01.dot(dd) * cross.recip();
56+
let x2 = d12.dot(dd) * cross.recip();
57+
let scale = (cross / (dd.hypot() * (x2 - x0))).abs();
58+
59+
// Compute number of subdivisions needed.
60+
let a0 = approx_parabola_integral(x0);
61+
let a2 = approx_parabola_integral(x2);
62+
let val = if scale.is_finite() {
63+
let da = (a2 - a0).abs();
64+
let sqrt_scale = scale.sqrt();
65+
if x0.signum() == x2.signum() {
66+
da * sqrt_scale
67+
} else {
68+
// Handle cusp case (segment contains curvature maximum)
69+
let xmin = sqrt_tol / sqrt_scale;
70+
sqrt_tol * da / approx_parabola_integral(xmin)
71+
}
72+
} else {
73+
0.0
74+
};
75+
let u0 = approx_parabola_inv_integral(a0);
76+
let u2 = approx_parabola_inv_integral(a2);
77+
let uscale = (u2 - u0).recip();
78+
FlattenParams {
79+
a0,
80+
a2,
81+
u0,
82+
uscale,
83+
val,
84+
}
85+
}
86+
87+
// Maps a value from 0..1 to 0..1.
88+
pub(crate) fn determine_subdiv_t(&self, params: &FlattenParams, x: f64) -> f64 {
89+
let a = params.a0 + (params.a2 - params.a0) * x;
90+
let u = approx_parabola_inv_integral(a);
91+
(u - params.u0) * params.uscale
92+
}
93+
}
94+
95+
pub(crate) struct FlattenParams {
96+
a0: f64,
97+
a2: f64,
98+
u0: f64,
99+
uscale: f64,
100+
/// The number of subdivisions * 2 * sqrt_tol.
101+
pub(crate) val: f64,
102+
}
103+
104+
/// An approximation to $\int (1 + 4x^2) ^ -0.25 dx$
105+
///
106+
/// This is used for flattening curves.
107+
fn approx_parabola_integral(x: f64) -> f64 {
108+
const D: f64 = 0.67;
109+
x / (1.0 - D + (D.powi(4) + 0.25 * x * x).sqrt().sqrt())
110+
}
111+
112+
/// An approximation to the inverse parabola integral.
113+
fn approx_parabola_inv_integral(x: f64) -> f64 {
114+
const B: f64 = 0.39;
115+
x * (1.0 - B + (B * B + 0.25 * x * x).sqrt())
47116
}
48117

49118
impl ParamCurve for QuadBez {

0 commit comments

Comments
 (0)