1// Copyright 2022 the Kurbo Authors
2// SPDX-License-Identifier: Apache-2.0 OR MIT
3
4//! Simplification of a Bézier path.
5//!
6//! This module is currently experimental.
7//!
8//! The methods in this module create a `SimplifyBezPath` object, which can then
9//! be fed to [`fit_to_bezpath`] or [`fit_to_bezpath_opt`] depending on the degree
10//! of optimization desired.
11//!
12//! The implementation uses a number of techniques to achieve high performance and
13//! accuracy. The parameter (generally written `t`) evenly divides the curve segments
14//! in the original, so sampling can be done in constant time. The derivatives are
15//! computed analytically, as that is straightforward with Béziers.
16//!
17//! The areas and moments are computed analytically (using Green's theorem), and
18//! the prefix sum is stored. Thus, it is possible to analytically compute the area
19//! and moment of any subdivision of the curve, also in constant time, by taking
20//! the difference of two stored prefix sum values, then fixing up the subsegments.
21//!
22//! A current limitation (hoped to be addressed in the future) is that non-regular
23//! cubic segments may have tangents computed incorrectly. This can easily happen,
24//! for example when setting a control point equal to an endpoint.
25//!
26//! In addition, this method does not report corners (adjoining segments where the
27//! tangents are not continuous). It is not clear whether it's best to handle such
28//! cases here, or in client code.
29//!
30//! [`fit_to_bezpath`]: crate::fit_to_bezpath
31//! [`fit_to_bezpath_opt`]: crate::fit_to_bezpath_opt
32
33use alloc::vec::Vec;
34
35use core::ops::Range;
36
37#[cfg(not(feature = "std"))]
38use crate::common::FloatFuncs;
39
40use crate::{
41 fit_to_bezpath, fit_to_bezpath_opt, BezPath, CubicBez, CurveFitSample, Line, ParamCurve,
42 ParamCurveDeriv, ParamCurveFit, PathEl, PathSeg, Point, QuadBez, Vec2,
43};
44
45/// A Bézier path which has been prepared for simplification.
46///
47/// See the [module-level documentation] for a bit more discussion of the approach,
48/// and how this struct is to be used.
49///
50/// [module-level documentation]: crate::simplify
51pub struct SimplifyBezPath(Vec<SimplifyCubic>);
52
53struct SimplifyCubic {
54 c: CubicBez,
55 // The inclusive prefix sum of the moment integrals
56 moments: (f64, f64, f64),
57}
58
59/// Options for simplifying paths.
60pub struct SimplifyOptions {
61 /// The tangent of the minimum angle below which the path is considered smooth.
62 angle_thresh: f64,
63 opt_level: SimplifyOptLevel,
64}
65
66/// Optimization level for simplification.
67pub enum SimplifyOptLevel {
68 /// Subdivide; faster but not as optimized results.
69 Subdivide,
70 /// Optimize subdivision points.
71 Optimize,
72}
73
74impl Default for SimplifyOptions {
75 fn default() -> Self {
76 let opt_level: SimplifyOptLevel = SimplifyOptLevel::Subdivide;
77 SimplifyOptions {
78 angle_thresh: 1e-3,
79 opt_level,
80 }
81 }
82}
83
84#[doc(hidden)]
85/// Compute moment integrals.
86///
87/// This is exposed for testing purposes but is an internal detail. We can
88/// add to the public, documented interface if there is a use case.
89pub fn moment_integrals(c: CubicBez) -> (f64, f64, f64) {
90 let (x0, y0) = (c.p0.x, c.p0.y);
91 let (x1, y1) = (c.p1.x - x0, c.p1.y - y0);
92 let (x2, y2) = (c.p2.x - x0, c.p2.y - y0);
93 let (x3, y3) = (c.p3.x - x0, c.p3.y - y0);
94
95 let r0 = 3. * x1;
96 let r1 = 3. * y1;
97 let r2 = x2 * y3;
98 let r3 = x3 * y2;
99 let r4 = x3 * y3;
100 let r5 = 27. * y1;
101 let r6 = x1 * x2;
102 let r7 = 27. * y2;
103 let r8 = 45. * r2;
104 let r9 = 18. * x3;
105 let r10 = x1 * y1;
106 let r11 = 30. * x1;
107 let r12 = 45. * x3;
108 let r13 = x2 * y1;
109 let r14 = 45. * r3;
110 let r15 = x1.powi(2);
111 let r16 = 18. * y3;
112 let r17 = x2.powi(2);
113 let r18 = 45. * y3;
114 let r19 = x3.powi(2);
115 let r20 = 30. * y1;
116 let r21 = y2.powi(2);
117 let r22 = y3.powi(2);
118 let r23 = y1.powi(2);
119 let a = -r0 * y2 - r0 * y3 + r1 * x2 + r1 * x3 - 6. * r2 + 6. * r3 + 10. * r4;
120
121 // Scale and add chord
122 let lift = x3 * y0;
123 let area = a * 0.05 + lift;
124 let x = r10 * r9 - r11 * r4 + r12 * r13 + r14 * x2 - r15 * r16 - r15 * r7 - r17 * r18
125 + r17 * r5
126 + r19 * r20
127 + 105. * r19 * y2
128 + 280. * r19 * y3
129 - 105. * r2 * x3
130 + r5 * r6
131 - r6 * r7
132 - r8 * x1;
133 let y = -r10 * r16 - r10 * r7 - r11 * r22 + r12 * r21 + r13 * r7 + r14 * y1 - r18 * x1 * y2
134 + r20 * r4
135 - 27. * r21 * x1
136 - 105. * r22 * x2
137 + 140. * r22 * x3
138 + r23 * r9
139 + 27. * r23 * x2
140 + 105. * r3 * y3
141 - r8 * y2;
142
143 let mx = x * (1. / 840.) + x0 * area + 0.5 * x3 * lift;
144 let my = y * (1. / 420.) + y0 * a * 0.1 + y0 * lift;
145
146 (area, mx, my)
147}
148
149impl SimplifyBezPath {
150 /// Set up a new Bézier path for simplification.
151 ///
152 /// Currently this is not dealing with discontinuities at all, but it
153 /// could be extended to do so.
154 pub fn new(path: impl IntoIterator<Item = PathEl>) -> Self {
155 let (mut a, mut x, mut y) = (0.0, 0.0, 0.0);
156 let els = crate::segments(path)
157 .map(|seg| {
158 let c = seg.to_cubic();
159 let (ai, xi, yi) = moment_integrals(c);
160 a += ai;
161 x += xi;
162 y += yi;
163 SimplifyCubic {
164 c,
165 moments: (a, x, y),
166 }
167 })
168 .collect();
169 SimplifyBezPath(els)
170 }
171
172 /// Resolve a `t` value to a cubic.
173 ///
174 /// Also return the resulting `t` value for the selected cubic.
175 fn scale(&self, t: f64) -> (usize, f64) {
176 let t_scale = t * self.0.len() as f64;
177 let t_floor = t_scale.floor();
178 (t_floor as usize, t_scale - t_floor)
179 }
180
181 fn moment_integrals(&self, i: usize, range: Range<f64>) -> (f64, f64, f64) {
182 if range.end == range.start {
183 (0.0, 0.0, 0.0)
184 } else {
185 moment_integrals(self.0[i].c.subsegment(range))
186 }
187 }
188}
189
190impl ParamCurveFit for SimplifyBezPath {
191 fn sample_pt_deriv(&self, t: f64) -> (Point, Vec2) {
192 let (mut i, mut t0) = self.scale(t);
193 let n = self.0.len();
194 if i == n {
195 i -= 1;
196 t0 = 1.0;
197 }
198 let c = self.0[i].c;
199 (c.eval(t0), c.deriv().eval(t0).to_vec2() * n as f64)
200 }
201
202 fn sample_pt_tangent(&self, t: f64, _: f64) -> CurveFitSample {
203 let (mut i, mut t0) = self.scale(t);
204 if i == self.0.len() {
205 i -= 1;
206 t0 = 1.0;
207 }
208 let c = self.0[i].c;
209 let p = c.eval(t0);
210 let tangent = c.deriv().eval(t0).to_vec2();
211 CurveFitSample { p, tangent }
212 }
213
214 // We could use the default implementation, but provide our own, mostly
215 // because it is possible to efficiently provide an analytically accurate
216 // answer.
217 fn moment_integrals(&self, range: Range<f64>) -> (f64, f64, f64) {
218 let (i0, t0) = self.scale(range.start);
219 let (i1, t1) = self.scale(range.end);
220 if i0 == i1 {
221 self.moment_integrals(i0, t0..t1)
222 } else {
223 let (a0, x0, y0) = self.moment_integrals(i0, t0..1.0);
224 let (a1, x1, y1) = self.moment_integrals(i1, 0.0..t1);
225 let (mut a, mut x, mut y) = (a0 + a1, x0 + x1, y0 + y1);
226 if i1 > i0 + 1 {
227 let (a2, x2, y2) = self.0[i0].moments;
228 let (a3, x3, y3) = self.0[i1 - 1].moments;
229 a += a3 - a2;
230 x += x3 - x2;
231 y += y3 - y2;
232 }
233 (a, x, y)
234 }
235 }
236
237 fn break_cusp(&self, _: Range<f64>) -> Option<f64> {
238 None
239 }
240}
241
242#[derive(Default)]
243struct SimplifyState {
244 queue: BezPath,
245 result: BezPath,
246 needs_moveto: bool,
247}
248
249impl SimplifyState {
250 fn add_seg(&mut self, seg: PathSeg) {
251 if self.queue.is_empty() {
252 self.queue.move_to(seg.start());
253 }
254 match seg {
255 PathSeg::Line(l) => self.queue.line_to(l.p1),
256 PathSeg::Quad(q) => self.queue.quad_to(q.p1, q.p2),
257 PathSeg::Cubic(c) => self.queue.curve_to(c.p1, c.p2, c.p3),
258 }
259 }
260
261 fn flush(&mut self, accuracy: f64, options: &SimplifyOptions) {
262 if self.queue.is_empty() {
263 return;
264 }
265 if self.queue.elements().len() == 2 {
266 // Queue is just one segment (count is moveto + primitive)
267 // Just output the segment, no simplification is possible.
268 self.result
269 .extend(self.queue.iter().skip(!self.needs_moveto as usize));
270 } else {
271 let s = SimplifyBezPath::new(&self.queue);
272 let b = match options.opt_level {
273 SimplifyOptLevel::Subdivide => fit_to_bezpath(&s, accuracy),
274 SimplifyOptLevel::Optimize => fit_to_bezpath_opt(&s, accuracy),
275 };
276 self.result
277 .extend(b.iter().skip(!self.needs_moveto as usize));
278 }
279 self.needs_moveto = false;
280 self.queue.truncate(0);
281 }
282}
283
284/// Simplify a Bézier path.
285///
286/// This function simplifies an arbitrary Bézier path; it is designed to handle
287/// multiple subpaths and also corners.
288///
289/// The underlying curve-fitting approach works best if the source path is very
290/// smooth. If it contains higher frequency noise, then results may be poor, as
291/// the resulting curve matches the original with G1 continuity at each subdivision
292/// point, and also preserves the area. For such inputs, consider some form of
293/// smoothing or low-pass filtering before simplification. In particular, if the
294/// input is derived from a sequence of points, consider fitting a smooth spline.
295///
296/// We may add such capabilities in the future, possibly as opt-in smoothing
297/// specified through the options.
298pub fn simplify_bezpath(
299 path: impl IntoIterator<Item = PathEl>,
300 accuracy: f64,
301 options: &SimplifyOptions,
302) -> BezPath {
303 let mut last_pt = None;
304 let mut last_seg: Option<PathSeg> = None;
305 let mut state = SimplifyState::default();
306 for el in path {
307 let mut this_seg = None;
308 match el {
309 PathEl::MoveTo(p) => {
310 state.flush(accuracy, options);
311 state.needs_moveto = true;
312 last_pt = Some(p);
313 }
314 PathEl::LineTo(p) => {
315 let last = last_pt.unwrap();
316 if last == p {
317 continue;
318 }
319 this_seg = Some(PathSeg::Line(Line::new(last, p)));
320 }
321 PathEl::QuadTo(p1, p2) => {
322 let last = last_pt.unwrap();
323 if last == p1 && last == p2 {
324 continue;
325 }
326 this_seg = Some(PathSeg::Quad(QuadBez::new(last, p1, p2)));
327 }
328 PathEl::CurveTo(p1, p2, p3) => {
329 let last = last_pt.unwrap();
330 if last == p1 && last == p2 && last == p3 {
331 continue;
332 }
333 this_seg = Some(PathSeg::Cubic(CubicBez::new(last, p1, p2, p3)));
334 }
335 PathEl::ClosePath => {
336 state.flush(accuracy, options);
337 state.result.close_path();
338 state.needs_moveto = true;
339 last_seg = None;
340 continue;
341 }
342 }
343 if let Some(seg) = this_seg {
344 if let Some(last) = last_seg {
345 let last_tan = last.tangents().1;
346 let this_tan = seg.tangents().0;
347 if last_tan.cross(this_tan).abs()
348 > last_tan.dot(this_tan).abs() * options.angle_thresh
349 {
350 state.flush(accuracy, options);
351 }
352 }
353 last_pt = Some(seg.end());
354 state.add_seg(seg);
355 }
356 last_seg = this_seg;
357 }
358 state.flush(accuracy, options);
359 state.result
360}
361
362impl SimplifyOptions {
363 /// Set optimization level.
364 pub fn opt_level(mut self, level: SimplifyOptLevel) -> Self {
365 self.opt_level = level;
366 self
367 }
368
369 /// Set angle threshold.
370 ///
371 /// The tangent of the angle below which joins are considered smooth and
372 /// not corners. The default is approximately 1 milliradian.
373 pub fn angle_thresh(mut self, thresh: f64) -> Self {
374 self.angle_thresh = thresh;
375 self
376 }
377}
378
379#[cfg(test)]
380mod tests {
381 use crate::BezPath;
382
383 use super::{simplify_bezpath, SimplifyOptions};
384
385 #[test]
386 fn simplify_lines_corner() {
387 // Make sure lines are passed through unchanged if there is a corner.
388 let mut path = BezPath::new();
389 path.move_to((1., 2.));
390 path.line_to((3., 4.));
391 path.line_to((10., 5.));
392 let options = SimplifyOptions::default();
393 let simplified = simplify_bezpath(path.clone(), 1.0, &options);
394 assert_eq!(path, simplified);
395 }
396}
397