bezierTools.py revision 05b4b4a27160e90307372f85dd99be69a9d972ff
1"""fontTools.misc.bezierTools.py -- tools for working with bezier path segments.""" 2 3 4__all__ = ["calcQuadraticBounds", "calcCubicBounds", "splitLine", "splitQuadratic", 5 "splitCubic", "solveQuadratic", "solveCubic"] 6 7 8from fontTools.misc.arrayTools import calcBounds 9import Numeric 10 11 12def calcQuadraticBounds(pt1, pt2, pt3): 13 """Return the bounding rectangle for a qudratic bezier segment. 14 pt1 and pt3 are the "anchor" points, pt2 is the "handle".""" 15 # convert points to Numeric arrays 16 pt1, pt2, pt3 = Numeric.array((pt1, pt2, pt3)) 17 18 # calc quadratic parameters 19 c = pt1 20 b = (pt2 - c) * 2.0 21 a = pt3 - c - b 22 23 # calc first derivative 24 ax, ay = a * 2 25 bx, by = b 26 roots = [] 27 if ax != 0: 28 roots.append(-bx/ax) 29 if ay != 0: 30 roots.append(-by/ay) 31 points = [a*t*t + b*t + c for t in roots if 0 <= t < 1] + [pt1, pt3] 32 return calcBounds(points) 33 34 35def calcCubicBounds(pt1, pt2, pt3, pt4): 36 """Return the bounding rectangle for a cubic bezier segment. 37 pt1 and pt4 are the "anchor" points, pt2 and pt3 are the "handles".""" 38 # convert points to Numeric arrays 39 pt1, pt2, pt3, pt4 = Numeric.array((pt1, pt2, pt3, pt4)) 40 41 # calc cubic parameters 42 d = pt1 43 c = (pt2 - d) * 3.0 44 b = (pt3 - pt2) * 3.0 - c 45 a = pt4 - d - c - b 46 47 # calc first derivative 48 ax, ay = a * 3.0 49 bx, by = b * 2.0 50 cx, cy = c 51 xRoots = [t for t in solveQuadratic(ax, bx, cx) if 0 <= t < 1] 52 yRoots = [t for t in solveQuadratic(ay, by, cy) if 0 <= t < 1] 53 roots = xRoots + yRoots 54 55 points = [(a*t*t*t + b*t*t + c * t + d) for t in roots] + [pt1, pt4] 56 return calcBounds(points) 57 58 59def splitLine(pt1, pt2, where, isHorizontal): 60 """Split the line between pt1 and pt2 at position 'where', which 61 is an x coordinate if isHorizontal is False, a y coordinate if 62 isHorizontal is True. Return a list of two line segments if the 63 line was successfully split, or a list containing the original 64 line.""" 65 pt1, pt2 = Numeric.array((pt1, pt2)) 66 a = (pt2 - pt1) 67 b = pt1 68 ax = a[isHorizontal] 69 if ax == 0: 70 return [(pt1, pt2)] 71 t = float(where - b[isHorizontal]) / ax 72 midPt = a * t + b 73 return [(pt1, midPt), (midPt, pt2)] 74 75 76def splitQuadratic(pt1, pt2, pt3, where, isHorizontal): 77 """Split the quadratic curve between pt1, pt2 and pt3 at position 'where', 78 which is an x coordinate if isHorizontal is False, a y coordinate if 79 isHorizontal is True. Return a list of curve segments.""" 80 pt1, pt2, pt3 = Numeric.array((pt1, pt2, pt3)) 81 c = pt1 82 b = (pt2 - c) * 2.0 83 a = pt3 - c - b 84 solutions = solveQuadratic(a[isHorizontal], b[isHorizontal], 85 c[isHorizontal] - where) 86 solutions = [t for t in solutions if 0 <= t < 1] 87 solutions.sort() 88 if not solutions: 89 return [(pt1, pt2, pt3)] 90 91 segments = [] 92 solutions.insert(0, 0.0) 93 solutions.append(1.0) 94 for i in range(len(solutions) - 1): 95 t1 = solutions[i] 96 t2 = solutions[i+1] 97 delta = (t2 - t1) 98 # calc new a, b and c 99 a1 = a * delta**2 100 b1 = (2*a*t1 + b) * delta 101 c1 = a*t1**2 + b*t1 + c 102 # calc new points 103 pt1 = c1 104 pt2 = (b1 * 0.5) + c1 105 pt3 = a1 + b1 + c1 106 segments.append((pt1, pt2, pt3)) 107 return segments 108 109 110def splitCubic(pt1, pt2, pt3, pt4, where, isHorizontal): 111 """Split the cubic curve between pt1, pt2, pt3 and pt4 at position 'where', 112 which is an x coordinate if isHorizontal is False, a y coordinate if 113 isHorizontal is True. Return a list of curve segments.""" 114 pt1, pt2, pt3, pt4 = Numeric.array((pt1, pt2, pt3, pt4)) 115 d = pt1 116 c = (pt2 - d) * 3.0 117 b = (pt3 - pt2) * 3.0 - c 118 a = pt4 - d - c - b 119 120 solutions = solveCubic(a[isHorizontal], b[isHorizontal], c[isHorizontal], 121 d[isHorizontal] - where) 122 solutions = [t for t in solutions if 0 <= t < 1] 123 solutions.sort() 124 if not solutions: 125 return [(pt1, pt2, pt3, pt4)] 126 127 segments = [] 128 solutions.insert(0, 0.0) 129 solutions.append(1.0) 130 for i in range(len(solutions) - 1): 131 t1 = solutions[i] 132 t2 = solutions[i+1] 133 delta = (t2 - t1) 134 # calc new a, b, c and d 135 a1 = a * delta**3 136 b1 = (3*a*t1 + b) * delta**2 137 c1 = (2*b*t1 + c + 3*a*t1**2) * delta 138 d1 = a*t1**3 + b*t1**2 + c*t1 + d 139 # calc new points 140 pt1 = d1 141 pt2 = (c1 / 3.0) + d1 142 pt3 = (b1 + c1) / 3.0 + pt2 143 pt4 = a1 + d1 + c1 + b1 144 segments.append((pt1, pt2, pt3, pt4)) 145 return segments 146 147 148# 149# Equation solvers. 150# 151 152from math import sqrt, acos, cos, pi 153 154 155def solveQuadratic(a, b, c, 156 sqrt=sqrt): 157 """Solve a quadratic equation where a, b and c are real. 158 a*x*x + b*x + c = 0 159 This function returns a list of roots. 160 """ 161 if a == 0.0: 162 if b == 0.0: 163 # We have a non-equation; therefore, we have no valid solution 164 roots = [] 165 else: 166 # We have a linear equation with 1 root. 167 roots = [-c/b] 168 else: 169 # We have a true quadratic equation. Apply the quadratic formula to find two roots. 170 DD = b*b - 4.0*a*c 171 if DD >= 0.0: 172 roots = [(-b+sqrt(DD))/2.0/a, (-b-sqrt(DD))/2.0/a] 173 else: 174 # complex roots, ignore 175 roots = [] 176 return roots 177 178 179def solveCubic(a, b, c, d, 180 abs=abs, pow=pow, sqrt=sqrt, cos=cos, acos=acos, pi=pi): 181 """Solve a cubic equation where a, b, c and d are real. 182 a*x*x*x + b*x*x + c*x + d = 0 183 This function returns a list of roots. 184 """ 185 # 186 # adapted from: 187 # CUBIC.C - Solve a cubic polynomial 188 # public domain by Ross Cottrell 189 # found at: http://www.strangecreations.com/library/snippets/Cubic.C 190 # 191 if abs(a) < 1e-6: 192 # don't just test for zero; for very small values of 'a' solveCubic() 193 # returns unreliable results, so we fall back to quad. 194 return solveQuadratic(b, c, d) 195 a1 = b/a 196 a2 = c/a 197 a3 = d/a 198 199 Q = (a1*a1 - 3.0*a2)/9.0 200 R = (2.0*a1*a1*a1 - 9.0*a1*a2 + 27.0*a3)/54.0 201 R2_Q3 = R*R - Q*Q*Q 202 203 if R2_Q3 <= 0: 204 theta = acos(R/sqrt(Q*Q*Q)) 205 x0 = -2.0*sqrt(Q)*cos(theta/3.0) - a1/3.0 206 x1 = -2.0*sqrt(Q)*cos((theta+2.0*pi)/3.0) - a1/3.0 207 x2 = -2.0*sqrt(Q)*cos((theta+4.0*pi)/3.0) - a1/3.0 208 return [x0, x1, x2] 209 else: 210 x = pow(sqrt(R2_Q3)+abs(R), 1/3.0) 211 x = x + Q/x 212 if R >= 0.0: 213 x = -x 214 x = x - a1/3.0 215 return [x] 216