1
2#include "edge-segments.h"
3
4#include "arithmetics.hpp"
5#include "equation-solver.h"
6
7namespace msdfgen {
8
9void EdgeSegment::distanceToPseudoDistance(SignedDistance &distance, Point2 origin, double param) const {
10 if (param < 0) {
11 Vector2 dir = direction(0).normalize();
12 Vector2 aq = origin-point(0);
13 double ts = dotProduct(aq, dir);
14 if (ts < 0) {
15 double pseudoDistance = crossProduct(aq, dir);
16 if (fabs(pseudoDistance) <= fabs(distance.distance)) {
17 distance.distance = pseudoDistance;
18 distance.dot = 0;
19 }
20 }
21 } else if (param > 1) {
22 Vector2 dir = direction(1).normalize();
23 Vector2 bq = origin-point(1);
24 double ts = dotProduct(bq, dir);
25 if (ts > 0) {
26 double pseudoDistance = crossProduct(bq, dir);
27 if (fabs(pseudoDistance) <= fabs(distance.distance)) {
28 distance.distance = pseudoDistance;
29 distance.dot = 0;
30 }
31 }
32 }
33}
34
35LinearSegment::LinearSegment(Point2 p0, Point2 p1, EdgeColor edgeColor) : EdgeSegment(edgeColor) {
36 p[0] = p0;
37 p[1] = p1;
38}
39
40QuadraticSegment::QuadraticSegment(Point2 p0, Point2 p1, Point2 p2, EdgeColor edgeColor) : EdgeSegment(edgeColor) {
41 if (p1 == p0 || p1 == p2)
42 p1 = 0.5*(p0+p2);
43 p[0] = p0;
44 p[1] = p1;
45 p[2] = p2;
46}
47
48CubicSegment::CubicSegment(Point2 p0, Point2 p1, Point2 p2, Point2 p3, EdgeColor edgeColor) : EdgeSegment(edgeColor) {
49 if ((p1 == p0 || p1 == p3) && (p2 == p0 || p2 == p3)) {
50 p1 = mix(p0, p3, 1/3.);
51 p2 = mix(p0, p3, 2/3.);
52 }
53 p[0] = p0;
54 p[1] = p1;
55 p[2] = p2;
56 p[3] = p3;
57}
58
59LinearSegment * LinearSegment::clone() const {
60 return new LinearSegment(p[0], p[1], color);
61}
62
63QuadraticSegment * QuadraticSegment::clone() const {
64 return new QuadraticSegment(p[0], p[1], p[2], color);
65}
66
67CubicSegment * CubicSegment::clone() const {
68 return new CubicSegment(p[0], p[1], p[2], p[3], color);
69}
70
71Point2 LinearSegment::point(double param) const {
72 return mix(p[0], p[1], param);
73}
74
75Point2 QuadraticSegment::point(double param) const {
76 return mix(mix(p[0], p[1], param), mix(p[1], p[2], param), param);
77}
78
79Point2 CubicSegment::point(double param) const {
80 Vector2 p12 = mix(p[1], p[2], param);
81 return mix(mix(mix(p[0], p[1], param), p12, param), mix(p12, mix(p[2], p[3], param), param), param);
82}
83
84Vector2 LinearSegment::direction(double param) const {
85 return p[1]-p[0];
86}
87
88Vector2 QuadraticSegment::direction(double param) const {
89 Vector2 tangent = mix(p[1]-p[0], p[2]-p[1], param);
90 if (!tangent)
91 return p[2]-p[0];
92 return tangent;
93}
94
95Vector2 CubicSegment::direction(double param) const {
96 Vector2 tangent = mix(mix(p[1]-p[0], p[2]-p[1], param), mix(p[2]-p[1], p[3]-p[2], param), param);
97 if (!tangent) {
98 if (param == 0) return p[2]-p[0];
99 if (param == 1) return p[3]-p[1];
100 }
101 return tangent;
102}
103
104Vector2 LinearSegment::directionChange(double param) const {
105 return Vector2();
106}
107
108Vector2 QuadraticSegment::directionChange(double param) const {
109 return (p[2]-p[1])-(p[1]-p[0]);
110}
111
112Vector2 CubicSegment::directionChange(double param) const {
113 return mix((p[2]-p[1])-(p[1]-p[0]), (p[3]-p[2])-(p[2]-p[1]), param);
114}
115
116double LinearSegment::length() const {
117 return (p[1]-p[0]).length();
118}
119
120double QuadraticSegment::length() const {
121 Vector2 ab = p[1]-p[0];
122 Vector2 br = p[2]-p[1]-ab;
123 double abab = dotProduct(ab, ab);
124 double abbr = dotProduct(ab, br);
125 double brbr = dotProduct(br, br);
126 double abLen = sqrt(abab);
127 double brLen = sqrt(brbr);
128 double crs = crossProduct(ab, br);
129 double h = sqrt(abab+abbr+abbr+brbr);
130 return (
131 brLen*((abbr+brbr)*h-abbr*abLen)+
132 crs*crs*log((brLen*h+abbr+brbr)/(brLen*abLen+abbr))
133 )/(brbr*brLen);
134}
135
136SignedDistance LinearSegment::signedDistance(Point2 origin, double &param) const {
137 Vector2 aq = origin-p[0];
138 Vector2 ab = p[1]-p[0];
139 param = dotProduct(aq, ab)/dotProduct(ab, ab);
140 Vector2 eq = p[param > .5]-origin;
141 double endpointDistance = eq.length();
142 if (param > 0 && param < 1) {
143 double orthoDistance = dotProduct(ab.getOrthonormal(false), aq);
144 if (fabs(orthoDistance) < endpointDistance)
145 return SignedDistance(orthoDistance, 0);
146 }
147 return SignedDistance(nonZeroSign(crossProduct(aq, ab))*endpointDistance, fabs(dotProduct(ab.normalize(), eq.normalize())));
148}
149
150SignedDistance QuadraticSegment::signedDistance(Point2 origin, double &param) const {
151 Vector2 qa = p[0]-origin;
152 Vector2 ab = p[1]-p[0];
153 Vector2 br = p[2]-p[1]-ab;
154 double a = dotProduct(br, br);
155 double b = 3*dotProduct(ab, br);
156 double c = 2*dotProduct(ab, ab)+dotProduct(qa, br);
157 double d = dotProduct(qa, ab);
158 double t[3];
159 int solutions = solveCubic(t, a, b, c, d);
160
161 Vector2 epDir = direction(0);
162 double minDistance = nonZeroSign(crossProduct(epDir, qa))*qa.length(); // distance from A
163 param = -dotProduct(qa, epDir)/dotProduct(epDir, epDir);
164 {
165 epDir = direction(1);
166 double distance = (p[2]-origin).length(); // distance from B
167 if (distance < fabs(minDistance)) {
168 minDistance = nonZeroSign(crossProduct(epDir, p[2]-origin))*distance;
169 param = dotProduct(origin-p[1], epDir)/dotProduct(epDir, epDir);
170 }
171 }
172 for (int i = 0; i < solutions; ++i) {
173 if (t[i] > 0 && t[i] < 1) {
174 Point2 qe = qa+2*t[i]*ab+t[i]*t[i]*br;
175 double distance = qe.length();
176 if (distance <= fabs(minDistance)) {
177 minDistance = nonZeroSign(crossProduct(ab+t[i]*br, qe))*distance;
178 param = t[i];
179 }
180 }
181 }
182
183 if (param >= 0 && param <= 1)
184 return SignedDistance(minDistance, 0);
185 if (param < .5)
186 return SignedDistance(minDistance, fabs(dotProduct(direction(0).normalize(), qa.normalize())));
187 else
188 return SignedDistance(minDistance, fabs(dotProduct(direction(1).normalize(), (p[2]-origin).normalize())));
189}
190
191SignedDistance CubicSegment::signedDistance(Point2 origin, double &param) const {
192 Vector2 qa = p[0]-origin;
193 Vector2 ab = p[1]-p[0];
194 Vector2 br = p[2]-p[1]-ab;
195 Vector2 as = (p[3]-p[2])-(p[2]-p[1])-br;
196
197 Vector2 epDir = direction(0);
198 double minDistance = nonZeroSign(crossProduct(epDir, qa))*qa.length(); // distance from A
199 param = -dotProduct(qa, epDir)/dotProduct(epDir, epDir);
200 {
201 epDir = direction(1);
202 double distance = (p[3]-origin).length(); // distance from B
203 if (distance < fabs(minDistance)) {
204 minDistance = nonZeroSign(crossProduct(epDir, p[3]-origin))*distance;
205 param = dotProduct(epDir-(p[3]-origin), epDir)/dotProduct(epDir, epDir);
206 }
207 }
208 // Iterative minimum distance search
209 for (int i = 0; i <= MSDFGEN_CUBIC_SEARCH_STARTS; ++i) {
210 double t = (double) i/MSDFGEN_CUBIC_SEARCH_STARTS;
211 Vector2 qe = qa+3*t*ab+3*t*t*br+t*t*t*as;
212 for (int step = 0; step < MSDFGEN_CUBIC_SEARCH_STEPS; ++step) {
213 // Improve t
214 Vector2 d1 = 3*ab+6*t*br+3*t*t*as;
215 Vector2 d2 = 6*br+6*t*as;
216 t -= dotProduct(qe, d1)/(dotProduct(d1, d1)+dotProduct(qe, d2));
217 if (t <= 0 || t >= 1)
218 break;
219 qe = qa+3*t*ab+3*t*t*br+t*t*t*as;
220 double distance = qe.length();
221 if (distance < fabs(minDistance)) {
222 minDistance = nonZeroSign(crossProduct(d1, qe))*distance;
223 param = t;
224 }
225 }
226 }
227
228 if (param >= 0 && param <= 1)
229 return SignedDistance(minDistance, 0);
230 if (param < .5)
231 return SignedDistance(minDistance, fabs(dotProduct(direction(0).normalize(), qa.normalize())));
232 else
233 return SignedDistance(minDistance, fabs(dotProduct(direction(1).normalize(), (p[3]-origin).normalize())));
234}
235
236int LinearSegment::scanlineIntersections(double x[3], int dy[3], double y) const {
237 if ((y >= p[0].y && y < p[1].y) || (y >= p[1].y && y < p[0].y)) {
238 double param = (y-p[0].y)/(p[1].y-p[0].y);
239 x[0] = mix(p[0].x, p[1].x, param);
240 dy[0] = sign(p[1].y-p[0].y);
241 return 1;
242 }
243 return 0;
244}
245
246int QuadraticSegment::scanlineIntersections(double x[3], int dy[3], double y) const {
247 int total = 0;
248 int nextDY = y > p[0].y ? 1 : -1;
249 x[total] = p[0].x;
250 if (p[0].y == y) {
251 if (p[0].y < p[1].y || (p[0].y == p[1].y && p[0].y < p[2].y))
252 dy[total++] = 1;
253 else
254 nextDY = 1;
255 }
256 {
257 Vector2 ab = p[1]-p[0];
258 Vector2 br = p[2]-p[1]-ab;
259 double t[2];
260 int solutions = solveQuadratic(t, br.y, 2*ab.y, p[0].y-y);
261 // Sort solutions
262 double tmp;
263 if (solutions >= 2 && t[0] > t[1])
264 tmp = t[0], t[0] = t[1], t[1] = tmp;
265 for (int i = 0; i < solutions && total < 2; ++i) {
266 if (t[i] >= 0 && t[i] <= 1) {
267 x[total] = p[0].x+2*t[i]*ab.x+t[i]*t[i]*br.x;
268 if (nextDY*(ab.y+t[i]*br.y) >= 0) {
269 dy[total++] = nextDY;
270 nextDY = -nextDY;
271 }
272 }
273 }
274 }
275 if (p[2].y == y) {
276 if (nextDY > 0 && total > 0) {
277 --total;
278 nextDY = -1;
279 }
280 if ((p[2].y < p[1].y || (p[2].y == p[1].y && p[2].y < p[0].y)) && total < 2) {
281 x[total] = p[2].x;
282 if (nextDY < 0) {
283 dy[total++] = -1;
284 nextDY = 1;
285 }
286 }
287 }
288 if (nextDY != (y >= p[2].y ? 1 : -1)) {
289 if (total > 0)
290 --total;
291 else {
292 if (fabs(p[2].y-y) < fabs(p[0].y-y))
293 x[total] = p[2].x;
294 dy[total++] = nextDY;
295 }
296 }
297 return total;
298}
299
300int CubicSegment::scanlineIntersections(double x[3], int dy[3], double y) const {
301 int total = 0;
302 int nextDY = y > p[0].y ? 1 : -1;
303 x[total] = p[0].x;
304 if (p[0].y == y) {
305 if (p[0].y < p[1].y || (p[0].y == p[1].y && (p[0].y < p[2].y || (p[0].y == p[2].y && p[0].y < p[3].y))))
306 dy[total++] = 1;
307 else
308 nextDY = 1;
309 }
310 {
311 Vector2 ab = p[1]-p[0];
312 Vector2 br = p[2]-p[1]-ab;
313 Vector2 as = (p[3]-p[2])-(p[2]-p[1])-br;
314 double t[3];
315 int solutions = solveCubic(t, as.y, 3*br.y, 3*ab.y, p[0].y-y);
316 // Sort solutions
317 double tmp;
318 if (solutions >= 2) {
319 if (t[0] > t[1])
320 tmp = t[0], t[0] = t[1], t[1] = tmp;
321 if (solutions >= 3 && t[1] > t[2]) {
322 tmp = t[1], t[1] = t[2], t[2] = tmp;
323 if (t[0] > t[1])
324 tmp = t[0], t[0] = t[1], t[1] = tmp;
325 }
326 }
327 for (int i = 0; i < solutions && total < 3; ++i) {
328 if (t[i] >= 0 && t[i] <= 1) {
329 x[total] = p[0].x+3*t[i]*ab.x+3*t[i]*t[i]*br.x+t[i]*t[i]*t[i]*as.x;
330 if (nextDY*(ab.y+2*t[i]*br.y+t[i]*t[i]*as.y) >= 0) {
331 dy[total++] = nextDY;
332 nextDY = -nextDY;
333 }
334 }
335 }
336 }
337 if (p[3].y == y) {
338 if (nextDY > 0 && total > 0) {
339 --total;
340 nextDY = -1;
341 }
342 if ((p[3].y < p[2].y || (p[3].y == p[2].y && (p[3].y < p[1].y || (p[3].y == p[1].y && p[3].y < p[0].y)))) && total < 3) {
343 x[total] = p[3].x;
344 if (nextDY < 0) {
345 dy[total++] = -1;
346 nextDY = 1;
347 }
348 }
349 }
350 if (nextDY != (y >= p[3].y ? 1 : -1)) {
351 if (total > 0)
352 --total;
353 else {
354 if (fabs(p[3].y-y) < fabs(p[0].y-y))
355 x[total] = p[3].x;
356 dy[total++] = nextDY;
357 }
358 }
359 return total;
360}
361
362static void pointBounds(Point2 p, double &l, double &b, double &r, double &t) {
363 if (p.x < l) l = p.x;
364 if (p.y < b) b = p.y;
365 if (p.x > r) r = p.x;
366 if (p.y > t) t = p.y;
367}
368
369void LinearSegment::bound(double &l, double &b, double &r, double &t) const {
370 pointBounds(p[0], l, b, r, t);
371 pointBounds(p[1], l, b, r, t);
372}
373
374void QuadraticSegment::bound(double &l, double &b, double &r, double &t) const {
375 pointBounds(p[0], l, b, r, t);
376 pointBounds(p[2], l, b, r, t);
377 Vector2 bot = (p[1]-p[0])-(p[2]-p[1]);
378 if (bot.x) {
379 double param = (p[1].x-p[0].x)/bot.x;
380 if (param > 0 && param < 1)
381 pointBounds(point(param), l, b, r, t);
382 }
383 if (bot.y) {
384 double param = (p[1].y-p[0].y)/bot.y;
385 if (param > 0 && param < 1)
386 pointBounds(point(param), l, b, r, t);
387 }
388}
389
390void CubicSegment::bound(double &l, double &b, double &r, double &t) const {
391 pointBounds(p[0], l, b, r, t);
392 pointBounds(p[3], l, b, r, t);
393 Vector2 a0 = p[1]-p[0];
394 Vector2 a1 = 2*(p[2]-p[1]-a0);
395 Vector2 a2 = p[3]-3*p[2]+3*p[1]-p[0];
396 double params[2];
397 int solutions;
398 solutions = solveQuadratic(params, a2.x, a1.x, a0.x);
399 for (int i = 0; i < solutions; ++i)
400 if (params[i] > 0 && params[i] < 1)
401 pointBounds(point(params[i]), l, b, r, t);
402 solutions = solveQuadratic(params, a2.y, a1.y, a0.y);
403 for (int i = 0; i < solutions; ++i)
404 if (params[i] > 0 && params[i] < 1)
405 pointBounds(point(params[i]), l, b, r, t);
406}
407
408void LinearSegment::reverse() {
409 Point2 tmp = p[0];
410 p[0] = p[1];
411 p[1] = tmp;
412}
413
414void QuadraticSegment::reverse() {
415 Point2 tmp = p[0];
416 p[0] = p[2];
417 p[2] = tmp;
418}
419
420void CubicSegment::reverse() {
421 Point2 tmp = p[0];
422 p[0] = p[3];
423 p[3] = tmp;
424 tmp = p[1];
425 p[1] = p[2];
426 p[2] = tmp;
427}
428
429void LinearSegment::moveStartPoint(Point2 to) {
430 p[0] = to;
431}
432
433void QuadraticSegment::moveStartPoint(Point2 to) {
434 Vector2 origSDir = p[0]-p[1];
435 Point2 origP1 = p[1];
436 p[1] += crossProduct(p[0]-p[1], to-p[0])/crossProduct(p[0]-p[1], p[2]-p[1])*(p[2]-p[1]);
437 p[0] = to;
438 if (dotProduct(origSDir, p[0]-p[1]) < 0)
439 p[1] = origP1;
440}
441
442void CubicSegment::moveStartPoint(Point2 to) {
443 p[1] += to-p[0];
444 p[0] = to;
445}
446
447void LinearSegment::moveEndPoint(Point2 to) {
448 p[1] = to;
449}
450
451void QuadraticSegment::moveEndPoint(Point2 to) {
452 Vector2 origEDir = p[2]-p[1];
453 Point2 origP1 = p[1];
454 p[1] += crossProduct(p[2]-p[1], to-p[2])/crossProduct(p[2]-p[1], p[0]-p[1])*(p[0]-p[1]);
455 p[2] = to;
456 if (dotProduct(origEDir, p[2]-p[1]) < 0)
457 p[1] = origP1;
458}
459
460void CubicSegment::moveEndPoint(Point2 to) {
461 p[2] += to-p[3];
462 p[3] = to;
463}
464
465void LinearSegment::splitInThirds(EdgeSegment *&part1, EdgeSegment *&part2, EdgeSegment *&part3) const {
466 part1 = new LinearSegment(p[0], point(1/3.), color);
467 part2 = new LinearSegment(point(1/3.), point(2/3.), color);
468 part3 = new LinearSegment(point(2/3.), p[1], color);
469}
470
471void QuadraticSegment::splitInThirds(EdgeSegment *&part1, EdgeSegment *&part2, EdgeSegment *&part3) const {
472 part1 = new QuadraticSegment(p[0], mix(p[0], p[1], 1/3.), point(1/3.), color);
473 part2 = new QuadraticSegment(point(1/3.), mix(mix(p[0], p[1], 5/9.), mix(p[1], p[2], 4/9.), .5), point(2/3.), color);
474 part3 = new QuadraticSegment(point(2/3.), mix(p[1], p[2], 2/3.), p[2], color);
475}
476
477void CubicSegment::splitInThirds(EdgeSegment *&part1, EdgeSegment *&part2, EdgeSegment *&part3) const {
478 part1 = new CubicSegment(p[0], p[0] == p[1] ? p[0] : mix(p[0], p[1], 1/3.), mix(mix(p[0], p[1], 1/3.), mix(p[1], p[2], 1/3.), 1/3.), point(1/3.), color);
479 part2 = new CubicSegment(point(1/3.),
480 mix(mix(mix(p[0], p[1], 1/3.), mix(p[1], p[2], 1/3.), 1/3.), mix(mix(p[1], p[2], 1/3.), mix(p[2], p[3], 1/3.), 1/3.), 2/3.),
481 mix(mix(mix(p[0], p[1], 2/3.), mix(p[1], p[2], 2/3.), 2/3.), mix(mix(p[1], p[2], 2/3.), mix(p[2], p[3], 2/3.), 2/3.), 1/3.),
482 point(2/3.), color);
483 part3 = new CubicSegment(point(2/3.), mix(mix(p[1], p[2], 2/3.), mix(p[2], p[3], 2/3.), 2/3.), p[2] == p[3] ? p[3] : mix(p[2], p[3], 2/3.), p[3], color);
484}
485
486EdgeSegment * QuadraticSegment::convertToCubic() const {
487 return new CubicSegment(p[0], mix(p[0], p[1], 2/3.), mix(p[1], p[2], 1/3.), p[2], color);
488}
489
490void CubicSegment::deconverge(int param, double amount) {
491 Vector2 dir = direction(param);
492 Vector2 normal = dir.getOrthonormal();
493 double h = dotProduct(directionChange(param)-dir, normal);
494 switch (param) {
495 case 0:
496 p[1] += amount*(dir+sign(h)*sqrt(fabs(h))*normal);
497 break;
498 case 1:
499 p[2] -= amount*(dir-sign(h)*sqrt(fabs(h))*normal);
500 break;
501 }
502}
503
504}
505