1 /*
2 * Copyright 2020 Google Inc.
3 *
4 * Use of this source code is governed by a BSD-style license that can be
5 * found in the LICENSE file.
6 */
7
8 #include "include/core/SkCanvas.h"
9 #include "include/core/SkFont.h"
10 #include "include/core/SkPaint.h"
11 #include "include/core/SkPath.h"
12 #include "include/private/base/SkTArray.h"
13 #include "tools/fonts/FontToolUtils.h"
14 #include "tools/viewer/ClickHandlerSlide.h"
15
16 #include <tuple>
17
18 using namespace skia_private;
19
20 // Math constants are not always defined.
21 #ifndef M_PI
22 #define M_PI 3.14159265358979323846264338327950288
23 #endif
24
25 #ifndef M_SQRT2
26 #define M_SQRT2 1.41421356237309504880168872420969808
27 #endif
28
29 constexpr static int kCenterX = 300;
30 constexpr static int kCenterY = 325;
31 constexpr static int kRadius = 250;
32
33 // This sample fits a cubic to the arc between two interactive points on a circle. It also finds the
34 // T-coordinate of max error, and outputs it and its value in pixels. (It turns out that max error
35 // always occurs at T=0.21132486540519.)
36 //
37 // Press 'E' to iteratively cut the arc in half and report the improvement in max error after each
38 // halving. (It turns out that max error improves by exactly 64x on every halving.)
39 class SampleFitCubicToCircle : public ClickHandlerSlide {
40 public:
SampleFitCubicToCircle()41 SampleFitCubicToCircle() { fName = "FitCubicToCircle"; }
load(SkScalar w,SkScalar h)42 void load(SkScalar w, SkScalar h) override { this->fitCubic(); }
43 void draw(SkCanvas*) override;
44 bool onChar(SkUnichar) override;
45
46 protected:
47 Click* onFindClickHandler(SkScalar x, SkScalar y, skui::ModifierKey) override;
48 bool onClick(Click*) override;
49
50 private:
51 void fitCubic();
52 // Coordinates of two points on the unit circle. These are the two endpoints of the arc we fit.
53 double fEndptsX[2] = {0, 1};
54 double fEndptsY[2] = {-1, 0};
55
56 // Fitted cubic and info, set by fitCubic().
57 double fControlLength; // Length of (p1 - p0) and/or (p3 - p2) in unit circle space.
58 double fMaxErrorT; // T value where the cubic diverges most from the true arc.
59 std::array<double, 4> fCubicX; // Screen space cubic control points.
60 std::array<double, 4> fCubicY;
61 double fMaxError; // Max error (in pixels) between the cubic and the screen-space arc.
62 double fTheta; // Angle of the arc. This is only used for informational purposes.
63 TArray<SkString> fInfoStrings;
64
65 class Click;
66 };
67
68 // Fits a cubic to an arc on the unit circle with endpoints (x0, y0) and (x1, y1). Using the
69 // following 3 constraints, we arrive at the formula used in the method:
70 //
71 // 1) The endpoints and tangent directions at the endpoints must match the arc.
72 // 2) The cubic must be symmetric (i.e., length(p1 - p0) == length(p3 - p2)).
73 // 3) The height of the cubic must match the height of the arc.
74 //
75 // Returns the "control length", or length of (p1 - p0) and/or (p3 - p2).
fit_cubic_to_unit_circle(double x0,double y0,double x1,double y1,std::array<double,4> * X,std::array<double,4> * Y)76 static float fit_cubic_to_unit_circle(double x0, double y0, double x1, double y1,
77 std::array<double, 4>* X, std::array<double, 4>* Y) {
78 constexpr static double kM = -4.0/3;
79 constexpr static double kA = 4*M_SQRT2/3;
80 double d = x0*x1 + y0*y1;
81 double c = (std::sqrt(1 + d) * kM + kA) / std::sqrt(1 - d);
82 *X = {x0, x0 - y0*c, x1 + y1*c, x1};
83 *Y = {y0, y0 + x0*c, y1 - x1*c, y1};
84 return c;
85 }
86
lerp(double x,double y,double T)87 static double lerp(double x, double y, double T) {
88 return x + T*(y - x);
89 }
90
91 // Evaluates the cubic and 1st and 2nd derivatives at T.
eval_cubic(double x[],double T)92 static std::tuple<double, double, double> eval_cubic(double x[], double T) {
93 // Use De Casteljau's algorithm for better accuracy and stability.
94 double ab = lerp(x[0], x[1], T);
95 double bc = lerp(x[1], x[2], T);
96 double cd = lerp(x[2], x[3], T);
97 double abc = lerp(ab, bc, T);
98 double bcd = lerp(bc, cd, T);
99 double abcd = lerp(abc, bcd, T);
100 return {abcd, 3 * (bcd - abc) /*1st derivative.*/, 6 * (cd - 2*bc + ab) /*2nd derivative.*/};
101 }
102
103 // Uses newton-raphson convergence to find the point where the provided cubic diverges most from the
104 // unit circle. i.e., the point where the derivative of error == 0. For error we use:
105 //
106 // error = x^2 + y^2 - 1
107 // error' = 2xx' + 2yy'
108 // error'' = 2xx'' + 2yy'' + 2x'^2 + 2y'^2
109 //
find_max_error_T(double cubicX[4],double cubicY[4])110 double find_max_error_T(double cubicX[4], double cubicY[4]) {
111 constexpr static double kInitialT = .25;
112 double T = kInitialT;
113 for (int i = 0; i < 64; ++i) {
114 auto [x, dx, ddx] = eval_cubic(cubicX, T);
115 auto [y, dy, ddy] = eval_cubic(cubicY, T);
116 double dError = 2*(x*dx + y*dy);
117 double ddError = 2*(x*ddx + y*ddy + dx*dx + dy*dy);
118 T -= dError / ddError;
119 }
120 return T;
121 }
122
fitCubic()123 void SampleFitCubicToCircle::fitCubic() {
124 fInfoStrings.clear();
125
126 std::array<double, 4> X, Y;
127 // "Control length" is the length of (p1 - p0) and/or (p3 - p2) in unit circle space.
128 fControlLength = fit_cubic_to_unit_circle(fEndptsX[0], fEndptsY[0], fEndptsX[1], fEndptsY[1],
129 &X, &Y);
130 fInfoStrings.push_back().printf("control length=%0.14f", fControlLength);
131
132 fMaxErrorT = find_max_error_T(X.data(), Y.data());
133 fInfoStrings.push_back().printf("max error T=%0.14f", fMaxErrorT);
134
135 for (int i = 0; i < 4; ++i) {
136 fCubicX[i] = X[i] * kRadius + kCenterX;
137 fCubicY[i] = Y[i] * kRadius + kCenterY;
138 }
139 double errX = std::get<0>(eval_cubic(fCubicX.data(), fMaxErrorT)) - kCenterX;
140 double errY = std::get<0>(eval_cubic(fCubicY.data(), fMaxErrorT)) - kCenterY;
141 fMaxError = std::sqrt(errX*errX + errY*errY) - kRadius;
142 fInfoStrings.push_back().printf("max error=%.5gpx", fMaxError);
143
144 fTheta = std::atan2(fEndptsY[1], fEndptsX[1]) - std::atan2(fEndptsY[0], fEndptsX[0]);
145 fTheta = std::abs(fTheta * 180/M_PI);
146 if (fTheta > 180) {
147 fTheta = 360 - fTheta;
148 }
149 fInfoStrings.push_back().printf("(theta=%.2f)", fTheta);
150
151 SkDebugf("\n");
152 for (const SkString& infoString : fInfoStrings) {
153 SkDebugf("%s\n", infoString.c_str());
154 }
155 }
156
draw(SkCanvas * canvas)157 void SampleFitCubicToCircle::draw(SkCanvas* canvas) {
158 canvas->clear(SK_ColorBLACK);
159
160 SkPaint circlePaint;
161 circlePaint.setColor(0x80ffffff);
162 circlePaint.setStyle(SkPaint::kStroke_Style);
163 circlePaint.setStrokeWidth(0);
164 circlePaint.setAntiAlias(true);
165 canvas->drawArc(SkRect::MakeXYWH(kCenterX - kRadius, kCenterY - kRadius, kRadius * 2,
166 kRadius * 2), 0, 360, false, circlePaint);
167
168 SkPaint cubicPaint;
169 cubicPaint.setColor(SK_ColorGREEN);
170 cubicPaint.setStyle(SkPaint::kStroke_Style);
171 cubicPaint.setStrokeWidth(10);
172 cubicPaint.setAntiAlias(true);
173 SkPath cubicPath;
174 cubicPath.moveTo(fCubicX[0], fCubicY[0]);
175 cubicPath.cubicTo(fCubicX[1], fCubicY[1], fCubicX[2], fCubicY[2], fCubicX[3], fCubicY[3]);
176 canvas->drawPath(cubicPath, cubicPaint);
177
178 SkPaint endpointsPaint;
179 endpointsPaint.setColor(SK_ColorBLUE);
180 endpointsPaint.setStrokeWidth(8);
181 endpointsPaint.setAntiAlias(true);
182 SkPoint points[2] = {{(float)fCubicX[0], (float)fCubicY[0]},
183 {(float)fCubicX[3], (float)fCubicY[3]}};
184 canvas->drawPoints(SkCanvas::kPoints_PointMode, 2, points, endpointsPaint);
185
186 SkPaint textPaint;
187 textPaint.setColor(SK_ColorWHITE);
188 constexpr static float kInfoTextSize = 16;
189 SkFont font(ToolUtils::DefaultTypeface(), kInfoTextSize);
190 int infoY = 10 + kInfoTextSize;
191 for (const SkString& infoString : fInfoStrings) {
192 canvas->drawString(infoString.c_str(), 10, infoY, font, textPaint);
193 infoY += kInfoTextSize * 3/2;
194 }
195 }
196
197 class SampleFitCubicToCircle::Click : public ClickHandlerSlide::Click {
198 public:
Click(int ptIdx)199 Click(int ptIdx) : fPtIdx(ptIdx) {}
200
doClick(SampleFitCubicToCircle * that)201 void doClick(SampleFitCubicToCircle* that) {
202 double dx = fCurr.fX - kCenterX;
203 double dy = fCurr.fY - kCenterY;
204 double l = std::sqrt(dx*dx + dy*dy);
205 that->fEndptsX[fPtIdx] = dx/l;
206 that->fEndptsY[fPtIdx] = dy/l;
207 if (that->fEndptsX[0] * that->fEndptsY[1] - that->fEndptsY[0] * that->fEndptsX[1] < 0) {
208 std::swap(that->fEndptsX[0], that->fEndptsX[1]);
209 std::swap(that->fEndptsY[0], that->fEndptsY[1]);
210 fPtIdx = 1 - fPtIdx;
211 }
212 that->fitCubic();
213 }
214
215 private:
216 int fPtIdx;
217 };
218
onFindClickHandler(SkScalar x,SkScalar y,skui::ModifierKey)219 ClickHandlerSlide::Click* SampleFitCubicToCircle::onFindClickHandler(SkScalar x, SkScalar y,
220 skui::ModifierKey) {
221 double dx0 = x - fCubicX[0];
222 double dy0 = y - fCubicY[0];
223 double dx3 = x - fCubicX[3];
224 double dy3 = y - fCubicY[3];
225 if (dx0*dx0 + dy0*dy0 < dx3*dx3 + dy3*dy3) {
226 return new Click(0);
227 } else {
228 return new Click(1);
229 }
230 }
231
onClick(ClickHandlerSlide::Click * click)232 bool SampleFitCubicToCircle::onClick(ClickHandlerSlide::Click* click) {
233 Click* myClick = (Click*)click;
234 myClick->doClick(this);
235 return true;
236 }
237
onChar(SkUnichar unichar)238 bool SampleFitCubicToCircle::onChar(SkUnichar unichar) {
239 if (unichar == 'E') {
240 constexpr static double kMaxErrorT = 0.21132486540519; // Always the same.
241 // Split the arc in half until error =~0, and report the improvement after each halving.
242 double lastError = -1;
243 for (double theta = fTheta; lastError != 0; theta /= 2) {
244 double rads = theta * M_PI/180;
245 std::array<double, 4> X, Y;
246 fit_cubic_to_unit_circle(1, 0, std::cos(rads), std::sin(rads), &X, &Y);
247 auto [x, dx, ddx] = eval_cubic(X.data(), kMaxErrorT);
248 auto [y, dy, ddy] = eval_cubic(Y.data(), kMaxErrorT);
249 double error = std::sqrt(x*x + y*y) * kRadius - kRadius;
250 if ((float)error <= 0) {
251 error = 0;
252 }
253 SkDebugf("%6.2f degrees: error= %10.5gpx", theta, error);
254 if (lastError > 0) {
255 SkDebugf(" (%17.14fx improvement)", lastError / error);
256 }
257 SkDebugf("\n");
258 lastError = error;
259 }
260 return true;
261 }
262 return false;
263 }
264
265 DEF_SLIDE(return new SampleFitCubicToCircle;)
266