1#cython: language_level=3 2#distutils: define_macros=CYTHON_TRACE_NOGIL=1 3 4# Copyright 2015 Google Inc. All Rights Reserved. 5# 6# Licensed under the Apache License, Version 2.0 (the "License"); 7# you may not use this file except in compliance with the License. 8# You may obtain a copy of the License at 9# 10# http://www.apache.org/licenses/LICENSE-2.0 11# 12# Unless required by applicable law or agreed to in writing, software 13# distributed under the License is distributed on an "AS IS" BASIS, 14# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 15# See the License for the specific language governing permissions and 16# limitations under the License. 17 18try: 19 import cython 20except ImportError: 21 # if cython not installed, use mock module with no-op decorators and types 22 from fontTools.misc import cython 23 24import math 25 26from .errors import Error as Cu2QuError, ApproxNotFoundError 27 28 29__all__ = ['curve_to_quadratic', 'curves_to_quadratic'] 30 31MAX_N = 100 32 33NAN = float("NaN") 34 35 36if cython.compiled: 37 # Yep, I'm compiled. 38 COMPILED = True 39else: 40 # Just a lowly interpreted script. 41 COMPILED = False 42 43 44@cython.cfunc 45@cython.inline 46@cython.returns(cython.double) 47@cython.locals(v1=cython.complex, v2=cython.complex) 48def dot(v1, v2): 49 """Return the dot product of two vectors. 50 51 Args: 52 v1 (complex): First vector. 53 v2 (complex): Second vector. 54 55 Returns: 56 double: Dot product. 57 """ 58 return (v1 * v2.conjugate()).real 59 60 61@cython.cfunc 62@cython.inline 63@cython.locals(a=cython.complex, b=cython.complex, c=cython.complex, d=cython.complex) 64@cython.locals(_1=cython.complex, _2=cython.complex, _3=cython.complex, _4=cython.complex) 65def calc_cubic_points(a, b, c, d): 66 _1 = d 67 _2 = (c / 3.0) + d 68 _3 = (b + c) / 3.0 + _2 69 _4 = a + d + c + b 70 return _1, _2, _3, _4 71 72 73@cython.cfunc 74@cython.inline 75@cython.locals(p0=cython.complex, p1=cython.complex, p2=cython.complex, p3=cython.complex) 76@cython.locals(a=cython.complex, b=cython.complex, c=cython.complex, d=cython.complex) 77def calc_cubic_parameters(p0, p1, p2, p3): 78 c = (p1 - p0) * 3.0 79 b = (p2 - p1) * 3.0 - c 80 d = p0 81 a = p3 - d - c - b 82 return a, b, c, d 83 84 85@cython.cfunc 86@cython.locals(p0=cython.complex, p1=cython.complex, p2=cython.complex, p3=cython.complex) 87def split_cubic_into_n_iter(p0, p1, p2, p3, n): 88 """Split a cubic Bezier into n equal parts. 89 90 Splits the curve into `n` equal parts by curve time. 91 (t=0..1/n, t=1/n..2/n, ...) 92 93 Args: 94 p0 (complex): Start point of curve. 95 p1 (complex): First handle of curve. 96 p2 (complex): Second handle of curve. 97 p3 (complex): End point of curve. 98 99 Returns: 100 An iterator yielding the control points (four complex values) of the 101 subcurves. 102 """ 103 # Hand-coded special-cases 104 if n == 2: 105 return iter(split_cubic_into_two(p0, p1, p2, p3)) 106 if n == 3: 107 return iter(split_cubic_into_three(p0, p1, p2, p3)) 108 if n == 4: 109 a, b = split_cubic_into_two(p0, p1, p2, p3) 110 return iter(split_cubic_into_two(*a) + split_cubic_into_two(*b)) 111 if n == 6: 112 a, b = split_cubic_into_two(p0, p1, p2, p3) 113 return iter(split_cubic_into_three(*a) + split_cubic_into_three(*b)) 114 115 return _split_cubic_into_n_gen(p0,p1,p2,p3,n) 116 117 118@cython.locals(p0=cython.complex, p1=cython.complex, p2=cython.complex, p3=cython.complex, n=cython.int) 119@cython.locals(a=cython.complex, b=cython.complex, c=cython.complex, d=cython.complex) 120@cython.locals(dt=cython.double, delta_2=cython.double, delta_3=cython.double, i=cython.int) 121@cython.locals(a1=cython.complex, b1=cython.complex, c1=cython.complex, d1=cython.complex) 122def _split_cubic_into_n_gen(p0, p1, p2, p3, n): 123 a, b, c, d = calc_cubic_parameters(p0, p1, p2, p3) 124 dt = 1 / n 125 delta_2 = dt * dt 126 delta_3 = dt * delta_2 127 for i in range(n): 128 t1 = i * dt 129 t1_2 = t1 * t1 130 # calc new a, b, c and d 131 a1 = a * delta_3 132 b1 = (3*a*t1 + b) * delta_2 133 c1 = (2*b*t1 + c + 3*a*t1_2) * dt 134 d1 = a*t1*t1_2 + b*t1_2 + c*t1 + d 135 yield calc_cubic_points(a1, b1, c1, d1) 136 137 138@cython.locals(p0=cython.complex, p1=cython.complex, p2=cython.complex, p3=cython.complex) 139@cython.locals(mid=cython.complex, deriv3=cython.complex) 140def split_cubic_into_two(p0, p1, p2, p3): 141 """Split a cubic Bezier into two equal parts. 142 143 Splits the curve into two equal parts at t = 0.5 144 145 Args: 146 p0 (complex): Start point of curve. 147 p1 (complex): First handle of curve. 148 p2 (complex): Second handle of curve. 149 p3 (complex): End point of curve. 150 151 Returns: 152 tuple: Two cubic Beziers (each expressed as a tuple of four complex 153 values). 154 """ 155 mid = (p0 + 3 * (p1 + p2) + p3) * .125 156 deriv3 = (p3 + p2 - p1 - p0) * .125 157 return ((p0, (p0 + p1) * .5, mid - deriv3, mid), 158 (mid, mid + deriv3, (p2 + p3) * .5, p3)) 159 160 161@cython.locals(p0=cython.complex, p1=cython.complex, p2=cython.complex, p3=cython.complex, _27=cython.double) 162@cython.locals(mid1=cython.complex, deriv1=cython.complex, mid2=cython.complex, deriv2=cython.complex) 163def split_cubic_into_three(p0, p1, p2, p3, _27=1/27): 164 """Split a cubic Bezier into three equal parts. 165 166 Splits the curve into three equal parts at t = 1/3 and t = 2/3 167 168 Args: 169 p0 (complex): Start point of curve. 170 p1 (complex): First handle of curve. 171 p2 (complex): Second handle of curve. 172 p3 (complex): End point of curve. 173 174 Returns: 175 tuple: Three cubic Beziers (each expressed as a tuple of four complex 176 values). 177 """ 178 # we define 1/27 as a keyword argument so that it will be evaluated only 179 # once but still in the scope of this function 180 mid1 = (8*p0 + 12*p1 + 6*p2 + p3) * _27 181 deriv1 = (p3 + 3*p2 - 4*p0) * _27 182 mid2 = (p0 + 6*p1 + 12*p2 + 8*p3) * _27 183 deriv2 = (4*p3 - 3*p1 - p0) * _27 184 return ((p0, (2*p0 + p1) / 3.0, mid1 - deriv1, mid1), 185 (mid1, mid1 + deriv1, mid2 - deriv2, mid2), 186 (mid2, mid2 + deriv2, (p2 + 2*p3) / 3.0, p3)) 187 188 189@cython.returns(cython.complex) 190@cython.locals(t=cython.double, p0=cython.complex, p1=cython.complex, p2=cython.complex, p3=cython.complex) 191@cython.locals(_p1=cython.complex, _p2=cython.complex) 192def cubic_approx_control(t, p0, p1, p2, p3): 193 """Approximate a cubic Bezier using a quadratic one. 194 195 Args: 196 t (double): Position of control point. 197 p0 (complex): Start point of curve. 198 p1 (complex): First handle of curve. 199 p2 (complex): Second handle of curve. 200 p3 (complex): End point of curve. 201 202 Returns: 203 complex: Location of candidate control point on quadratic curve. 204 """ 205 _p1 = p0 + (p1 - p0) * 1.5 206 _p2 = p3 + (p2 - p3) * 1.5 207 return _p1 + (_p2 - _p1) * t 208 209 210@cython.returns(cython.complex) 211@cython.locals(a=cython.complex, b=cython.complex, c=cython.complex, d=cython.complex) 212@cython.locals(ab=cython.complex, cd=cython.complex, p=cython.complex, h=cython.double) 213def calc_intersect(a, b, c, d): 214 """Calculate the intersection of two lines. 215 216 Args: 217 a (complex): Start point of first line. 218 b (complex): End point of first line. 219 c (complex): Start point of second line. 220 d (complex): End point of second line. 221 222 Returns: 223 complex: Location of intersection if one present, ``complex(NaN,NaN)`` 224 if no intersection was found. 225 """ 226 ab = b - a 227 cd = d - c 228 p = ab * 1j 229 try: 230 h = dot(p, a - c) / dot(p, cd) 231 except ZeroDivisionError: 232 return complex(NAN, NAN) 233 return c + cd * h 234 235 236@cython.cfunc 237@cython.returns(cython.int) 238@cython.locals(tolerance=cython.double, p0=cython.complex, p1=cython.complex, p2=cython.complex, p3=cython.complex) 239@cython.locals(mid=cython.complex, deriv3=cython.complex) 240def cubic_farthest_fit_inside(p0, p1, p2, p3, tolerance): 241 """Check if a cubic Bezier lies within a given distance of the origin. 242 243 "Origin" means *the* origin (0,0), not the start of the curve. Note that no 244 checks are made on the start and end positions of the curve; this function 245 only checks the inside of the curve. 246 247 Args: 248 p0 (complex): Start point of curve. 249 p1 (complex): First handle of curve. 250 p2 (complex): Second handle of curve. 251 p3 (complex): End point of curve. 252 tolerance (double): Distance from origin. 253 254 Returns: 255 bool: True if the cubic Bezier ``p`` entirely lies within a distance 256 ``tolerance`` of the origin, False otherwise. 257 """ 258 # First check p2 then p1, as p2 has higher error early on. 259 if abs(p2) <= tolerance and abs(p1) <= tolerance: 260 return True 261 262 # Split. 263 mid = (p0 + 3 * (p1 + p2) + p3) * .125 264 if abs(mid) > tolerance: 265 return False 266 deriv3 = (p3 + p2 - p1 - p0) * .125 267 return (cubic_farthest_fit_inside(p0, (p0+p1)*.5, mid-deriv3, mid, tolerance) and 268 cubic_farthest_fit_inside(mid, mid+deriv3, (p2+p3)*.5, p3, tolerance)) 269 270 271@cython.cfunc 272@cython.locals(tolerance=cython.double, _2_3=cython.double) 273@cython.locals(q1=cython.complex, c0=cython.complex, c1=cython.complex, c2=cython.complex, c3=cython.complex) 274def cubic_approx_quadratic(cubic, tolerance, _2_3=2/3): 275 """Approximate a cubic Bezier with a single quadratic within a given tolerance. 276 277 Args: 278 cubic (sequence): Four complex numbers representing control points of 279 the cubic Bezier curve. 280 tolerance (double): Permitted deviation from the original curve. 281 282 Returns: 283 Three complex numbers representing control points of the quadratic 284 curve if it fits within the given tolerance, or ``None`` if no suitable 285 curve could be calculated. 286 """ 287 # we define 2/3 as a keyword argument so that it will be evaluated only 288 # once but still in the scope of this function 289 290 q1 = calc_intersect(*cubic) 291 if math.isnan(q1.imag): 292 return None 293 c0 = cubic[0] 294 c3 = cubic[3] 295 c1 = c0 + (q1 - c0) * _2_3 296 c2 = c3 + (q1 - c3) * _2_3 297 if not cubic_farthest_fit_inside(0, 298 c1 - cubic[1], 299 c2 - cubic[2], 300 0, tolerance): 301 return None 302 return c0, q1, c3 303 304 305@cython.cfunc 306@cython.locals(n=cython.int, tolerance=cython.double, _2_3=cython.double) 307@cython.locals(i=cython.int) 308@cython.locals(c0=cython.complex, c1=cython.complex, c2=cython.complex, c3=cython.complex) 309@cython.locals(q0=cython.complex, q1=cython.complex, next_q1=cython.complex, q2=cython.complex, d1=cython.complex) 310def cubic_approx_spline(cubic, n, tolerance, _2_3=2/3): 311 """Approximate a cubic Bezier curve with a spline of n quadratics. 312 313 Args: 314 cubic (sequence): Four complex numbers representing control points of 315 the cubic Bezier curve. 316 n (int): Number of quadratic Bezier curves in the spline. 317 tolerance (double): Permitted deviation from the original curve. 318 319 Returns: 320 A list of ``n+2`` complex numbers, representing control points of the 321 quadratic spline if it fits within the given tolerance, or ``None`` if 322 no suitable spline could be calculated. 323 """ 324 # we define 2/3 as a keyword argument so that it will be evaluated only 325 # once but still in the scope of this function 326 327 if n == 1: 328 return cubic_approx_quadratic(cubic, tolerance) 329 330 cubics = split_cubic_into_n_iter(cubic[0], cubic[1], cubic[2], cubic[3], n) 331 332 # calculate the spline of quadratics and check errors at the same time. 333 next_cubic = next(cubics) 334 next_q1 = cubic_approx_control(0, *next_cubic) 335 q2 = cubic[0] 336 d1 = 0j 337 spline = [cubic[0], next_q1] 338 for i in range(1, n+1): 339 340 # Current cubic to convert 341 c0, c1, c2, c3 = next_cubic 342 343 # Current quadratic approximation of current cubic 344 q0 = q2 345 q1 = next_q1 346 if i < n: 347 next_cubic = next(cubics) 348 next_q1 = cubic_approx_control(i / (n-1), *next_cubic) 349 spline.append(next_q1) 350 q2 = (q1 + next_q1) * .5 351 else: 352 q2 = c3 353 354 # End-point deltas 355 d0 = d1 356 d1 = q2 - c3 357 358 if (abs(d1) > tolerance or 359 not cubic_farthest_fit_inside(d0, 360 q0 + (q1 - q0) * _2_3 - c1, 361 q2 + (q1 - q2) * _2_3 - c2, 362 d1, 363 tolerance)): 364 return None 365 spline.append(cubic[3]) 366 367 return spline 368 369 370@cython.locals(max_err=cython.double) 371@cython.locals(n=cython.int) 372def curve_to_quadratic(curve, max_err): 373 """Approximate a cubic Bezier curve with a spline of n quadratics. 374 375 Args: 376 cubic (sequence): Four 2D tuples representing control points of 377 the cubic Bezier curve. 378 max_err (double): Permitted deviation from the original curve. 379 380 Returns: 381 A list of 2D tuples, representing control points of the quadratic 382 spline if it fits within the given tolerance, or ``None`` if no 383 suitable spline could be calculated. 384 """ 385 386 curve = [complex(*p) for p in curve] 387 388 for n in range(1, MAX_N + 1): 389 spline = cubic_approx_spline(curve, n, max_err) 390 if spline is not None: 391 # done. go home 392 return [(s.real, s.imag) for s in spline] 393 394 raise ApproxNotFoundError(curve) 395 396 397 398@cython.locals(l=cython.int, last_i=cython.int, i=cython.int) 399def curves_to_quadratic(curves, max_errors): 400 """Return quadratic Bezier splines approximating the input cubic Beziers. 401 402 Args: 403 curves: A sequence of *n* curves, each curve being a sequence of four 404 2D tuples. 405 max_errors: A sequence of *n* floats representing the maximum permissible 406 deviation from each of the cubic Bezier curves. 407 408 Example:: 409 410 >>> curves_to_quadratic( [ 411 ... [ (50,50), (100,100), (150,100), (200,50) ], 412 ... [ (75,50), (120,100), (150,75), (200,60) ] 413 ... ], [1,1] ) 414 [[(50.0, 50.0), (75.0, 75.0), (125.0, 91.66666666666666), (175.0, 75.0), (200.0, 50.0)], [(75.0, 50.0), (97.5, 75.0), (135.41666666666666, 82.08333333333333), (175.0, 67.5), (200.0, 60.0)]] 415 416 The returned splines have "implied oncurve points" suitable for use in 417 TrueType ``glif`` outlines - i.e. in the first spline returned above, 418 the first quadratic segment runs from (50,50) to 419 ( (75 + 125)/2 , (120 + 91.666..)/2 ) = (100, 83.333...). 420 421 Returns: 422 A list of splines, each spline being a list of 2D tuples. 423 424 Raises: 425 fontTools.cu2qu.Errors.ApproxNotFoundError: if no suitable approximation 426 can be found for all curves with the given parameters. 427 """ 428 429 curves = [[complex(*p) for p in curve] for curve in curves] 430 assert len(max_errors) == len(curves) 431 432 l = len(curves) 433 splines = [None] * l 434 last_i = i = 0 435 n = 1 436 while True: 437 spline = cubic_approx_spline(curves[i], n, max_errors[i]) 438 if spline is None: 439 if n == MAX_N: 440 break 441 n += 1 442 last_i = i 443 continue 444 splines[i] = spline 445 i = (i + 1) % l 446 if i == last_i: 447 # done. go home 448 return [[(s.real, s.imag) for s in spline] for spline in splines] 449 450 raise ApproxNotFoundError(curves) 451 452 453if __name__ == '__main__': 454 import random 455 import timeit 456 457 MAX_ERR = 5 458 459 def generate_curve(): 460 return [ 461 tuple(float(random.randint(0, 2048)) for coord in range(2)) 462 for point in range(4)] 463 464 def setup_curve_to_quadratic(): 465 return generate_curve(), MAX_ERR 466 467 def setup_curves_to_quadratic(): 468 num_curves = 3 469 return ( 470 [generate_curve() for curve in range(num_curves)], 471 [MAX_ERR] * num_curves) 472 473 def run_benchmark( 474 benchmark_module, module, function, setup_suffix='', repeat=5, number=1000): 475 setup_func = 'setup_' + function 476 if setup_suffix: 477 print('%s with %s:' % (function, setup_suffix), end='') 478 setup_func += '_' + setup_suffix 479 else: 480 print('%s:' % function, end='') 481 482 def wrapper(function, setup_func): 483 function = globals()[function] 484 setup_func = globals()[setup_func] 485 def wrapped(): 486 return function(*setup_func()) 487 return wrapped 488 results = timeit.repeat(wrapper(function, setup_func), repeat=repeat, number=number) 489 print('\t%5.1fus' % (min(results) * 1000000. / number)) 490 491 def main(): 492 run_benchmark('cu2qu.benchmark', 'cu2qu', 'curve_to_quadratic') 493 run_benchmark('cu2qu.benchmark', 'cu2qu', 'curves_to_quadratic') 494 495 random.seed(1) 496 main() 497