1// Copyright 2018 the Kurbo Authors
2// SPDX-License-Identifier: Apache-2.0 OR MIT
3
4//! A trait for curves parametrized by a scalar.
5
6use core::ops::Range;
7
8use arrayvec::ArrayVec;
9
10use crate::{common, Point, Rect};
11
12#[cfg(not(feature = "std"))]
13use crate::common::FloatFuncs;
14
15/// A default value for methods that take an 'accuracy' argument.
16///
17/// This value is intended to be suitable for general-purpose use, such as
18/// 2d graphics.
19pub const DEFAULT_ACCURACY: f64 = 1e-6;
20
21/// A curve parametrized by a scalar.
22///
23/// If the result is interpreted as a point, this represents a curve.
24/// But the result can be interpreted as a vector as well.
25pub trait ParamCurve: Sized {
26 /// Evaluate the curve at parameter `t`.
27 ///
28 /// Generally `t` is in the range [0..1].
29 fn eval(&self, t: f64) -> Point;
30
31 /// Get a subsegment of the curve for the given parameter range.
32 fn subsegment(&self, range: Range<f64>) -> Self;
33
34 /// Subdivide into (roughly) halves.
35 #[inline]
36 fn subdivide(&self) -> (Self, Self) {
37 (self.subsegment(0.0..0.5), self.subsegment(0.5..1.0))
38 }
39
40 /// The start point.
41 fn start(&self) -> Point {
42 self.eval(0.0)
43 }
44
45 /// The end point.
46 fn end(&self) -> Point {
47 self.eval(1.0)
48 }
49}
50
51// TODO: I might not want to have separate traits for all these.
52
53/// A differentiable parametrized curve.
54pub trait ParamCurveDeriv {
55 /// The parametric curve obtained by taking the derivative of this one.
56 type DerivResult: ParamCurve;
57
58 /// The derivative of the curve.
59 ///
60 /// Note that the type of the return value is somewhat inaccurate, as
61 /// the derivative of a curve (mapping of param to point) is a mapping
62 /// of param to vector. We choose to accept this rather than have a
63 /// more complex type scheme.
64 fn deriv(&self) -> Self::DerivResult;
65
66 /// Estimate arclength using Gaussian quadrature.
67 ///
68 /// The coefficients are assumed to cover the range (-1..1), which is
69 /// traditional.
70 #[inline]
71 fn gauss_arclen(&self, coeffs: &[(f64, f64)]) -> f64 {
72 let d = self.deriv();
73 coeffs
74 .iter()
75 .map(|(wi, xi)| wi * d.eval(0.5 * (xi + 1.0)).to_vec2().hypot())
76 .sum::<f64>()
77 * 0.5
78 }
79}
80
81/// A parametrized curve that can have its arc length measured.
82pub trait ParamCurveArclen: ParamCurve {
83 /// The arc length of the curve.
84 ///
85 /// The result is accurate to the given accuracy (subject to
86 /// roundoff errors for ridiculously low values). Compute time
87 /// may vary with accuracy, if the curve needs to be subdivided.
88 fn arclen(&self, accuracy: f64) -> f64;
89
90 /// Solve for the parameter that has the given arc length from the start.
91 ///
92 /// This implementation uses the IPT method, as provided by
93 /// [`common::solve_itp`]. This is as robust as bisection but
94 /// typically converges faster. In addition, the method takes
95 /// care to compute arc lengths of increasingly smaller segments
96 /// of the curve, as that is likely faster than repeatedly
97 /// computing the arc length of the segment starting at t=0.
98 fn inv_arclen(&self, arclen: f64, accuracy: f64) -> f64 {
99 if arclen <= 0.0 {
100 return 0.0;
101 }
102 let total_arclen = self.arclen(accuracy);
103 if arclen >= total_arclen {
104 return 1.0;
105 }
106 let mut t_last = 0.0;
107 let mut arclen_last = 0.0;
108 let epsilon = accuracy / total_arclen;
109 let n = 1.0 - epsilon.log2().ceil().min(0.0);
110 let inner_accuracy = accuracy / n;
111 let f = |t: f64| {
112 let (range, dir) = if t > t_last {
113 (t_last..t, 1.0)
114 } else {
115 (t..t_last, -1.0)
116 };
117 let arc = self.subsegment(range).arclen(inner_accuracy);
118 arclen_last += arc * dir;
119 t_last = t;
120 arclen_last - arclen
121 };
122 common::solve_itp(f, 0.0, 1.0, epsilon, 1, 0.2, -arclen, total_arclen - arclen)
123 }
124}
125
126/// A parametrized curve that can have its signed area measured.
127pub trait ParamCurveArea {
128 /// Compute the signed area under the curve.
129 ///
130 /// For a closed path, the signed area of the path is the sum of signed
131 /// areas of the segments. This is a variant of the "shoelace formula."
132 /// See:
133 /// <https://github.com/Pomax/bezierinfo/issues/44> and
134 /// <http://ich.deanmcnamee.com/graphics/2016/03/30/CurveArea.html>
135 ///
136 /// This can be computed exactly for Béziers thanks to Green's theorem,
137 /// and also for simple curves such as circular arcs. For more exotic
138 /// curves, it's probably best to subdivide to cubics. We leave that
139 /// to the caller, which is why we don't give an accuracy param here.
140 fn signed_area(&self) -> f64;
141}
142
143/// The nearest position on a curve to some point.
144///
145/// This is returned by [`ParamCurveNearest::nearest`]
146#[derive(Debug, Clone, Copy)]
147pub struct Nearest {
148 /// The square of the distance from the nearest position on the curve
149 /// to the given point.
150 pub distance_sq: f64,
151 /// The position on the curve of the nearest point, as a parameter.
152 ///
153 /// To resolve this to a [`Point`], use [`ParamCurve::eval`].
154 pub t: f64,
155}
156
157/// A parametrized curve that reports the nearest point.
158pub trait ParamCurveNearest {
159 /// Find the position on the curve that is nearest to the given point.
160 ///
161 /// This returns a [`Nearest`] struct that contains information about
162 /// the position.
163 fn nearest(&self, p: Point, accuracy: f64) -> Nearest;
164}
165
166/// A parametrized curve that reports its curvature.
167pub trait ParamCurveCurvature: ParamCurveDeriv
168where
169 Self::DerivResult: ParamCurveDeriv,
170{
171 /// Compute the signed curvature at parameter `t`.
172 #[inline]
173 fn curvature(&self, t: f64) -> f64 {
174 let deriv: impl ParamCurveDeriv = self.deriv();
175 let deriv2: ::DerivResult = deriv.deriv();
176 let d: Vec2 = deriv.eval(t).to_vec2();
177 let d2: Vec2 = deriv2.eval(t).to_vec2();
178 // TODO: What's the convention for sign? I think it should match signed
179 // area - a positive area curve should have positive curvature.
180 d2.cross(d) * d.hypot2().powf(-1.5)
181 }
182}
183
184/// The maximum number of extrema that can be reported in the `ParamCurveExtrema` trait.
185///
186/// This is 4 to support cubic Béziers. If other curves are used, they should be
187/// subdivided to limit the number of extrema.
188pub const MAX_EXTREMA: usize = 4;
189
190/// A parametrized curve that reports its extrema.
191pub trait ParamCurveExtrema: ParamCurve {
192 /// Compute the extrema of the curve.
193 ///
194 /// Only extrema within the interior of the curve count.
195 /// At most four extrema can be reported, which is sufficient for
196 /// cubic Béziers.
197 ///
198 /// The extrema should be reported in increasing parameter order.
199 fn extrema(&self) -> ArrayVec<f64, MAX_EXTREMA>;
200
201 /// Return parameter ranges, each of which is monotonic within the range.
202 fn extrema_ranges(&self) -> ArrayVec<Range<f64>, { MAX_EXTREMA + 1 }> {
203 let mut result = ArrayVec::new();
204 let mut t0 = 0.0;
205 for t in self.extrema() {
206 result.push(t0..t);
207 t0 = t;
208 }
209 result.push(t0..1.0);
210 result
211 }
212
213 /// The smallest rectangle that encloses the curve in the range (0..1).
214 fn bounding_box(&self) -> Rect {
215 let mut bbox = Rect::from_points(self.start(), self.end());
216 for t in self.extrema() {
217 bbox = bbox.union_pt(self.eval(t));
218 }
219 bbox
220 }
221}
222