qu2cu.py 12 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408
  1. # cython: language_level=3
  2. # distutils: define_macros=CYTHON_TRACE_NOGIL=1
  3. # Copyright 2023 Google Inc. All Rights Reserved.
  4. # Copyright 2023 Behdad Esfahbod. 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. try:
  18. import cython
  19. COMPILED = cython.compiled
  20. except (AttributeError, ImportError):
  21. # if cython not installed, use mock module with no-op decorators and types
  22. from fontTools.misc import cython
  23. COMPILED = False
  24. from fontTools.misc.bezierTools import splitCubicAtTC
  25. from collections import namedtuple
  26. import math
  27. from typing import (
  28. List,
  29. Tuple,
  30. Union,
  31. )
  32. __all__ = ["quadratic_to_curves"]
  33. # Copied from cu2qu
  34. @cython.cfunc
  35. @cython.returns(cython.int)
  36. @cython.locals(
  37. tolerance=cython.double,
  38. p0=cython.complex,
  39. p1=cython.complex,
  40. p2=cython.complex,
  41. p3=cython.complex,
  42. )
  43. @cython.locals(mid=cython.complex, deriv3=cython.complex)
  44. def cubic_farthest_fit_inside(p0, p1, p2, p3, tolerance):
  45. """Check if a cubic Bezier lies within a given distance of the origin.
  46. "Origin" means *the* origin (0,0), not the start of the curve. Note that no
  47. checks are made on the start and end positions of the curve; this function
  48. only checks the inside of the curve.
  49. Args:
  50. p0 (complex): Start point of curve.
  51. p1 (complex): First handle of curve.
  52. p2 (complex): Second handle of curve.
  53. p3 (complex): End point of curve.
  54. tolerance (double): Distance from origin.
  55. Returns:
  56. bool: True if the cubic Bezier ``p`` entirely lies within a distance
  57. ``tolerance`` of the origin, False otherwise.
  58. """
  59. # First check p2 then p1, as p2 has higher error early on.
  60. if abs(p2) <= tolerance and abs(p1) <= tolerance:
  61. return True
  62. # Split.
  63. mid = (p0 + 3 * (p1 + p2) + p3) * 0.125
  64. if abs(mid) > tolerance:
  65. return False
  66. deriv3 = (p3 + p2 - p1 - p0) * 0.125
  67. return cubic_farthest_fit_inside(
  68. p0, (p0 + p1) * 0.5, mid - deriv3, mid, tolerance
  69. ) and cubic_farthest_fit_inside(mid, mid + deriv3, (p2 + p3) * 0.5, p3, tolerance)
  70. @cython.locals(
  71. p0=cython.complex,
  72. p1=cython.complex,
  73. p2=cython.complex,
  74. p1_2_3=cython.complex,
  75. )
  76. def elevate_quadratic(p0, p1, p2):
  77. """Given a quadratic bezier curve, return its degree-elevated cubic."""
  78. # https://pomax.github.io/bezierinfo/#reordering
  79. p1_2_3 = p1 * (2 / 3)
  80. return (
  81. p0,
  82. (p0 * (1 / 3) + p1_2_3),
  83. (p2 * (1 / 3) + p1_2_3),
  84. p2,
  85. )
  86. @cython.cfunc
  87. @cython.locals(
  88. start=cython.int,
  89. n=cython.int,
  90. k=cython.int,
  91. prod_ratio=cython.double,
  92. sum_ratio=cython.double,
  93. ratio=cython.double,
  94. t=cython.double,
  95. p0=cython.complex,
  96. p1=cython.complex,
  97. p2=cython.complex,
  98. p3=cython.complex,
  99. )
  100. def merge_curves(curves, start, n):
  101. """Give a cubic-Bezier spline, reconstruct one cubic-Bezier
  102. that has the same endpoints and tangents and approxmates
  103. the spline."""
  104. # Reconstruct the t values of the cut segments
  105. prod_ratio = 1.0
  106. sum_ratio = 1.0
  107. ts = [1]
  108. for k in range(1, n):
  109. ck = curves[start + k]
  110. c_before = curves[start + k - 1]
  111. # |t_(k+1) - t_k| / |t_k - t_(k - 1)| = ratio
  112. assert ck[0] == c_before[3]
  113. ratio = abs(ck[1] - ck[0]) / abs(c_before[3] - c_before[2])
  114. prod_ratio *= ratio
  115. sum_ratio += prod_ratio
  116. ts.append(sum_ratio)
  117. # (t(n) - t(n - 1)) / (t_(1) - t(0)) = prod_ratio
  118. ts = [t / sum_ratio for t in ts[:-1]]
  119. p0 = curves[start][0]
  120. p1 = curves[start][1]
  121. p2 = curves[start + n - 1][2]
  122. p3 = curves[start + n - 1][3]
  123. # Build the curve by scaling the control-points.
  124. p1 = p0 + (p1 - p0) / (ts[0] if ts else 1)
  125. p2 = p3 + (p2 - p3) / ((1 - ts[-1]) if ts else 1)
  126. curve = (p0, p1, p2, p3)
  127. return curve, ts
  128. @cython.locals(
  129. count=cython.int,
  130. num_offcurves=cython.int,
  131. i=cython.int,
  132. off1=cython.complex,
  133. off2=cython.complex,
  134. on=cython.complex,
  135. )
  136. def add_implicit_on_curves(p):
  137. q = list(p)
  138. count = 0
  139. num_offcurves = len(p) - 2
  140. for i in range(1, num_offcurves):
  141. off1 = p[i]
  142. off2 = p[i + 1]
  143. on = off1 + (off2 - off1) * 0.5
  144. q.insert(i + 1 + count, on)
  145. count += 1
  146. return q
  147. Point = Union[Tuple[float, float], complex]
  148. @cython.locals(
  149. cost=cython.int,
  150. is_complex=cython.int,
  151. )
  152. def quadratic_to_curves(
  153. quads: List[List[Point]],
  154. max_err: float = 0.5,
  155. all_cubic: bool = False,
  156. ) -> List[Tuple[Point, ...]]:
  157. """Converts a connecting list of quadratic splines to a list of quadratic
  158. and cubic curves.
  159. A quadratic spline is specified as a list of points. Either each point is
  160. a 2-tuple of X,Y coordinates, or each point is a complex number with
  161. real/imaginary components representing X,Y coordinates.
  162. The first and last points are on-curve points and the rest are off-curve
  163. points, with an implied on-curve point in the middle between every two
  164. consequtive off-curve points.
  165. Returns:
  166. The output is a list of tuples of points. Points are represented
  167. in the same format as the input, either as 2-tuples or complex numbers.
  168. Each tuple is either of length three, for a quadratic curve, or four,
  169. for a cubic curve. Each curve's last point is the same as the next
  170. curve's first point.
  171. Args:
  172. quads: quadratic splines
  173. max_err: absolute error tolerance; defaults to 0.5
  174. all_cubic: if True, only cubic curves are generated; defaults to False
  175. """
  176. is_complex = type(quads[0][0]) is complex
  177. if not is_complex:
  178. quads = [[complex(x, y) for (x, y) in p] for p in quads]
  179. q = [quads[0][0]]
  180. costs = [1]
  181. cost = 1
  182. for p in quads:
  183. assert q[-1] == p[0]
  184. for i in range(len(p) - 2):
  185. cost += 1
  186. costs.append(cost)
  187. costs.append(cost)
  188. qq = add_implicit_on_curves(p)[1:]
  189. costs.pop()
  190. q.extend(qq)
  191. cost += 1
  192. costs.append(cost)
  193. curves = spline_to_curves(q, costs, max_err, all_cubic)
  194. if not is_complex:
  195. curves = [tuple((c.real, c.imag) for c in curve) for curve in curves]
  196. return curves
  197. Solution = namedtuple("Solution", ["num_points", "error", "start_index", "is_cubic"])
  198. @cython.locals(
  199. i=cython.int,
  200. j=cython.int,
  201. k=cython.int,
  202. start=cython.int,
  203. i_sol_count=cython.int,
  204. j_sol_count=cython.int,
  205. this_sol_count=cython.int,
  206. tolerance=cython.double,
  207. err=cython.double,
  208. error=cython.double,
  209. i_sol_error=cython.double,
  210. j_sol_error=cython.double,
  211. all_cubic=cython.int,
  212. is_cubic=cython.int,
  213. count=cython.int,
  214. p0=cython.complex,
  215. p1=cython.complex,
  216. p2=cython.complex,
  217. p3=cython.complex,
  218. v=cython.complex,
  219. u=cython.complex,
  220. )
  221. def spline_to_curves(q, costs, tolerance=0.5, all_cubic=False):
  222. """
  223. q: quadratic spline with alternating on-curve / off-curve points.
  224. costs: cumulative list of encoding cost of q in terms of number of
  225. points that need to be encoded. Implied on-curve points do not
  226. contribute to the cost. If all points need to be encoded, then
  227. costs will be range(1, len(q)+1).
  228. """
  229. assert len(q) >= 3, "quadratic spline requires at least 3 points"
  230. # Elevate quadratic segments to cubic
  231. elevated_quadratics = [
  232. elevate_quadratic(*q[i : i + 3]) for i in range(0, len(q) - 2, 2)
  233. ]
  234. # Find sharp corners; they have to be oncurves for sure.
  235. forced = set()
  236. for i in range(1, len(elevated_quadratics)):
  237. p0 = elevated_quadratics[i - 1][2]
  238. p1 = elevated_quadratics[i][0]
  239. p2 = elevated_quadratics[i][1]
  240. if abs(p1 - p0) + abs(p2 - p1) > tolerance + abs(p2 - p0):
  241. forced.add(i)
  242. # Dynamic-Programming to find the solution with fewest number of
  243. # cubic curves, and within those the one with smallest error.
  244. sols = [Solution(0, 0, 0, False)]
  245. impossible = Solution(len(elevated_quadratics) * 3 + 1, 0, 1, False)
  246. start = 0
  247. for i in range(1, len(elevated_quadratics) + 1):
  248. best_sol = impossible
  249. for j in range(start, i):
  250. j_sol_count, j_sol_error = sols[j].num_points, sols[j].error
  251. if not all_cubic:
  252. # Solution with quadratics between j:i
  253. this_count = costs[2 * i - 1] - costs[2 * j] + 1
  254. i_sol_count = j_sol_count + this_count
  255. i_sol_error = j_sol_error
  256. i_sol = Solution(i_sol_count, i_sol_error, i - j, False)
  257. if i_sol < best_sol:
  258. best_sol = i_sol
  259. if this_count <= 3:
  260. # Can't get any better than this in the path below
  261. continue
  262. # Fit elevated_quadratics[j:i] into one cubic
  263. try:
  264. curve, ts = merge_curves(elevated_quadratics, j, i - j)
  265. except ZeroDivisionError:
  266. continue
  267. # Now reconstruct the segments from the fitted curve
  268. reconstructed_iter = splitCubicAtTC(*curve, *ts)
  269. reconstructed = []
  270. # Knot errors
  271. error = 0
  272. for k, reconst in enumerate(reconstructed_iter):
  273. orig = elevated_quadratics[j + k]
  274. err = abs(reconst[3] - orig[3])
  275. error = max(error, err)
  276. if error > tolerance:
  277. break
  278. reconstructed.append(reconst)
  279. if error > tolerance:
  280. # Not feasible
  281. continue
  282. # Interior errors
  283. for k, reconst in enumerate(reconstructed):
  284. orig = elevated_quadratics[j + k]
  285. p0, p1, p2, p3 = tuple(v - u for v, u in zip(reconst, orig))
  286. if not cubic_farthest_fit_inside(p0, p1, p2, p3, tolerance):
  287. error = tolerance + 1
  288. break
  289. if error > tolerance:
  290. # Not feasible
  291. continue
  292. # Save best solution
  293. i_sol_count = j_sol_count + 3
  294. i_sol_error = max(j_sol_error, error)
  295. i_sol = Solution(i_sol_count, i_sol_error, i - j, True)
  296. if i_sol < best_sol:
  297. best_sol = i_sol
  298. if i_sol_count == 3:
  299. # Can't get any better than this
  300. break
  301. sols.append(best_sol)
  302. if i in forced:
  303. start = i
  304. # Reconstruct solution
  305. splits = []
  306. cubic = []
  307. i = len(sols) - 1
  308. while i:
  309. count, is_cubic = sols[i].start_index, sols[i].is_cubic
  310. splits.append(i)
  311. cubic.append(is_cubic)
  312. i -= count
  313. curves = []
  314. j = 0
  315. for i, is_cubic in reversed(list(zip(splits, cubic))):
  316. if is_cubic:
  317. curves.append(merge_curves(elevated_quadratics, j, i - j)[0])
  318. else:
  319. for k in range(j, i):
  320. curves.append(q[k * 2 : k * 2 + 3])
  321. j = i
  322. return curves
  323. def main():
  324. from fontTools.cu2qu.benchmark import generate_curve
  325. from fontTools.cu2qu import curve_to_quadratic
  326. tolerance = 0.05
  327. reconstruct_tolerance = tolerance * 1
  328. curve = generate_curve()
  329. quadratics = curve_to_quadratic(curve, tolerance)
  330. print(
  331. "cu2qu tolerance %g. qu2cu tolerance %g." % (tolerance, reconstruct_tolerance)
  332. )
  333. print("One random cubic turned into %d quadratics." % len(quadratics))
  334. curves = quadratic_to_curves([quadratics], reconstruct_tolerance)
  335. print("Those quadratics turned back into %d cubics. " % len(curves))
  336. print("Original curve:", curve)
  337. print("Reconstructed curve(s):", curves)
  338. if __name__ == "__main__":
  339. main()