pyx/path.py

   1 #!/usr/bin/env python
   2 # -*- coding: ISO-8859-1 -*-
   3 #
   4 #
   5 # Copyright (C) 2002-2004 Jörg Lehmann <joergl@users.sourceforge.net>
   6 # Copyright (C) 2003-2004 Michael Schindler <m-schindler@users.sourceforge.net>
   7 # Copyright (C) 2002-2004 André Wobst <wobsta@users.sourceforge.net>
   8 #
   9 # This file is part of PyX (http://pyx.sourceforge.net/).
  10 #
  11 # PyX is free software; you can redistribute it and/or modify
  12 # it under the terms of the GNU General Public License as published by
  13 # the Free Software Foundation; either version 2 of the License, or
  14 # (at your option) any later version.
  15 #
  16 # PyX is distributed in the hope that it will be useful,
  17 # but WITHOUT ANY WARRANTY; without even the implied warranty of
  18 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  19 # GNU General Public License for more details.
  20 #
  21 # You should have received a copy of the GNU General Public License
  22 # along with PyX; if not, write to the Free Software
  23 # Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
  24
  25 #       - exceptions: nocurrentpoint, paramrange
  26 #       - correct bbox for curveto and normcurve
  27 #         (maybe we still need the current bbox implementation (then maybe called
  28 #          cbox = control box) for normcurve for the use during the
  29 #          intersection of bpaths)
  30
  31 import copy, math, bisect
  32 from math import cos, sin, pi
  33 try:
  34     from math import radians, degrees
  35 except ImportError:
  36     # fallback implementation for Python 2.1 and below
  37     def radians(x): return x*pi/180
  38     def degrees(x): return x*180/pi
  39 import base, bbox, trafo, unit, helper
  40
  41 try:
  42     sum([])
  43 except NameError:
  44     # fallback implementation for Python 2.2. and below
  45     def sum(list):
  46         return reduce(lambda x, y: x+y, list, 0)
  47
  48 try:
  49     enumerate([])
  50 except NameError:
  51     # fallback implementation for Python 2.2. and below
  52     def enumerate(list):
  53         return zip(xrange(len(list)), list)
  54
  55 # use new style classes when possible
  56 __metaclass__ = type
  57
  58 ################################################################################
  59 # Bezier helper functions
  60 ################################################################################
  61
  62 def _arctobcurve(x_pt, y_pt, r_pt, phi1, phi2):
  63     """generate the best bpathel corresponding to an arc segment"""
  64
  65     dphi = phi2-phi1
  66
  67     if dphi==0: return None
  68
  69     # the two endpoints should be clear
  70     x0_pt, y0_pt = x_pt+r_pt*cos(phi1), y_pt+r_pt*sin(phi1)
  71     x3_pt, y3_pt = x_pt+r_pt*cos(phi2), y_pt+r_pt*sin(phi2)
  72
  73     # optimal relative distance along tangent for second and third
  74     # control point
  75     l = r_pt*4*(1-cos(dphi/2))/(3*sin(dphi/2))
  76
  77     x1_pt, y1_pt = x0_pt-l*sin(phi1), y0_pt+l*cos(phi1)
  78     x2_pt, y2_pt = x3_pt+l*sin(phi2), y3_pt-l*cos(phi2)
  79
  80     return normcurve(x0_pt, y0_pt, x1_pt, y1_pt, x2_pt, y2_pt, x3_pt, y3_pt)
  81
  82
  83 def _arctobezierpath(x_pt, y_pt, r_pt, phi1, phi2, dphimax=45):
  84     apath = []
  85
  86     phi1 = radians(phi1)
  87     phi2 = radians(phi2)
  88     dphimax = radians(dphimax)
  89
  90     if phi2<phi1:
  91         # guarantee that phi2>phi1 ...
  92         phi2 = phi2 + (math.floor((phi1-phi2)/(2*pi))+1)*2*pi
  93     elif phi2>phi1+2*pi:
  94         # ... or remove unnecessary multiples of 2*pi
  95         phi2 = phi2 - (math.floor((phi2-phi1)/(2*pi))-1)*2*pi
  96
  97     if r_pt == 0 or phi1-phi2 == 0: return []
  98
  99     subdivisions = abs(int((1.0*(phi1-phi2))/dphimax))+1
 100
 101     dphi = (1.0*(phi2-phi1))/subdivisions
 102
 103     for i in range(subdivisions):
 104         apath.append(_arctobcurve(x_pt, y_pt, r_pt, phi1+i*dphi, phi1+(i+1)*dphi))
 105
 106     return apath
 107
 108 #
 109 # we define one exception
 110 #
 111
 112 class PathException(Exception): pass
 113
 114 ################################################################################
 115 # _pathcontext: context during walk along path
 116 ################################################################################
 117
 118 class _pathcontext:
 119
 120     """context during walk along path"""
 121
 122     __slots__ = "currentpoint", "currentsubpath"
 123
 124     def __init__(self, currentpoint=None, currentsubpath=None):
 125         """ initialize context
 126
 127         currentpoint:   position of current point
 128         currentsubpath: position of first point of current subpath
 129
 130         """
 131
 132         self.currentpoint = currentpoint
 133         self.currentsubpath = currentsubpath
 134
 135 ################################################################################
 136 # pathel: element of a PS style path
 137 ################################################################################
 138
 139 class pathel(base.PSOp):
 140
 141     """element of a PS style path"""
 142
 143     def _updatecontext(self, context):
 144         """update context of during walk along pathel
 145
 146         changes context in place
 147         """
 148         pass
 149
 150
 151     def _bbox(self, context):
 152         """calculate bounding box of pathel
 153
 154         context: context of pathel
 155
 156         returns bounding box of pathel (in given context)
 157
 158         Important note: all coordinates in bbox, currentpoint, and
 159         currrentsubpath have to be floats (in unit.topt)
 160
 161         """
 162         pass
 163
 164     def _normalized(self, context):
 165         """returns list of normalized version of pathel
 166
 167         context: context of pathel
 168
 169         Returns the path converted into a list of closepath, moveto_pt,
 170         normline, or normcurve instances.
 171
 172         """
 173         pass
 174
 175     def outputPS(self, file):
 176         """write PS code corresponding to pathel to file"""
 177         pass
 178
 179     def outputPDF(self, file):
 180         """write PDF code corresponding to pathel to file"""
 181         pass
 182
 183 #
 184 # various pathels
 185 #
 186 # Each one comes in two variants:
 187 #  - one which requires the coordinates to be already in pts (mainly
 188 #    used for internal purposes)
 189 #  - another which accepts arbitrary units
 190
 191 class closepath(pathel):
 192
 193     """Connect subpath back to its starting point"""
 194
 195     __slots__ = ()
 196
 197     def __str__(self):
 198         return "closepath"
 199
 200     def _updatecontext(self, context):
 201         context.currentpoint = None
 202         context.currentsubpath = None
 203
 204     def _bbox(self, context):
 205         x0_pt, y0_pt = context.currentpoint
 206         x1_pt, y1_pt = context.currentsubpath
 207
 208         return bbox.bbox_pt(min(x0_pt, x1_pt), min(y0_pt, y1_pt),
 209                           max(x0_pt, x1_pt), max(y0_pt, y1_pt))
 210
 211     def _normalized(self, context):
 212         return [closepath()]
 213
 214     def outputPS(self, file):
 215         file.write("closepath\n")
 216
 217     def outputPDF(self, file):
 218         file.write("h\n")
 219
 220
 221 class moveto_pt(pathel):
 222
 223     """Set current point to (x_pt, y_pt) (coordinates in pts)"""
 224
 225     __slots__ = "x_pt", "y_pt"
 226
 227     def __init__(self, x_pt, y_pt):
 228          self.x_pt = x_pt
 229          self.y_pt = y_pt
 230
 231     def __str__(self):
 232         return "%g %g moveto" % (self.x_pt, self.y_pt)
 233
 234     def _updatecontext(self, context):
 235         context.currentpoint = self.x_pt, self.y_pt
 236         context.currentsubpath = self.x_pt, self.y_pt
 237
 238     def _bbox(self, context):
 239         return None
 240
 241     def _normalized(self, context):
 242         return [moveto_pt(self.x_pt, self.y_pt)]
 243
 244     def outputPS(self, file):
 245         file.write("%g %g moveto\n" % (self.x_pt, self.y_pt) )
 246
 247     def outputPDF(self, file):
 248         file.write("%g %g m\n" % (self.x_pt, self.y_pt) )
 249
 250
 251 class lineto_pt(pathel):
 252
 253     """Append straight line to (x_pt, y_pt) (coordinates in pts)"""
 254
 255     __slots__ = "x_pt", "y_pt"
 256
 257     def __init__(self, x_pt, y_pt):
 258          self.x_pt = x_pt
 259          self.y_pt = y_pt
 260
 261     def __str__(self):
 262         return "%g %g lineto" % (self.x_pt, self.y_pt)
 263
 264     def _updatecontext(self, context):
 265         context.currentsubpath = context.currentsubpath or context.currentpoint
 266         context.currentpoint = self.x_pt, self.y_pt
 267
 268     def _bbox(self, context):
 269         return bbox.bbox_pt(min(context.currentpoint[0], self.x_pt),
 270                           min(context.currentpoint[1], self.y_pt),
 271                           max(context.currentpoint[0], self.x_pt),
 272                           max(context.currentpoint[1], self.y_pt))
 273
 274     def _normalized(self, context):
 275         return [normline(context.currentpoint[0], context.currentpoint[1], self.x_pt, self.y_pt)]
 276
 277     def outputPS(self, file):
 278         file.write("%g %g lineto\n" % (self.x_pt, self.y_pt) )
 279
 280     def outputPDF(self, file):
 281         file.write("%g %g l\n" % (self.x_pt, self.y_pt) )
 282
 283
 284 class curveto_pt(pathel):
 285
 286     """Append curveto (coordinates in pts)"""
 287
 288     __slots__ = "x1_pt", "y1_pt", "x2_pt", "y2_pt", "x3_pt", "y3_pt"
 289
 290     def __init__(self, x1_pt, y1_pt, x2_pt, y2_pt, x3_pt, y3_pt):
 291         self.x1_pt = x1_pt
 292         self.y1_pt = y1_pt
 293         self.x2_pt = x2_pt
 294         self.y2_pt = y2_pt
 295         self.x3_pt = x3_pt
 296         self.y3_pt = y3_pt
 297
 298     def __str__(self):
 299         return "%g %g %g %g %g %g curveto" % (self.x1_pt, self.y1_pt,
 300                                               self.x2_pt, self.y2_pt,
 301                                               self.x3_pt, self.y3_pt)
 302
 303     def _updatecontext(self, context):
 304         context.currentsubpath = context.currentsubpath or context.currentpoint
 305         context.currentpoint = self.x3_pt, self.y3_pt
 306
 307     def _bbox(self, context):
 308         return bbox.bbox_pt(min(context.currentpoint[0], self.x1_pt, self.x2_pt, self.x3_pt),
 309                           min(context.currentpoint[1], self.y1_pt, self.y2_pt, self.y3_pt),
 310                           max(context.currentpoint[0], self.x1_pt, self.x2_pt, self.x3_pt),
 311                           max(context.currentpoint[1], self.y1_pt, self.y2_pt, self.y3_pt))
 312
 313     def _normalized(self, context):
 314         return [normcurve(context.currentpoint[0], context.currentpoint[1],
 315                           self.x1_pt, self.y1_pt,
 316                           self.x2_pt, self.y2_pt,
 317                           self.x3_pt, self.y3_pt)]
 318
 319     def outputPS(self, file):
 320         file.write("%g %g %g %g %g %g curveto\n" % ( self.x1_pt, self.y1_pt,
 321                                                      self.x2_pt, self.y2_pt,
 322                                                      self.x3_pt, self.y3_pt ) )
 323
 324     def outputPDF(self, file):
 325         file.write("%f %f %f %f %f %f c\n" % ( self.x1_pt, self.y1_pt,
 326                                                self.x2_pt, self.y2_pt,
 327                                                self.x3_pt, self.y3_pt ) )
 328
 329
 330 class rmoveto_pt(pathel):
 331
 332     """Perform relative moveto (coordinates in pts)"""
 333
 334     __slots__ = "dx_pt", "dy_pt"
 335
 336     def __init__(self, dx_pt, dy_pt):
 337          self.dx_pt = dx_pt
 338          self.dy_pt = dy_pt
 339
 340     def _updatecontext(self, context):
 341         context.currentpoint = (context.currentpoint[0] + self.dx_pt,
 342                                 context.currentpoint[1] + self.dy_pt)
 343         context.currentsubpath = context.currentpoint
 344
 345     def _bbox(self, context):
 346         return None
 347
 348     def _normalized(self, context):
 349         x_pt = context.currentpoint[0]+self.dx_pt
 350         y_pt = context.currentpoint[1]+self.dy_pt
 351         return [moveto_pt(x_pt, y_pt)]
 352
 353     def outputPS(self, file):
 354         file.write("%g %g rmoveto\n" % (self.dx_pt, self.dy_pt) )
 355
 356
 357 class rlineto_pt(pathel):
 358
 359     """Perform relative lineto (coordinates in pts)"""
 360
 361     __slots__ = "dx_pt", "dy_pt"
 362
 363     def __init__(self, dx_pt, dy_pt):
 364          self.dx_pt = dx_pt
 365          self.dy_pt = dy_pt
 366
 367     def _updatecontext(self, context):
 368         context.currentsubpath = context.currentsubpath or context.currentpoint
 369         context.currentpoint = (context.currentpoint[0]+self.dx_pt,
 370                                 context.currentpoint[1]+self.dy_pt)
 371
 372     def _bbox(self, context):
 373         x = context.currentpoint[0] + self.dx_pt
 374         y = context.currentpoint[1] + self.dy_pt
 375         return bbox.bbox_pt(min(context.currentpoint[0], x),
 376                           min(context.currentpoint[1], y),
 377                           max(context.currentpoint[0], x),
 378                           max(context.currentpoint[1], y))
 379
 380     def _normalized(self, context):
 381         x0_pt = context.currentpoint[0]
 382         y0_pt = context.currentpoint[1]
 383         return [normline(x0_pt, y0_pt, x0_pt+self.dx_pt, y0_pt+self.dy_pt)]
 384
 385     def outputPS(self, file):
 386         file.write("%g %g rlineto\n" % (self.dx_pt, self.dy_pt) )
 387
 388
 389 class rcurveto_pt(pathel):
 390
 391     """Append rcurveto (coordinates in pts)"""
 392
 393     __slots__ = "dx1_pt", "dy1_pt", "dx2_pt", "dy2_pt", "dx3_pt", "dy3_pt"
 394
 395     def __init__(self, dx1_pt, dy1_pt, dx2_pt, dy2_pt, dx3_pt, dy3_pt):
 396         self.dx1_pt = dx1_pt
 397         self.dy1_pt = dy1_pt
 398         self.dx2_pt = dx2_pt
 399         self.dy2_pt = dy2_pt
 400         self.dx3_pt = dx3_pt
 401         self.dy3_pt = dy3_pt
 402
 403     def outputPS(self, file):
 404         file.write("%g %g %g %g %g %g rcurveto\n" % ( self.dx1_pt, self.dy1_pt,
 405                                                     self.dx2_pt, self.dy2_pt,
 406                                                     self.dx3_pt, self.dy3_pt ) )
 407
 408     def _updatecontext(self, context):
 409         x3_pt = context.currentpoint[0]+self.dx3_pt
 410         y3_pt = context.currentpoint[1]+self.dy3_pt
 411
 412         context.currentsubpath = context.currentsubpath or context.currentpoint
 413         context.currentpoint = x3_pt, y3_pt
 414
 415
 416     def _bbox(self, context):
 417         x1_pt = context.currentpoint[0]+self.dx1_pt
 418         y1_pt = context.currentpoint[1]+self.dy1_pt
 419         x2_pt = context.currentpoint[0]+self.dx2_pt
 420         y2_pt = context.currentpoint[1]+self.dy2_pt
 421         x3_pt = context.currentpoint[0]+self.dx3_pt
 422         y3_pt = context.currentpoint[1]+self.dy3_pt
 423         return bbox.bbox_pt(min(context.currentpoint[0], x1_pt, x2_pt, x3_pt),
 424                           min(context.currentpoint[1], y1_pt, y2_pt, y3_pt),
 425                           max(context.currentpoint[0], x1_pt, x2_pt, x3_pt),
 426                           max(context.currentpoint[1], y1_pt, y2_pt, y3_pt))
 427
 428     def _normalized(self, context):
 429         x0_pt = context.currentpoint[0]
 430         y0_pt = context.currentpoint[1]
 431         return [normcurve(x0_pt, y0_pt, x0_pt+self.dx1_pt, y0_pt+self.dy1_pt, x0_pt+self.dx2_pt, y0_pt+self.dy2_pt, x0_pt+self.dx3_pt, y0_pt+self.dy3_pt)]
 432
 433
 434 class arc_pt(pathel):
 435
 436     """Append counterclockwise arc (coordinates in pts)"""
 437
 438     __slots__ = "x_pt", "y_pt", "r_pt", "angle1", "angle2"
 439
 440     def __init__(self, x_pt, y_pt, r_pt, angle1, angle2):
 441         self.x_pt = x_pt
 442         self.y_pt = y_pt
 443         self.r_pt = r_pt
 444         self.angle1 = angle1
 445         self.angle2 = angle2
 446
 447     def _sarc(self):
 448         """Return starting point of arc segment"""
 449         return (self.x_pt+self.r_pt*cos(radians(self.angle1)),
 450                 self.y_pt+self.r_pt*sin(radians(self.angle1)))
 451
 452     def _earc(self):
 453         """Return end point of arc segment"""
 454         return (self.x_pt+self.r_pt*cos(radians(self.angle2)),
 455                 self.y_pt+self.r_pt*sin(radians(self.angle2)))
 456
 457     def _updatecontext(self, context):
 458         if context.currentpoint:
 459             context.currentsubpath = context.currentsubpath or context.currentpoint
 460         else:
 461             # we assert that currentsubpath is also None
 462             context.currentsubpath = self._sarc()
 463
 464         context.currentpoint = self._earc()
 465
 466     def _bbox(self, context):
 467         phi1 = radians(self.angle1)
 468         phi2 = radians(self.angle2)
 469
 470         # starting end end point of arc segment
 471         sarcx_pt, sarcy_pt = self._sarc()
 472         earcx_pt, earcy_pt = self._earc()
 473
 474         # Now, we have to determine the corners of the bbox for the
 475         # arc segment, i.e. global maxima/mimima of cos(phi) and sin(phi)
 476         # in the interval [phi1, phi2]. These can either be located
 477         # on the borders of this interval or in the interior.
 478
 479         if phi2 < phi1:
 480             # guarantee that phi2>phi1
 481             phi2 = phi2 + (math.floor((phi1-phi2)/(2*pi))+1)*2*pi
 482
 483         # next minimum of cos(phi) looking from phi1 in counterclockwise
 484         # direction: 2*pi*floor((phi1-pi)/(2*pi)) + 3*pi
 485
 486         if phi2 < (2*math.floor((phi1-pi)/(2*pi))+3)*pi:
 487             minarcx_pt = min(sarcx_pt, earcx_pt)
 488         else:
 489             minarcx_pt = self.x_pt-self.r_pt
 490
 491         # next minimum of sin(phi) looking from phi1 in counterclockwise
 492         # direction: 2*pi*floor((phi1-3*pi/2)/(2*pi)) + 7/2*pi
 493
 494         if phi2 < (2*math.floor((phi1-3.0*pi/2)/(2*pi))+7.0/2)*pi:
 495             minarcy_pt = min(sarcy_pt, earcy_pt)
 496         else:
 497             minarcy_pt = self.y_pt-self.r_pt
 498
 499         # next maximum of cos(phi) looking from phi1 in counterclockwise
 500         # direction: 2*pi*floor((phi1)/(2*pi))+2*pi
 501
 502         if phi2 < (2*math.floor((phi1)/(2*pi))+2)*pi:
 503             maxarcx_pt = max(sarcx_pt, earcx_pt)
 504         else:
 505             maxarcx_pt = self.x_pt+self.r_pt
 506
 507         # next maximum of sin(phi) looking from phi1 in counterclockwise
 508         # direction: 2*pi*floor((phi1-pi/2)/(2*pi)) + 1/2*pi
 509
 510         if phi2 < (2*math.floor((phi1-pi/2)/(2*pi))+5.0/2)*pi:
 511             maxarcy_pt = max(sarcy_pt, earcy_pt)
 512         else:
 513             maxarcy_pt = self.y_pt+self.r_pt
 514
 515         # Finally, we are able to construct the bbox for the arc segment.
 516         # Note that if there is a currentpoint defined, we also
 517         # have to include the straight line from this point
 518         # to the first point of the arc segment
 519
 520         if context.currentpoint:
 521             return (bbox.bbox_pt(min(context.currentpoint[0], sarcx_pt),
 522                                min(context.currentpoint[1], sarcy_pt),
 523                                max(context.currentpoint[0], sarcx_pt),
 524                                max(context.currentpoint[1], sarcy_pt)) +
 525                     bbox.bbox_pt(minarcx_pt, minarcy_pt, maxarcx_pt, maxarcy_pt)
 526                     )
 527         else:
 528             return  bbox.bbox_pt(minarcx_pt, minarcy_pt, maxarcx_pt, maxarcy_pt)
 529
 530     def _normalized(self, context):
 531         # get starting and end point of arc segment and bpath corresponding to arc
 532         sarcx_pt, sarcy_pt = self._sarc()
 533         earcx_pt, earcy_pt = self._earc()
 534         barc = _arctobezierpath(self.x_pt, self.y_pt, self.r_pt, self.angle1, self.angle2)
 535
 536         # convert to list of curvetos omitting movetos
 537         nbarc = []
 538
 539         for bpathel in barc:
 540             nbarc.append(normcurve(bpathel.x0_pt, bpathel.y0_pt,
 541                                    bpathel.x1_pt, bpathel.y1_pt,
 542                                    bpathel.x2_pt, bpathel.y2_pt,
 543                                    bpathel.x3_pt, bpathel.y3_pt))
 544
 545         # Note that if there is a currentpoint defined, we also
 546         # have to include the straight line from this point
 547         # to the first point of the arc segment.
 548         # Otherwise, we have to add a moveto at the beginning
 549         if context.currentpoint:
 550             return [normline(context.currentpoint[0], context.currentpoint[1], sarcx_pt, sarcy_pt)] + nbarc
 551         else:
 552             return [moveto_pt(sarcx_pt, sarcy_pt)] + nbarc
 553
 554     def outputPS(self, file):
 555         file.write("%g %g %g %g %g arc\n" % ( self.x_pt, self.y_pt,
 556                                             self.r_pt,
 557                                             self.angle1,
 558                                             self.angle2 ) )
 559
 560
 561 class arcn_pt(pathel):
 562
 563     """Append clockwise arc (coordinates in pts)"""
 564
 565     __slots__ = "x_pt", "y_pt", "r_pt", "angle1", "angle2"
 566
 567     def __init__(self, x_pt, y_pt, r_pt, angle1, angle2):
 568         self.x_pt = x_pt
 569         self.y_pt = y_pt
 570         self.r_pt = r_pt
 571         self.angle1 = angle1
 572         self.angle2 = angle2
 573
 574     def _sarc(self):
 575         """Return starting point of arc segment"""
 576         return (self.x_pt+self.r_pt*cos(radians(self.angle1)),
 577                 self.y_pt+self.r_pt*sin(radians(self.angle1)))
 578
 579     def _earc(self):
 580         """Return end point of arc segment"""
 581         return (self.x_pt+self.r_pt*cos(radians(self.angle2)),
 582                 self.y_pt+self.r_pt*sin(radians(self.angle2)))
 583
 584     def _updatecontext(self, context):
 585         if context.currentpoint:
 586             context.currentsubpath = context.currentsubpath or context.currentpoint
 587         else:  # we assert that currentsubpath is also None
 588             context.currentsubpath = self._sarc()
 589
 590         context.currentpoint = self._earc()
 591
 592     def _bbox(self, context):
 593         # in principle, we obtain bbox of an arcn element from
 594         # the bounding box of the corrsponding arc element with
 595         # angle1 and angle2 interchanged. Though, we have to be carefull
 596         # with the straight line segment, which is added if currentpoint
 597         # is defined.
 598
 599         # Hence, we first compute the bbox of the arc without this line:
 600
 601         a = arc_pt(self.x_pt, self.y_pt, self.r_pt,
 602                  self.angle2,
 603                  self.angle1)
 604
 605         sarcx_pt, sarcy_pt = self._sarc()
 606         arcbb = a._bbox(_pathcontext())
 607
 608         # Then, we repeat the logic from arc.bbox, but with interchanged
 609         # start and end points of the arc
 610
 611         if context.currentpoint:
 612             return  bbox.bbox_pt(min(context.currentpoint[0], sarcx_pt),
 613                                min(context.currentpoint[1], sarcy_pt),
 614                                max(context.currentpoint[0], sarcx_pt),
 615                                max(context.currentpoint[1], sarcy_pt))+ arcbb
 616         else:
 617             return arcbb
 618
 619     def _normalized(self, context):
 620         # get starting and end point of arc segment and bpath corresponding to arc
 621         sarcx_pt, sarcy_pt = self._sarc()
 622         earcx_pt, earcy_pt = self._earc()
 623         barc = _arctobezierpath(self.x_pt, self.y_pt, self.r_pt, self.angle2, self.angle1)
 624         barc.reverse()
 625
 626         # convert to list of curvetos omitting movetos
 627         nbarc = []
 628
 629         for bpathel in barc:
 630             nbarc.append(normcurve(bpathel.x3_pt, bpathel.y3_pt,
 631                                    bpathel.x2_pt, bpathel.y2_pt,
 632                                    bpathel.x1_pt, bpathel.y1_pt,
 633                                    bpathel.x0_pt, bpathel.y0_pt))
 634
 635         # Note that if there is a currentpoint defined, we also
 636         # have to include the straight line from this point
 637         # to the first point of the arc segment.
 638         # Otherwise, we have to add a moveto at the beginning
 639         if context.currentpoint:
 640             return [normline(context.currentpoint[0], context.currentpoint[1], sarcx_pt, sarcy_pt)] + nbarc
 641         else:
 642             return [moveto_pt(sarcx_pt, sarcy_pt)] + nbarc
 643
 644
 645     def outputPS(self, file):
 646         file.write("%g %g %g %g %g arcn\n" % ( self.x_pt, self.y_pt,
 647                                                self.r_pt,
 648                                                self.angle1,
 649                                                self.angle2 ) )
 650
 651
 652 class arct_pt(pathel):
 653
 654     """Append tangent arc (coordinates in pts)"""
 655
 656     __slots__ = "x1_pt", "y1_pt", "x2_pt", "y2_pt", "r_pt"
 657
 658     def __init__(self, x1_pt, y1_pt, x2_pt, y2_pt, r_pt):
 659         self.x1_pt = x1_pt
 660         self.y1_pt = y1_pt
 661         self.x2_pt = x2_pt
 662         self.y2_pt = y2_pt
 663         self.r_pt = r_pt
 664
 665     def _path(self, currentpoint, currentsubpath):
 666         """returns new currentpoint, currentsubpath and path consisting
 667         of arc and/or line which corresponds to arct
 668
 669         this is a helper routine for _bbox and _normalized, which both need
 670         this path. Note: we don't want to calculate the bbox from a bpath
 671
 672         """
 673
 674         # direction and length of tangent 1
 675         dx1_pt = currentpoint[0]-self.x1_pt
 676         dy1_pt = currentpoint[1]-self.y1_pt
 677         l1 = math.hypot(dx1_pt, dy1_pt)
 678
 679         # direction and length of tangent 2
 680         dx2_pt = self.x2_pt-self.x1_pt
 681         dy2_pt = self.y2_pt-self.y1_pt
 682         l2 = math.hypot(dx2_pt, dy2_pt)
 683
 684         # intersection angle between two tangents
 685         alpha = math.acos((dx1_pt*dx2_pt+dy1_pt*dy2_pt)/(l1*l2))
 686
 687         if math.fabs(sin(alpha)) >= 1e-15 and 1.0+self.r_pt != 1.0:
 688             cotalpha2 = 1.0/math.tan(alpha/2)
 689
 690             # two tangent points
 691             xt1_pt = self.x1_pt + dx1_pt*self.r_pt*cotalpha2/l1
 692             yt1_pt = self.y1_pt + dy1_pt*self.r_pt*cotalpha2/l1
 693             xt2_pt = self.x1_pt + dx2_pt*self.r_pt*cotalpha2/l2
 694             yt2_pt = self.y1_pt + dy2_pt*self.r_pt*cotalpha2/l2
 695
 696             # direction of center of arc
 697             rx_pt = self.x1_pt - 0.5*(xt1_pt+xt2_pt)
 698             ry_pt = self.y1_pt - 0.5*(yt1_pt+yt2_pt)
 699             lr = math.hypot(rx_pt, ry_pt)
 700
 701             # angle around which arc is centered
 702             if rx_pt >= 0:
 703                 phi = degrees(math.atan2(ry_pt, rx_pt))
 704             else:
 705                 # XXX why is rx_pt/ry_pt and not ry_pt/rx_pt used???
 706                 phi = degrees(math.atan(rx_pt/ry_pt))+180
 707
 708             # half angular width of arc
 709             deltaphi = 90*(1-alpha/pi)
 710
 711             # center position of arc
 712             mx_pt = self.x1_pt - rx_pt*self.r_pt/(lr*sin(alpha/2))
 713             my_pt = self.y1_pt - ry_pt*self.r_pt/(lr*sin(alpha/2))
 714
 715             # now we are in the position to construct the path
 716             p = path(moveto_pt(*currentpoint))
 717
 718             if phi<0:
 719                 p.append(arc_pt(mx_pt, my_pt, self.r_pt, phi-deltaphi, phi+deltaphi))
 720             else:
 721                 p.append(arcn_pt(mx_pt, my_pt, self.r_pt, phi+deltaphi, phi-deltaphi))
 722
 723             return ( (xt2_pt, yt2_pt),
 724                      currentsubpath or (xt2_pt, yt2_pt),
 725                      p )
 726
 727         else:
 728             # we need no arc, so just return a straight line to currentpoint to x1_pt, y1_pt
 729             return  ( (self.x1_pt, self.y1_pt),
 730                       currentsubpath or (self.x1_pt, self.y1_pt),
 731                       line_pt(currentpoint[0], currentpoint[1], self.x1_pt, self.y1_pt) )
 732
 733     def _updatecontext(self, context):
 734         result = self._path(context.currentpoint, context.currentsubpath)
 735         context.currentpoint, context.currentsubpath = result[:2]
 736
 737     def _bbox(self, context):
 738         return self._path(context.currentpoint, context.currentsubpath)[2].bbox()
 739
 740     def _normalized(self, context):
 741         # XXX TODO
 742         return normpath(self._path(context.currentpoint,
 743                                    context.currentsubpath)[2]).subpaths[0].normpathels
 744     def outputPS(self, file):
 745         file.write("%g %g %g %g %g arct\n" % ( self.x1_pt, self.y1_pt,
 746                                                self.x2_pt, self.y2_pt,
 747                                                self.r_pt ) )
 748
 749 #
 750 # now the pathels that convert from user coordinates to pts
 751 #
 752
 753 class moveto(moveto_pt):
 754
 755     """Set current point to (x, y)"""
 756
 757     __slots__ = "x_pt", "y_pt"
 758
 759     def __init__(self, x, y):
 760         moveto_pt.__init__(self, unit.topt(x), unit.topt(y))
 761
 762
 763 class lineto(lineto_pt):
 764
 765     """Append straight line to (x, y)"""
 766
 767     __slots__ = "x_pt", "y_pt"
 768
 769     def __init__(self, x, y):
 770         lineto_pt.__init__(self, unit.topt(x), unit.topt(y))
 771
 772
 773 class curveto(curveto_pt):
 774
 775     """Append curveto"""
 776
 777     __slots__ = "x1_pt", "y1_pt", "x2_pt", "y2_pt", "x3_pt", "y3_pt"
 778
 779     def __init__(self, x1, y1, x2, y2, x3, y3):
 780         curveto_pt.__init__(self,
 781                             unit.topt(x1), unit.topt(y1),
 782                             unit.topt(x2), unit.topt(y2),
 783                             unit.topt(x3), unit.topt(y3))
 784
 785 class rmoveto(rmoveto_pt):
 786
 787     """Perform relative moveto"""
 788
 789     __slots__ = "dx_pt", "dy_pt"
 790
 791     def __init__(self, dx, dy):
 792         rmoveto_pt.__init__(self, unit.topt(dx), unit.topt(dy))
 793
 794
 795 class rlineto(rlineto_pt):
 796
 797     """Perform relative lineto"""
 798
 799     __slots__ = "dx_pt", "dy_pt"
 800
 801     def __init__(self, dx, dy):
 802         rlineto_pt.__init__(self, unit.topt(dx), unit.topt(dy))
 803
 804
 805 class rcurveto(rcurveto_pt):
 806
 807     """Append rcurveto"""
 808
 809     __slots__ = "dx1_pt", "dy1_pt", "dx2_pt", "dy2_pt", "dx3_pt", "dy3_pt"
 810
 811     def __init__(self, dx1, dy1, dx2, dy2, dx3, dy3):
 812         rcurveto_pt.__init__(self,
 813                              unit.topt(dx1), unit.topt(dy1),
 814                              unit.topt(dx2), unit.topt(dy2),
 815                              unit.topt(dx3), unit.topt(dy3))
 816
 817
 818 class arcn(arcn_pt):
 819
 820     """Append clockwise arc"""
 821
 822     __slots__ = "x_pt", "y_pt", "r_pt", "angle1", "angle2"
 823
 824     def __init__(self, x, y, r, angle1, angle2):
 825         arcn_pt.__init__(self, unit.topt(x), unit.topt(y), unit.topt(r), angle1, angle2)
 826
 827
 828 class arc(arc_pt):
 829
 830     """Append counterclockwise arc"""
 831
 832     __slots__ = "x_pt", "y_pt", "r_pt", "angle1", "angle2"
 833
 834     def __init__(self, x, y, r, angle1, angle2):
 835         arc_pt.__init__(self, unit.topt(x), unit.topt(y), unit.topt(r), angle1, angle2)
 836
 837
 838 class arct(arct_pt):
 839
 840     """Append tangent arc"""
 841
 842     __slots__ = "x1_pt", "y1_pt", "x2_pt", "y2_pt", "r"
 843
 844     def __init__(self, x1, y1, x2, y2, r):
 845         arct_pt.__init__(self, unit.topt(x1), unit.topt(y1),
 846                          unit.topt(x2), unit.topt(y2), unit.topt(r))
 847
 848 #
 849 # "combined" pathels provided for performance reasons
 850 #
 851
 852 class multilineto_pt(pathel):
 853
 854     """Perform multiple linetos (coordinates in pts)"""
 855
 856     __slots__ = "points_pt"
 857
 858     def __init__(self, points_pt):
 859         self.points_pt = points_pt
 860
 861     def _updatecontext(self, context):
 862         context.currentsubpath = context.currentsubpath or context.currentpoint
 863         context.currentpoint = self.points_pt[-1]
 864
 865     def _bbox(self, context):
 866         xs_pt = [point[0] for point in self.points_pt]
 867         ys_pt = [point[1] for point in self.points_pt]
 868         return bbox.bbox_pt(min(context.currentpoint[0], *xs_pt),
 869                           min(context.currentpoint[1], *ys_pt),
 870                           max(context.currentpoint[0], *xs_pt),
 871                           max(context.currentpoint[1], *ys_pt))
 872
 873     def _normalized(self, context):
 874         result = []
 875         x0_pt, y0_pt = context.currentpoint
 876         for x_pt, y_pt in self.points_pt:
 877             result.append(normline(x0_pt, y0_pt, x_pt, y_pt))
 878             x0_pt, y0_pt = x_pt, y_pt
 879         return result
 880
 881     def outputPS(self, file):
 882         for point_pt in self.points_pt:
 883             file.write("%g %g lineto\n" % point_pt )
 884
 885     def outputPDF(self, file):
 886         for point_pt in self.points_pt:
 887             file.write("%f %f l\n" % point_pt )
 888
 889
 890 class multicurveto_pt(pathel):
 891
 892     """Perform multiple curvetos (coordinates in pts)"""
 893
 894     __slots__ = "points_pt"
 895
 896     def __init__(self, points_pt):
 897         self.points_pt = points_pt
 898
 899     def _updatecontext(self, context):
 900         context.currentsubpath = context.currentsubpath or context.currentpoint
 901         context.currentpoint = self.points_pt[-1]
 902
 903     def _bbox(self, context):
 904         xs = ( [point[0] for point in self.points_pt] +
 905                [point[2] for point in self.points_pt] +
 906                [point[4] for point in self.points_pt] )
 907         ys = ( [point[1] for point in self.points_pt] +
 908                [point[3] for point in self.points_pt] +
 909                [point[5] for point in self.points_pt] )
 910         return bbox.bbox_pt(min(context.currentpoint[0], *xs_pt),
 911                           min(context.currentpoint[1], *ys_pt),
 912                           max(context.currentpoint[0], *xs_pt),
 913                           max(context.currentpoint[1], *ys_pt))
 914
 915     def _normalized(self, context):
 916         result = []
 917         x0_pt, y0_pt = context.currentpoint
 918         for point_pt in self.points_pt:
 919             result.append(normcurve(x0_pt, y0_pt, *point_pt))
 920             x0_pt, y0_pt = point_pt[4:]
 921         return result
 922
 923     def outputPS(self, file):
 924         for point_pt in self.points_pt:
 925             file.write("%g %g %g %g %g %g curveto\n" % point_pt)
 926
 927     def outputPDF(self, file):
 928         for point_pt in self.points_pt:
 929             file.write("%f %f %f %f %f %f c\n" % point_pt)
 930
 931
 932 ################################################################################
 933 # path: PS style path
 934 ################################################################################
 935
 936 class path(base.PSCmd):
 937
 938     """PS style path"""
 939
 940     __slots__ = "path"
 941
 942     def __init__(self, *args):
 943         if len(args)==1 and isinstance(args[0], path):
 944             self.path = args[0].path
 945         else:
 946             self.path = list(args)
 947
 948     def __add__(self, other):
 949         return path(*(self.path+other.path))
 950
 951     def __iadd__(self, other):
 952         self.path += other.path
 953         return self
 954
 955     def __getitem__(self, i):
 956         return self.path[i]
 957
 958     def __len__(self):
 959         return len(self.path)
 960
 961     def append(self, pathel):
 962         self.path.append(pathel)
 963
 964     def arclen_pt(self):
 965         """returns total arc length of path in pts"""
 966         return normpath(self).arclen_pt()
 967
 968     def arclen(self):
 969         """returns total arc length of path"""
 970         return normpath(self).arclen()
 971
 972     def arclentoparam(self, lengths):
 973         """returns the parameter value(s) matching the given length(s)"""
 974         return normpath(self).arclentoparam(lengths)
 975
 976     def at_pt(self, param=None, arclen=None):
 977         """return coordinates of path in pts at either parameter value param
 978         or arc length arclen.
 979
 980         At discontinuities in the path, the limit from below is returned
 981         """
 982         return normpath(self).at_pt(param, arclen)
 983
 984     def at(self, param=None, arclen=None):
 985         """return coordinates of path at either parameter value param
 986         or arc length arclen.
 987
 988         At discontinuities in the path, the limit from below is returned
 989         """
 990         return normpath(self).at(param, arclen)
 991
 992     def bbox(self):
 993         context = _pathcontext()
 994         abbox = None
 995
 996         for pel in self.path:
 997             nbbox =  pel._bbox(context)
 998             pel._updatecontext(context)
 999             if abbox is None:
1000                 abbox = nbbox
1001             elif nbbox:
1002                 abbox += nbbox
1003
1004         return abbox
1005
1006     def begin_pt(self):
1007         """return coordinates of first point of first subpath in path (in pts)"""
1008         return normpath(self).begin_pt()
1009
1010     def begin(self):
1011         """return coordinates of first point of first subpath in path"""
1012         return normpath(self).begin()
1013
1014     def curvradius_pt(self, param=None, arclen=None):
1015         """Returns the curvature radius in pts (or None if infinite)
1016         at parameter param or arc length arclen.  This is the inverse
1017         of the curvature at this parameter
1018
1019         Please note that this radius can be negative or positive,
1020         depending on the sign of the curvature"""
1021         return normpath(self).curvradius_pt(param, arclen)
1022
1023     def curvradius(self, param=None, arclen=None):
1024         """Returns the curvature radius (or None if infinite) at
1025         parameter param or arc length arclen.  This is the inverse of
1026         the curvature at this parameter
1027
1028         Please note that this radius can be negative or positive,
1029         depending on the sign of the curvature"""
1030         return normpath(self).curvradius(param, arclen)
1031
1032     def end_pt(self):
1033         """return coordinates of last point of last subpath in path (in pts)"""
1034         return normpath(self).end_pt()
1035
1036     def end(self):
1037         """return coordinates of last point of last subpath in path"""
1038         return normpath(self).end()
1039
1040     def joined(self, other):
1041         """return path consisting of self and other joined together"""
1042         return normpath(self).joined(other)
1043
1044     # << operator also designates joining
1045     __lshift__ = joined
1046
1047     def intersect(self, other):
1048         """intersect normpath corresponding to self with other path"""
1049         return normpath(self).intersect(other)
1050
1051     def range(self):
1052         """return maximal value for parameter value t for corr. normpath"""
1053         return normpath(self).range()
1054
1055     def reversed(self):
1056         """return reversed path"""
1057         return normpath(self).reversed()
1058
1059     def split(self, params):
1060         """return corresponding normpaths split at parameter values params"""
1061         return normpath(self).split(params)
1062
1063     def tangent(self, param=None, arclen=None, length=None):
1064         """return tangent vector of path at either parameter value param
1065         or arc length arclen.
1066
1067         At discontinuities in the path, the limit from below is returned.
1068         If length is not None, the tangent vector will be scaled to
1069         the desired length.
1070         """
1071         return normpath(self).tangent(param, arclen, length)
1072
1073     def trafo(self, param=None, arclen=None):
1074         """return transformation at either parameter value param or arc length arclen"""
1075         return normpath(self).trafo(param, arclen)
1076
1077     def transformed(self, trafo):
1078         """return transformed path"""
1079         return normpath(self).transformed(trafo)
1080
1081     def outputPS(self, file):
1082         if not (isinstance(self.path[0], moveto_pt) or
1083                 isinstance(self.path[0], arc_pt) or
1084                 isinstance(self.path[0], arcn_pt)):
1085             raise PathException("first path element must be either moveto, arc, or arcn")
1086         for pel in self.path:
1087             pel.outputPS(file)
1088
1089     def outputPDF(self, file):
1090         if not (isinstance(self.path[0], moveto_pt) or
1091                 isinstance(self.path[0], arc_pt) or
1092                 isinstance(self.path[0], arcn_pt)):
1093             raise PathException("first path element must be either moveto, arc, or arcn")
1094         # PDF practically only supports normpathels
1095         context = _pathcontext()
1096         for pel in self.path:
1097             for npel in pel._normalized(context):
1098                 npel.outputPDF(file)
1099             pel._updatecontext(context)
1100
1101 ################################################################################
1102 # some special kinds of path, again in two variants
1103 ################################################################################
1104
1105 class line_pt(path):
1106
1107    """straight line from (x1_pt, y1_pt) to (x2_pt, y2_pt) (coordinates in pts)"""
1108
1109    def __init__(self, x1_pt, y1_pt, x2_pt, y2_pt):
1110        path.__init__(self, moveto_pt(x1_pt, y1_pt), lineto_pt(x2_pt, y2_pt))
1111
1112
1113 class curve_pt(path):
1114
1115    """Bezier curve with control points (x0_pt, y1_pt),..., (x3_pt, y3_pt)
1116    (coordinates in pts)"""
1117
1118    def __init__(self, x0_pt, y0_pt, x1_pt, y1_pt, x2_pt, y2_pt, x3_pt, y3_pt):
1119        path.__init__(self,
1120                      moveto_pt(x0_pt, y0_pt),
1121                      curveto_pt(x1_pt, y1_pt, x2_pt, y2_pt, x3_pt, y3_pt))
1122
1123
1124 class rect_pt(path):
1125
1126    """rectangle at position (x,y) with width and height (coordinates in pts)"""
1127
1128    def __init__(self, x, y, width, height):
1129        path.__init__(self, moveto_pt(x, y),
1130                            lineto_pt(x+width, y),
1131                            lineto_pt(x+width, y+height),
1132                            lineto_pt(x, y+height),
1133                            closepath())
1134
1135
1136 class circle_pt(path):
1137
1138    """circle with center (x,y) and radius"""
1139
1140    def __init__(self, x, y, radius):
1141        path.__init__(self, arc_pt(x, y, radius, 0, 360),
1142                      closepath())
1143
1144
1145 class line(line_pt):
1146
1147    """straight line from (x1, y1) to (x2, y2)"""
1148
1149    def __init__(self, x1, y1, x2, y2):
1150        line_pt.__init__(self,
1151                         unit.topt(x1), unit.topt(y1),
1152                         unit.topt(x2), unit.topt(y2))
1153
1154
1155 class curve(curve_pt):
1156
1157    """Bezier curve with control points (x0, y1),..., (x3, y3)"""
1158
1159    def __init__(self, x0, y0, x1, y1, x2, y2, x3, y3):
1160        curve_pt.__init__(self,
1161                          unit.topt(x0), unit.topt(y0),
1162                          unit.topt(x1), unit.topt(y1),
1163                          unit.topt(x2), unit.topt(y2),
1164                          unit.topt(x3), unit.topt(y3))
1165
1166
1167 class rect(rect_pt):
1168
1169    """rectangle at position (x,y) with width and height"""
1170
1171    def __init__(self, x, y, width, height):
1172        rect_pt.__init__(self,
1173                         unit.topt(x), unit.topt(y),
1174                         unit.topt(width), unit.topt(height))
1175
1176
1177 class circle(circle_pt):
1178
1179    """circle with center (x,y) and radius"""
1180
1181    def __init__(self, x, y, radius):
1182        circle_pt.__init__(self,
1183                         unit.topt(x), unit.topt(y),
1184                         unit.topt(radius))
1185
1186 ################################################################################
1187 # normpath and corresponding classes
1188 ################################################################################
1189
1190 # two helper functions for the intersection of normpathels
1191
1192 def _intersectnormcurves(a, a_t0, a_t1, b, b_t0, b_t1, epsilon=1e-5):
1193     """intersect two bpathels
1194
1195     a and b are bpathels with parameter ranges [a_t0, a_t1],
1196     respectively [b_t0, b_t1].
1197     epsilon determines when the bpathels are assumed to be straight
1198
1199     """
1200
1201     # intersection of bboxes is a necessary criterium for intersection
1202     if not a.bbox().intersects(b.bbox()): return []
1203
1204     if not a.isstraight(epsilon):
1205         (aa, ab) = a.midpointsplit()
1206         a_tm = 0.5*(a_t0+a_t1)
1207
1208         if not b.isstraight(epsilon):
1209             (ba, bb) = b.midpointsplit()
1210             b_tm = 0.5*(b_t0+b_t1)
1211
1212             return ( _intersectnormcurves(aa, a_t0, a_tm,
1213                                        ba, b_t0, b_tm, epsilon) +
1214                      _intersectnormcurves(ab, a_tm, a_t1,
1215                                        ba, b_t0, b_tm, epsilon) +
1216                      _intersectnormcurves(aa, a_t0, a_tm,
1217                                        bb, b_tm, b_t1, epsilon) +
1218                      _intersectnormcurves(ab, a_tm, a_t1,
1219                                        bb, b_tm, b_t1, epsilon) )
1220         else:
1221             return ( _intersectnormcurves(aa, a_t0, a_tm,
1222                                        b, b_t0, b_t1, epsilon) +
1223                      _intersectnormcurves(ab, a_tm, a_t1,
1224                                        b, b_t0, b_t1, epsilon) )
1225     else:
1226         if not b.isstraight(epsilon):
1227             (ba, bb) = b.midpointsplit()
1228             b_tm = 0.5*(b_t0+b_t1)
1229
1230             return  ( _intersectnormcurves(a, a_t0, a_t1,
1231                                        ba, b_t0, b_tm, epsilon) +
1232                       _intersectnormcurves(a, a_t0, a_t1,
1233                                        bb, b_tm, b_t1, epsilon) )
1234         else:
1235             # no more subdivisions of either a or b
1236             # => try to intersect a and b as straight line segments
1237
1238             a_deltax = a.x3_pt - a.x0_pt
1239             a_deltay = a.y3_pt - a.y0_pt
1240             b_deltax = b.x3_pt - b.x0_pt
1241             b_deltay = b.y3_pt - b.y0_pt
1242
1243             det = b_deltax*a_deltay - b_deltay*a_deltax
1244
1245             ba_deltax0_pt = b.x0_pt - a.x0_pt
1246             ba_deltay0_pt = b.y0_pt - a.y0_pt
1247
1248             try:
1249                 a_t = ( b_deltax*ba_deltay0_pt - b_deltay*ba_deltax0_pt)/det
1250                 b_t = ( a_deltax*ba_deltay0_pt - a_deltay*ba_deltax0_pt)/det
1251             except ArithmeticError:
1252                 return []
1253
1254             # check for intersections out of bound
1255             if not (0<=a_t<=1 and 0<=b_t<=1): return []
1256
1257             # return rescaled parameters of the intersection
1258             return [ ( a_t0 + a_t * (a_t1 - a_t0),
1259                        b_t0 + b_t * (b_t1 - b_t0) ) ]
1260
1261
1262 def _intersectnormlines(a, b):
1263     """return one-element list constisting either of tuple of
1264     parameters of the intersection point of the two normlines a and b
1265     or empty list if both normlines do not intersect each other"""
1266
1267     a_deltax_pt = a.x1_pt - a.x0_pt
1268     a_deltay_pt = a.y1_pt - a.y0_pt
1269     b_deltax_pt = b.x1_pt - b.x0_pt
1270     b_deltay_pt = b.y1_pt - b.y0_pt
1271
1272     det = b_deltax_pt * a_deltay_pt - b_deltay_pt * a_deltax_pt
1273
1274     ba_deltax0_pt = b.x0_pt - a.x0_pt
1275     ba_deltay0_pt = b.y0_pt - a.y0_pt
1276
1277     try:
1278         a_t = ( b_deltax_pt * ba_deltay0_pt - b_deltay_pt * ba_deltax0_pt)/det
1279         b_t = ( a_deltax_pt * ba_deltay0_pt - a_deltay_pt * ba_deltax0_pt)/det
1280     except ArithmeticError:
1281         return []
1282
1283     # check for intersections out of bound
1284     if not (0<=a_t<=1 and 0<=b_t<=1): return []
1285
1286     # return parameters of the intersection
1287     return [( a_t, b_t)]
1288
1289 #
1290 # normpathel: normalized element
1291 #
1292
1293 class normpathel:
1294
1295     """element of a normalized sub path"""
1296
1297     def at_pt(self, t):
1298         """returns coordinates of point in pts at parameter t (0<=t<=1) """
1299         pass
1300
1301     def arclen_pt(self, epsilon=1e-5):
1302         """returns arc length of normpathel in pts with given accuracy epsilon"""
1303         pass
1304
1305     def _arclentoparam_pt(self, lengths, epsilon=1e-5):
1306         """returns tuple (t,l) with
1307           t the parameter where the arclen of normpathel is length and
1308           l the total arclen
1309
1310         length:  length (in pts) to find the parameter for
1311         epsilon: epsilon controls the accuracy for calculation of the
1312                  length of the Bezier elements
1313         """
1314         # Note: _arclentoparam returns both, parameters and total lengths
1315         # while  arclentoparam returns only parameters
1316         pass
1317
1318     def bbox(self):
1319         """return bounding box of normpathel"""
1320         pass
1321
1322     def curvradius_pt(self, param):
1323         """Returns the curvature radius in pts at parameter param.
1324         This is the inverse of the curvature at this parameter
1325
1326         Please note that this radius can be negative or positive,
1327         depending on the sign of the curvature"""
1328         pass
1329
1330     def intersect(self, other, epsilon=1e-5):
1331         """intersect self with other normpathel"""
1332         pass
1333
1334     def reversed(self):
1335         """return reversed normpathel"""
1336         pass
1337
1338     def split(self, parameters):
1339         """splits normpathel
1340
1341         parameters: list of parameter values (0<=t<=1) at which to split
1342
1343         returns None or list of tuple of normpathels corresponding to
1344         the orginal normpathel.
1345
1346         """
1347         pass
1348
1349     def tangentvector_pt(self, t):
1350         """returns tangent vector of normpathel in pts at parameter t (0<=t<=1)"""
1351         pass
1352
1353     def transformed(self, trafo):
1354         """return transformed normpathel according to trafo"""
1355         pass
1356
1357     def outputPS(self, file):
1358         """write PS code corresponding to normpathel to file"""
1359         pass
1360
1361     def outputPS(self, file):
1362         """write PDF code corresponding to normpathel to file"""
1363         pass
1364
1365 #
1366 # there are only two normpathels: normline and normcurve
1367 #
1368
1369 class normline(normpathel):
1370
1371     """Straight line from (x0_pt, y0_pt) to (x1_pt, y1_pt) (coordinates in pts)"""
1372
1373     __slots__ = "x0_pt", "y0_pt", "x1_pt", "y1_pt"
1374
1375     def __init__(self, x0_pt, y0_pt, x1_pt, y1_pt):
1376          self.x0_pt = x0_pt
1377          self.y0_pt = y0_pt
1378          self.x1_pt = x1_pt
1379          self.y1_pt = y1_pt
1380
1381     def __str__(self):
1382         return "normline(%g, %g, %g, %g)" % (self.x0_pt, self.y0_pt, self.x1_pt, self.y1_pt)
1383
1384     def _arclentoparam_pt(self, lengths, epsilon=1e-5):
1385         l = self.arclen_pt(epsilon)
1386         return ([max(min(1.0 * length / l, 1), 0) for length in lengths], l)
1387
1388     def _normcurve(self):
1389         """ return self as equivalent normcurve """
1390         xa_pt = self.x0_pt+(self.x1_pt-self.x0_pt)/3.0
1391         ya_pt = self.y0_pt+(self.y1_pt-self.y0_pt)/3.0
1392         xb_pt = self.x0_pt+2.0*(self.x1_pt-self.x0_pt)/3.0
1393         yb_pt = self.y0_pt+2.0*(self.y1_pt-self.y0_pt)/3.0
1394         return normcurve(self.x0_pt, self.y0_pt, xa_pt, ya_pt, xb_pt, yb_pt, self.x1_pt, self.y1_pt)
1395
1396     def arclen_pt(self,  epsilon=1e-5):
1397         return math.hypot(self.x0_pt-self.x1_pt, self.y0_pt-self.y1_pt)
1398
1399     def at_pt(self, t):
1400         return self.x0_pt+(self.x1_pt-self.x0_pt)*t, self.y0_pt+(self.y1_pt-self.y0_pt)*t
1401
1402     def bbox(self):
1403         return bbox.bbox_pt(min(self.x0_pt, self.x1_pt), min(self.y0_pt, self.y1_pt),
1404                           max(self.x0_pt, self.x1_pt), max(self.y0_pt, self.y1_pt))
1405
1406     def begin_pt(self):
1407         return self.x0_pt, self.y0_pt
1408
1409     def curvradius_pt(self, param):
1410         return None
1411
1412     def end_pt(self):
1413         return self.x1_pt, self.y1_pt
1414
1415     def intersect(self, other, epsilon=1e-5):
1416         if isinstance(other, normline):
1417             return _intersectnormlines(self, other)
1418         else:
1419             return  _intersectnormcurves(self._normcurve(), 0, 1, other, 0, 1, epsilon)
1420
1421     def isstraight(self, epsilon):
1422         return 1
1423
1424     def reverse(self):
1425         self.x0_pt, self.y0_pt, self.x1_pt, self.y1_pt = self.x1_pt, self.y1_pt, self.x0_pt, self.y0_pt
1426
1427     def reversed(self):
1428         return normline(self.x1_pt, self.y1_pt, self.x0_pt, self.y0_pt)
1429
1430     def split(self, params):
1431         x0_pt, y0_pt = self.x0_pt, self.y0_pt
1432         x1_pt, y1_pt = self.x1_pt, self.y1_pt
1433         if params:
1434             xl_pt, yl_pt = x0_pt, y0_pt
1435             result = []
1436
1437             if params[0] == 0:
1438                 result.append(None)
1439                 params = params[1:]
1440
1441             if params:
1442                 for t in params:
1443                     xs_pt, ys_pt = x0_pt + (x1_pt-x0_pt)*t, y0_pt + (y1_pt-y0_pt)*t
1444                     result.append(normline(xl_pt, yl_pt, xs_pt, ys_pt))
1445                     xl_pt, yl_pt = xs_pt, ys_pt
1446
1447                 if params[-1]!=1:
1448                     result.append(normline(xs_pt, ys_pt, x1_pt, y1_pt))
1449                 else:
1450                     result.append(None)
1451             else:
1452                 result.append(normline(x0_pt, y0_pt, x1_pt, y1_pt))
1453         else:
1454             result = []
1455         return result
1456
1457     def tangentvector_pt(self, param):
1458         return self.x1_pt-self.x0_pt, self.y1_pt-self.y0_pt
1459
1460     def trafo(self, param):
1461         tx_pt, ty_pt = self.at_pt(param)
1462         tdx_pt, tdy_pt = self.x1_pt-self.x0_pt, self.y1_pt-self.y0_pt
1463         return trafo.translate_pt(tx_pt, ty_pt)*trafo.rotate(degrees(math.atan2(tdy_pt, tdx_pt)))
1464
1465     def transformed(self, trafo):
1466         return normline(*(trafo._apply(self.x0_pt, self.y0_pt) + trafo._apply(self.x1_pt, self.y1_pt)))
1467
1468     def outputPS(self, file):
1469         file.write("%g %g lineto\n" % (self.x1_pt, self.y1_pt))
1470
1471     def outputPDF(self, file):
1472         file.write("%f %f l\n" % (self.x1_pt, self.y1_pt))
1473
1474
1475 class normcurve(normpathel):
1476
1477     """Bezier curve with control points x0_pt, y0_pt, x1_pt, y1_pt, x2_pt, y2_pt, x3_pt, y3_pt (coordinates in pts)"""
1478
1479     __slots__ = "x0_pt", "y0_pt", "x1_pt", "y1_pt", "x2_pt", "y2_pt", "x3_pt", "y3_pt"
1480
1481     def __init__(self, x0_pt, y0_pt, x1_pt, y1_pt, x2_pt, y2_pt, x3_pt, y3_pt):
1482         self.x0_pt = x0_pt
1483         self.y0_pt = y0_pt
1484         self.x1_pt = x1_pt
1485         self.y1_pt = y1_pt
1486         self.x2_pt = x2_pt
1487         self.y2_pt = y2_pt
1488         self.x3_pt = x3_pt
1489         self.y3_pt = y3_pt
1490
1491     def __str__(self):
1492         return "normcurve(%g, %g, %g, %g, %g, %g, %g, %g)" % (self.x0_pt, self.y0_pt, self.x1_pt, self.y1_pt,
1493                                                               self.x2_pt, self.y2_pt, self.x3_pt, self.y3_pt)
1494
1495     def _arclentoparam_pt(self, lengths, epsilon=1e-5):
1496         """computes the parameters [t] of bpathel where the given lengths (in pts) are assumed
1497         returns ( [parameters], total arclen)
1498         A negative length gives a parameter 0"""
1499
1500         # create the list of accumulated lengths
1501         # and the length of the parameters
1502         seg = self.seglengths(1, epsilon)
1503         arclens = [seg[i][0] for i in range(len(seg))]
1504         Dparams = [seg[i][1] for i in range(len(seg))]
1505         l = len(arclens)
1506         for i in range(1,l):
1507             arclens[i] += arclens[i-1]
1508
1509         # create the list of parameters to be returned
1510         params = []
1511         for length in lengths:
1512             # find the last index that is smaller than length
1513             try:
1514                 lindex = bisect.bisect_left(arclens, length)
1515             except: # workaround for python 2.0
1516                 lindex = bisect.bisect(arclens, length)
1517                 while lindex and (lindex >= len(arclens) or
1518                                   arclens[lindex] >= length):
1519                     lindex -= 1
1520             if lindex == 0:
1521                 param = Dparams[0] * length * 1.0 / arclens[0]
1522             elif lindex < l-1:
1523                 param = Dparams[lindex+1] * (length - arclens[lindex]) * 1.0 / (arclens[lindex+1] - arclens[lindex])
1524                 for i in range(lindex+1):
1525                     param += Dparams[i]
1526             else:
1527                 param = 1 + Dparams[-1] * (length - arclens[-1]) * 1.0 / (arclens[-1] - arclens[-2])
1528
1529             param = max(min(param,1),0)
1530             params.append(param)
1531         return (params, arclens[-1])
1532
1533     def arclen_pt(self, epsilon=1e-5):
1534         """computes arclen of bpathel in pts using successive midpoint split"""
1535         if self.isstraight(epsilon):
1536             return math.hypot(self.x3_pt-self.x0_pt, self.y3_pt-self.y0_pt)
1537         else:
1538             a, b = self.midpointsplit()
1539             return a.arclen_pt(epsilon) + b.arclen_pt(epsilon)
1540
1541
1542     def at_pt(self, t):
1543         xt_pt = ( (-self.x0_pt+3*self.x1_pt-3*self.x2_pt+self.x3_pt)*t*t*t +
1544                   (3*self.x0_pt-6*self.x1_pt+3*self.x2_pt          )*t*t +
1545                   (-3*self.x0_pt+3*self.x1_pt                      )*t +
1546                   self.x0_pt )
1547         yt_pt = ( (-self.y0_pt+3*self.y1_pt-3*self.y2_pt+self.y3_pt)*t*t*t +
1548                   (3*self.y0_pt-6*self.y1_pt+3*self.y2_pt          )*t*t +
1549                   (-3*self.y0_pt+3*self.y1_pt                      )*t +
1550                   self.y0_pt )
1551         return xt_pt, yt_pt
1552
1553     def bbox(self):
1554         return bbox.bbox_pt(min(self.x0_pt, self.x1_pt, self.x2_pt, self.x3_pt),
1555                           min(self.y0_pt, self.y1_pt, self.y2_pt, self.y3_pt),
1556                           max(self.x0_pt, self.x1_pt, self.x2_pt, self.x3_pt),
1557                           max(self.y0_pt, self.y1_pt, self.y2_pt, self.y3_pt))
1558
1559     def begin_pt(self):
1560         return self.x0_pt, self.y0_pt
1561
1562     def curvradius_pt(self, param):
1563         xdot = ( 3 * (1-param)*(1-param) * (-self.x0_pt + self.x1_pt) +
1564                  6 * (1-param)*param * (-self.x1_pt + self.x2_pt) +
1565                  3 * param*param * (-self.x2_pt + self.x3_pt) )
1566         ydot = ( 3 * (1-param)*(1-param) * (-self.y0_pt + self.y1_pt) +
1567                  6 * (1-param)*param * (-self.y1_pt + self.y2_pt) +
1568                  3 * param*param * (-self.y2_pt + self.y3_pt) )
1569         xddot = ( 6 * (1-param) * (self.x0_pt - 2*self.x1_pt + self.x2_pt) +
1570                   6 * param * (self.x1_pt - 2*self.x2_pt + self.x3_pt) )
1571         yddot = ( 6 * (1-param) * (self.y0_pt - 2*self.y1_pt + self.y2_pt) +
1572                   6 * param * (self.y1_pt - 2*self.y2_pt + self.y3_pt) )
1573         return (xdot**2 + ydot**2)**1.5 / (xdot*yddot - ydot*xddot)
1574
1575     def end_pt(self):
1576         return self.x3_pt, self.y3_pt
1577
1578     def intersect(self, other, epsilon=1e-5):
1579         if isinstance(other, normline):
1580             return  _intersectnormcurves(self, 0, 1, other._normcurve(), 0, 1, epsilon)
1581         else:
1582             return  _intersectnormcurves(self, 0, 1, other, 0, 1, epsilon)
1583
1584     def isstraight(self, epsilon=1e-5):
1585         """check wheter the normcurve is approximately straight"""
1586
1587         # just check, whether the modulus of the difference between
1588         # the length of the control polygon
1589         # (i.e. |P1-P0|+|P2-P1|+|P3-P2|) and the length of the
1590         # straight line between starting and ending point of the
1591         # normcurve (i.e. |P3-P1|) is smaller the epsilon
1592         return abs(math.hypot(self.x1_pt-self.x0_pt, self.y1_pt-self.y0_pt)+
1593                    math.hypot(self.x2_pt-self.x1_pt, self.y2_pt-self.y1_pt)+
1594                    math.hypot(self.x3_pt-self.x2_pt, self.y3_pt-self.y2_pt)-
1595                    math.hypot(self.x3_pt-self.x0_pt, self.y3_pt-self.y0_pt))<epsilon
1596
1597     def midpointsplit(self):
1598         """splits bpathel at midpoint returning bpath with two bpathels"""
1599
1600         # for efficiency reason, we do not use self.split(0.5)!
1601
1602         # first, we have to calculate the  midpoints between adjacent
1603         # control points
1604         x01_pt = 0.5*(self.x0_pt + self.x1_pt)
1605         y01_pt = 0.5*(self.y0_pt + self.y1_pt)
1606         x12_pt = 0.5*(self.x1_pt + self.x2_pt)
1607         y12_pt = 0.5*(self.y1_pt + self.y2_pt)
1608         x23_pt = 0.5*(self.x2_pt + self.x3_pt)
1609         y23_pt = 0.5*(self.y2_pt + self.y3_pt)
1610
1611         # In the next iterative step, we need the midpoints between 01 and 12
1612         # and between 12 and 23
1613         x01_12_pt = 0.5*(x01_pt + x12_pt)
1614         y01_12_pt = 0.5*(y01_pt + y12_pt)
1615         x12_23_pt = 0.5*(x12_pt + x23_pt)
1616         y12_23_pt = 0.5*(y12_pt + y23_pt)
1617
1618         # Finally the midpoint is given by
1619         xmidpoint_pt = 0.5*(x01_12_pt + x12_23_pt)
1620         ymidpoint_pt = 0.5*(y01_12_pt + y12_23_pt)
1621
1622         return (normcurve(self.x0_pt, self.y0_pt,
1623                           x01_pt, y01_pt,
1624                           x01_12_pt, y01_12_pt,
1625                           xmidpoint_pt, ymidpoint_pt),
1626                 normcurve(xmidpoint_pt, ymidpoint_pt,
1627                           x12_23_pt, y12_23_pt,
1628                           x23_pt, y23_pt,
1629                           self.x3_pt, self.y3_pt))
1630
1631     def reverse(self):
1632         self.x0_pt, self.y0_pt, self.x1_pt, self.y1_pt, self.x2_pt, self.y2_pt, self.x3_pt, self.y3_pt = \
1633         self.x3_pt, self.y3_pt, self.x2_pt, self.y2_pt, self.x1_pt, self.y1_pt, self.x0_pt, self.y0_pt
1634
1635     def reversed(self):
1636         return normcurve(self.x3_pt, self.y3_pt, self.x2_pt, self.y2_pt, self.x1_pt, self.y1_pt, self.x0_pt, self.y0_pt)
1637
1638     def seglengths(self, paraminterval, epsilon=1e-5):
1639         """returns the list of segment line lengths (in pts) of the normcurve
1640            together with the length of the parameterinterval"""
1641
1642         # lower and upper bounds for the arclen
1643         lowerlen = math.hypot(self.x3_pt-self.x0_pt, self.y3_pt-self.y0_pt)
1644         upperlen = ( math.hypot(self.x1_pt-self.x0_pt, self.y1_pt-self.y0_pt) +
1645                      math.hypot(self.x2_pt-self.x1_pt, self.y2_pt-self.y1_pt) +
1646                      math.hypot(self.x3_pt-self.x2_pt, self.y3_pt-self.y2_pt) )
1647
1648         # instead of isstraight method:
1649         if abs(upperlen-lowerlen)<epsilon:
1650             return [( 0.5*(upperlen+lowerlen), paraminterval )]
1651         else:
1652             a, b = self.midpointsplit()
1653             return a.seglengths(0.5*paraminterval, epsilon) + b.seglengths(0.5*paraminterval, epsilon)
1654
1655     def _split(self, parameters):
1656         """return list of normcurve corresponding to split at parameters"""
1657
1658         # first, we calculate the coefficients corresponding to our
1659         # original bezier curve. These represent a useful starting
1660         # point for the following change of the polynomial parameter
1661         a0x_pt = self.x0_pt
1662         a0y_pt = self.y0_pt
1663         a1x_pt = 3*(-self.x0_pt+self.x1_pt)
1664         a1y_pt = 3*(-self.y0_pt+self.y1_pt)
1665         a2x_pt = 3*(self.x0_pt-2*self.x1_pt+self.x2_pt)
1666         a2y_pt = 3*(self.y0_pt-2*self.y1_pt+self.y2_pt)
1667         a3x_pt = -self.x0_pt+3*(self.x1_pt-self.x2_pt)+self.x3_pt
1668         a3y_pt = -self.y0_pt+3*(self.y1_pt-self.y2_pt)+self.y3_pt
1669
1670         if parameters[0]!=0:
1671             parameters = [0] + parameters
1672         if parameters[-1]!=1:
1673             parameters = parameters + [1]
1674
1675         result = []
1676
1677         for i in range(len(parameters)-1):
1678             t1 = parameters[i]
1679             dt = parameters[i+1]-t1
1680
1681             # [t1,t2] part
1682             #
1683             # the new coefficients of the [t1,t1+dt] part of the bezier curve
1684             # are then given by expanding
1685             #  a0 + a1*(t1+dt*u) + a2*(t1+dt*u)**2 +
1686             #  a3*(t1+dt*u)**3 in u, yielding
1687             #
1688             #   a0 + a1*t1 + a2*t1**2 + a3*t1**3        +
1689             #   ( a1 + 2*a2 + 3*a3*t1**2 )*dt    * u    +
1690             #   ( a2 + 3*a3*t1 )*dt**2           * u**2 +
1691             #   a3*dt**3                         * u**3
1692             #
1693             # from this values we obtain the new control points by inversion
1694             #
1695             # XXX: we could do this more efficiently by reusing for
1696             # (x0_pt, y0_pt) the control point (x3_pt, y3_pt) from the previous
1697             # Bezier curve
1698
1699             x0_pt = a0x_pt + a1x_pt*t1 + a2x_pt*t1*t1 + a3x_pt*t1*t1*t1
1700             y0_pt = a0y_pt + a1y_pt*t1 + a2y_pt*t1*t1 + a3y_pt*t1*t1*t1
1701             x1_pt = (a1x_pt+2*a2x_pt*t1+3*a3x_pt*t1*t1)*dt/3.0 + x0_pt
1702             y1_pt = (a1y_pt+2*a2y_pt*t1+3*a3y_pt*t1*t1)*dt/3.0 + y0_pt
1703             x2_pt = (a2x_pt+3*a3x_pt*t1)*dt*dt/3.0 - x0_pt + 2*x1_pt
1704             y2_pt = (a2y_pt+3*a3y_pt*t1)*dt*dt/3.0 - y0_pt + 2*y1_pt
1705             x3_pt = a3x_pt*dt*dt*dt + x0_pt - 3*x1_pt + 3*x2_pt
1706             y3_pt = a3y_pt*dt*dt*dt + y0_pt - 3*y1_pt + 3*y2_pt
1707
1708             result.append(normcurve(x0_pt, y0_pt, x1_pt, y1_pt, x2_pt, y2_pt, x3_pt, y3_pt))
1709
1710         return result
1711
1712     def split(self, params):
1713         if params:
1714             # we need to split
1715             bps = self._split(list(params))
1716
1717             if params[0]==0:
1718                 result = [None]
1719             else:
1720                 bp0 = bps[0]
1721                 result = [normcurve(self.x0_pt, self.y0_pt, bp0.x1_pt, bp0.y1_pt, bp0.x2_pt, bp0.y2_pt, bp0.x3_pt, bp0.y3_pt)]
1722                 bps = bps[1:]
1723
1724             for bp in bps:
1725                 result.append(normcurve(bp.x0_pt, bp.y0_pt, bp.x1_pt, bp.y1_pt, bp.x2_pt, bp.y2_pt, bp.x3_pt, bp.y3_pt))
1726
1727             if params[-1]==1:
1728                 result.append(None)
1729         else:
1730             result = []
1731         return result
1732
1733     def tangentvector_pt(self, param):
1734         tvectx = (3*(  -self.x0_pt+3*self.x1_pt-3*self.x2_pt+self.x3_pt)*param*param +
1735                   2*( 3*self.x0_pt-6*self.x1_pt+3*self.x2_pt        )*param +
1736                     (-3*self.x0_pt+3*self.x1_pt                  ))
1737         tvecty = (3*(  -self.y0_pt+3*self.y1_pt-3*self.y2_pt+self.y3_pt)*param*param +
1738                   2*( 3*self.y0_pt-6*self.y1_pt+3*self.y2_pt        )*param +
1739                     (-3*self.y0_pt+3*self.y1_pt                  ))
1740         return (tvectx, tvecty)
1741
1742     def trafo(self, param):
1743         tx_pt, ty_pt = self.at_pt(param)
1744         tdx_pt, tdy_pt = self.tangentvector_pt(param)
1745         return trafo.translate_pt(tx_pt, ty_pt)*trafo.rotate(degrees(math.atan2(tdy_pt, tdx_pt)))
1746
1747     def transform(self, trafo):
1748         self.x0_pt, self.y0_pt = trafo._apply(self.x0_pt, self.y0_pt)
1749         self.x1_pt, self.y1_pt = trafo._apply(self.x1_pt, self.y1_pt)
1750         self.x2_pt, self.y2_pt = trafo._apply(self.x2_pt, self.y2_pt)
1751         self.x3_pt, self.y3_pt = trafo._apply(self.x3_pt, self.y3_pt)
1752
1753     def transformed(self, trafo):
1754         return normcurve(*(trafo._apply(self.x0_pt, self.y0_pt)+
1755                            trafo._apply(self.x1_pt, self.y1_pt)+
1756                            trafo._apply(self.x2_pt, self.y2_pt)+
1757                            trafo._apply(self.x3_pt, self.y3_pt)))
1758
1759     def outputPS(self, file):
1760         file.write("%g %g %g %g %g %g curveto\n" % (self.x1_pt, self.y1_pt, self.x2_pt, self.y2_pt, self.x3_pt, self.y3_pt))
1761
1762     def outputPDF(self, file):
1763         file.write("%f %f %f %f %f %f c\n" % (self.x1_pt, self.y1_pt, self.x2_pt, self.y2_pt, self.x3_pt, self.y3_pt))
1764
1765 #
1766 # normpaths are made up of normsubpaths, which represent connected line segments
1767 #
1768
1769 class normsubpath:
1770
1771     """sub path of a normalized path
1772
1773     A subpath consists of a list of normpathels, i.e., lines and bcurves
1774     and can either be closed or not.
1775
1776     Some invariants, which have to be obeyed:
1777     - All normpathels have to be longer than epsilon pts.
1778     - The last point of a normpathel and the first point of the next
1779     element have to be equal.
1780     - When the path is closed, the last normpathel has to be a
1781     normline and the last point of this normline has to be equal
1782     to the first point of the first normpathel, except when
1783     this normline would be too short.
1784     """
1785
1786     __slots__ = "normpathels", "closed", "epsilon"
1787
1788     def __init__(self, normpathels, closed, epsilon=1e-5):
1789         self.normpathels = [npel for npel in normpathels if not npel.isstraight(epsilon) or npel.arclen_pt(epsilon)>epsilon]
1790         self.closed = closed
1791         self.epsilon = epsilon
1792
1793     def __str__(self):
1794         return "subpath(%s, [%s])" % (self.closed and "closed" or "open",
1795                                     ", ".join(map(str, self.normpathels)))
1796
1797     def arclen_pt(self):
1798         """returns total arc length of normsubpath in pts with accuracy epsilon"""
1799         return sum([npel.arclen_pt(self.epsilon) for npel in self.normpathels])
1800
1801     def _arclentoparam_pt(self, lengths):
1802         """returns [t, l] where t are parameter value(s) matching given length(s)
1803         and l is the total length of the normsubpath
1804         The parameters are with respect to the normsubpath: t in [0, self.range()]
1805         lengths that are < 0 give parameter 0"""
1806
1807         allarclen = 0
1808         allparams = [0] * len(lengths)
1809         rests = copy.copy(lengths)
1810
1811         for pel in self.normpathels:
1812             params, arclen = pel._arclentoparam_pt(rests, self.epsilon)
1813             allarclen += arclen
1814             for i in range(len(rests)):
1815                 if rests[i] >= 0:
1816                     rests[i] -= arclen
1817                     allparams[i] += params[i]
1818
1819         return (allparams, allarclen)
1820
1821     def at_pt(self, param):
1822         """return coordinates in pts of sub path at parameter value param
1823
1824         The parameter param must be smaller or equal to the number of
1825         segments in the normpath, otherwise None is returned.
1826         """
1827         try:
1828             return self.normpathels[int(param-self.epsilon)].at_pt(param-int(param-self.epsilon))
1829         except:
1830             raise PathException("parameter value param out of range")
1831
1832     def bbox(self):
1833         if self.normpathels:
1834             abbox = self.normpathels[0].bbox()
1835             for anormpathel in self.normpathels[1:]:
1836                 abbox += anormpathel.bbox()
1837             return abbox
1838         else:
1839             return None
1840
1841     def begin_pt(self):
1842         return self.normpathels[0].begin_pt()
1843
1844     def curvradius_pt(self, param):
1845         try:
1846             return self.normpathels[int(param-self.epsilon)].curvradius_pt(param-int(param-self.epsilon))
1847         except:
1848             raise PathException("parameter value param out of range")
1849
1850     def end_pt(self):
1851         return self.normpathels[-1].end_pt()
1852
1853     def intersect(self, other):
1854         """intersect self with other normsubpath
1855
1856         returns a tuple of lists consisting of the parameter values
1857         of the intersection points of the corresponding normsubpath
1858
1859         """
1860         intersections = ([], [])
1861         epsilon = min(self.epsilon, other.epsilon)
1862         # Intersect all subpaths of self with the subpaths of other
1863         for t_a, pel_a  in enumerate(self.normpathels):
1864             for t_b, pel_b in enumerate(other.normpathels):
1865                 for intersection in pel_a.intersect(pel_b, epsilon):
1866                     # check whether an intersection occurs at the end
1867                     # of a closed subpath. If yes, we don't include it
1868                     # in the list of intersections to prevent a
1869                     # duplication of intersection points
1870                     if not ((self.closed and self.range()-intersection[0]-t_a<epsilon) or
1871                             (other.closed and other.range()-intersection[1]-t_b<epsilon)):
1872                         intersections[0].append(intersection[0]+t_a)
1873                         intersections[1].append(intersection[1]+t_b)
1874         return intersections
1875
1876     def range(self):
1877         """return maximal parameter value, i.e. number of line/curve segments"""
1878         return len(self.normpathels)
1879
1880     def reverse(self):
1881         self.normpathels.reverse()
1882         for npel in self.normpathels:
1883             npel.reverse()
1884
1885     def reversed(self):
1886         nnormpathels = []
1887         for i in range(len(self.normpathels)):
1888             nnormpathels.append(self.normpathels[-(i+1)].reversed())
1889         return normsubpath(nnormpathels, self.closed)
1890
1891     def split(self, params):
1892         """split normsubpath at list of parameter values params and return list
1893         of normsubpaths
1894
1895         The parameter list params has to be sorted. Note that each element of
1896         the resulting list is an open normsubpath.
1897         """
1898
1899         if min(params) < -self.epsilon or max(params) > self.range()+self.epsilon:
1900             raise PathException("parameter for split of subpath out of range")
1901
1902         result = []
1903         npels = None
1904         for t, pel in enumerate(self.normpathels):
1905             # determine list of splitting parameters relevant for pel
1906             nparams = []
1907             for nt in params:
1908                 if t+1 >= nt:
1909                     nparams.append(nt-t)
1910                     params = params[1:]
1911
1912             # now we split the path at the filtered parameter values
1913             # This yields a list of normpathels and possibly empty
1914             # segments marked by None
1915             splitresult = pel.split(nparams)
1916             if splitresult:
1917                 # first split?
1918                 if npels is None:
1919                     if splitresult[0] is None:
1920                         # mark split at the beginning of the normsubpath
1921                         result = [None]
1922                     else:
1923                         result.append(normsubpath([splitresult[0]], 0))
1924                 else:
1925                     npels.append(splitresult[0])
1926                     result.append(normsubpath(npels, 0))
1927                 for npel in splitresult[1:-1]:
1928                     result.append(normsubpath([npel], 0))
1929                 if len(splitresult)>1 and splitresult[-1] is not None:
1930                     npels = [splitresult[-1]]
1931                 else:
1932                     npels = []
1933             else:
1934                 if npels is None:
1935                     npels = [pel]
1936                 else:
1937                     npels.append(pel)
1938
1939         if npels:
1940             result.append(normsubpath(npels, 0))
1941         else:
1942             # mark split at the end of the normsubpath
1943             result.append(None)
1944
1945         # join last and first segment together if the normsubpath was originally closed
1946         if self.closed:
1947             if result[0] is None:
1948                 result = result[1:]
1949             elif result[-1] is None:
1950                 result = result[:-1]
1951             else:
1952                 result[-1].normpathels.extend(result[0].normpathels)
1953                 result = result[1:]
1954         return result
1955
1956     def tangent(self, param, length=None):
1957         tx_pt, ty_pt = self.at_pt(param)
1958         try:
1959             tdx_pt, tdy_pt = self.normpathels[int(param-self.epsilon)].tangentvector_pt(param-int(param-self.epsilon))
1960         except:
1961             raise PathException("parameter value param out of range")
1962         tlen = math.hypot(tdx_pt, tdy_pt)
1963         if not (length is None or tlen==0):
1964             sfactor = unit.topt(length)/tlen
1965             tdx_pt *= sfactor
1966             tdy_pt *= sfactor
1967         return line_pt(tx_pt, ty_pt, tx_pt+tdx_pt, ty_pt+tdy_pt)
1968
1969     def trafo(self, param):
1970         try:
1971             return self.normpathels[int(param-self.epsilon)].trafo(param-int(param-self.epsilon))
1972         except:
1973             raise PathException("parameter value param out of range")
1974
1975     def transform(self, trafo):
1976         """transform sub path according to trafo"""
1977         for pel in self.normpathels:
1978             pel.transform(trafo)
1979
1980     def transformed(self, trafo):
1981         """return sub path transformed according to trafo"""
1982         nnormpathels = []
1983         for pel in self.normpathels:
1984             nnormpathels.append(pel.transformed(trafo))
1985         return normsubpath(nnormpathels, self.closed)
1986
1987     def outputPS(self, file):
1988         # if the normsubpath is closed, we must not output a normline at
1989         # the end
1990         if not self.normpathels:
1991             return
1992         if self.closed and isinstance(self.normpathels[-1], normline):
1993             normpathels = self.normpathels[:-1]
1994         else:
1995             normpathels = self.normpathels
1996         if normpathels:
1997             file.write("%g %g moveto\n" % self.begin_pt())
1998             for anormpathel in normpathels:
1999                 anormpathel.outputPS(file)
2000         if self.closed:
2001             file.write("closepath\n")
2002
2003     def outputPDF(self, file):
2004         # if the normsubpath is closed, we must not output a normline at
2005         # the end
2006         if not self.normpathels:
2007             return
2008         if self.closed and isinstance(self.normpathels[-1], normline):
2009             normpathels = self.normpathels[:-1]
2010         else:
2011             normpathels = self.normpathels
2012         if normpathels:
2013             file.write("%f %f m\n" % self.begin_pt())
2014             for anormpathel in normpathels:
2015                 anormpathel.outputPDF(file)
2016         if self.closed:
2017             file.write("h\n")
2018
2019 #
2020 # the normpath class
2021 #
2022
2023 class normpath(path):
2024
2025     """normalized path
2026
2027     A normalized path consists of a list of normalized sub paths.
2028
2029     """
2030
2031     def __init__(self, arg=[], epsilon=1e-5):
2032         """ construct a normpath from another normpath passed as arg,
2033         a path or a list of normsubpaths. An accuracy of epsilon pts
2034         is used for numerical calculations.
2035         """
2036
2037         self.epsilon = epsilon
2038         if isinstance(arg, normpath):
2039             self.subpaths = copy.copy(arg.subpaths)
2040             return
2041         elif isinstance(arg, path):
2042             # split path in sub paths
2043             self.subpaths = []
2044             currentsubpathels = []
2045             context = _pathcontext()
2046             for pel in arg.path:
2047                 for npel in pel._normalized(context):
2048                     if isinstance(npel, moveto_pt):
2049                         if currentsubpathels:
2050                             # append open sub path
2051                             self.subpaths.append(normsubpath(currentsubpathels, 0, epsilon))
2052                         # start new sub path
2053                         currentsubpathels = []
2054                     elif isinstance(npel, closepath):
2055                         if currentsubpathels:
2056                             # append closed sub path
2057                             currentsubpathels.append(normline(context.currentpoint[0], context.currentpoint[1],
2058                                                               context.currentsubpath[0], context.currentsubpath[1]))
2059                         self.subpaths.append(normsubpath(currentsubpathels, 1, epsilon))
2060                         currentsubpathels = []
2061                     else:
2062                         currentsubpathels.append(npel)
2063                 pel._updatecontext(context)
2064
2065             if currentsubpathels:
2066                 # append open sub path
2067                 self.subpaths.append(normsubpath(currentsubpathels, 0, epsilon))
2068         else:
2069             # we expect a list of normsubpaths
2070             self.subpaths = list(arg)
2071
2072     def __add__(self, other):
2073         result = normpath(other)
2074         result.subpaths = self.subpaths + result.subpaths
2075         return result
2076
2077     def __iadd__(self, other):
2078         self.subpaths += normpath(other).subpaths
2079         return self
2080
2081     def __nonzero__(self):
2082         return len(self.subpaths)>0
2083
2084     def __str__(self):
2085         return "normpath(%s)" % ", ".join(map(str, self.subpaths))
2086
2087     def _findsubpath(self, param, arclen):
2088         """return a tuple (subpath, rparam), where subpath is the subpath
2089         containing the position specified by either param or arclen and rparam
2090         is the corresponding parameter value in this subpath.
2091         """
2092
2093         if param is not None and arclen is not None:
2094             raise PathException("either param or arclen has to be specified, but not both")
2095         elif arclen is not None:
2096             param = self.arclentoparam(arclen)
2097
2098         spt = 0
2099         for sp in self.subpaths:
2100             sprange = sp.range()
2101             if spt <= param <= sprange+spt+self.epsilon:
2102                 return sp, param-spt
2103             spt += sprange
2104         raise PathException("parameter value out of range")
2105
2106     def append(self, pathel):
2107         # XXX factor parts of this code out
2108         if self.subpaths[-1].closed:
2109             context = _pathcontext(self.end_pt(), None)
2110             currentsubpathels = []
2111         else:
2112             context = _pathcontext(self.end_pt(), self.subpaths[-1].begin_pt())
2113             currentsubpathels = self.subpaths[-1].normpathels
2114             self.subpaths = self.subpaths[:-1]
2115         for npel in pathel._normalized(context):
2116             if isinstance(npel, moveto_pt):
2117                 if currentsubpathels:
2118                     # append open sub path
2119                     self.subpaths.append(normsubpath(currentsubpathels, 0, self.epsilon))
2120                 # start new sub path
2121                 currentsubpathels = []
2122             elif isinstance(npel, closepath):
2123                 if currentsubpathels:
2124                     # append closed sub path
2125                     currentsubpathels.append(normline(context.currentpoint[0], context.currentpoint[1],
2126                                                       context.currentsubpath[0], context.currentsubpath[1]))
2127                     self.subpaths.append(normsubpath(currentsubpathels, 1, self.epsilon))
2128                 currentsubpathels = []
2129             else:
2130                 currentsubpathels.append(npel)
2131
2132         if currentsubpathels:
2133             # append open sub path
2134             self.subpaths.append(normsubpath(currentsubpathels, 0, self.epsilon))
2135
2136     def arclen_pt(self):
2137         """returns total arc length of normpath in pts"""
2138         return sum([sp.arclen_pt() for sp in self.subpaths])
2139
2140     def arclen(self):
2141         """returns total arc length of normpath"""
2142         return self.arclen_pt() * unit.t_pt
2143
2144     def arclentoparam_pt(self, lengths):
2145         rests = copy.copy(lengths)
2146         allparams = [0] * len(lengths)
2147
2148         for sp in self.subpaths:
2149             # we need arclen for knowing when all the parameters are done
2150             # for lengths that are done: rests[i] is negative
2151             # sp._arclentoparam has to ignore such lengths
2152             params, arclen = sp._arclentoparam_pt(rests)
2153             finis = 0 # number of lengths that are done
2154             for i in range(len(rests)):
2155                 if rests[i] >= 0:
2156                   rests[i] -= arclen
2157                   allparams[i] += params[i]
2158                 else:
2159                     finis += 1
2160             if finis == len(rests): break
2161
2162         if len(lengths) == 1: allparams = allparams[0]
2163         return allparams
2164
2165     def arclentoparam(self, lengths):
2166         """returns the parameter value(s) matching the given length(s)
2167
2168         all given lengths must be positive.
2169         A length greater than the total arclength will give self.range()
2170         """
2171         l = [unit.topt(length) for length in helper.ensuresequence(lengths)]
2172         return self.arclentoparam_pt(l)
2173
2174     def at_pt(self, param=None, arclen=None):
2175         """return coordinates in pts of path at either parameter value param
2176         or arc length arclen.
2177
2178         At discontinuities in the path, the limit from below is returned.
2179         """
2180         sp, param = self._findsubpath(param, arclen)
2181         return sp.at_pt(param)
2182
2183     def at(self, param=None, arclen=None):
2184         """return coordinates of path at either parameter value param
2185         or arc length arclen.
2186
2187         At discontinuities in the path, the limit from below is returned
2188         """
2189         x, y = self.at_pt(param, arclen)
2190         return x * unit.t_pt, y * unit.t_pt
2191
2192     def bbox(self):
2193         abbox = None
2194         for sp in self.subpaths:
2195             nbbox =  sp.bbox()
2196             if abbox is None:
2197                 abbox = nbbox
2198             elif nbbox:
2199                 abbox += nbbox
2200         return abbox
2201
2202     def begin_pt(self):
2203         """return coordinates of first point of first subpath in path (in pts)"""
2204         if self.subpaths:
2205             return self.subpaths[0].begin_pt()
2206         else:
2207             raise PathException("cannot return first point of empty path")
2208
2209     def begin(self):
2210         """return coordinates of first point of first subpath in path"""
2211         x_pt, y_pt = self.begin_pt()
2212         return x_pt * unit.t_pt, y_pt * unit.t_pt
2213
2214     def curvradius_pt(self, param=None, arclen=None):
2215         """Returns the curvature radius in pts (or None if infinite)
2216         at parameter param or arc length arclen.  This is the inverse
2217         of the curvature at this parameter
2218
2219         Please note that this radius can be negative or positive,
2220         depending on the sign of the curvature"""
2221         sp, param = self._findsubpath(param, arclen)
2222         return sp.curvradius_pt(param)
2223
2224     def curvradius(self, param=None, arclen=None):
2225         """Returns the curvature radius (or None if infinite) at
2226         parameter param or arc length arclen.  This is the inverse of
2227         the curvature at this parameter
2228
2229         Please note that this radius can be negative or positive,
2230         depending on the sign of the curvature"""
2231         radius = self.curvradius_pt(param, arclen)
2232         if radius is not None:
2233             radius = radius * unit.t_pt
2234         return radius
2235
2236     def end_pt(self):
2237         """return coordinates of last point of last subpath in path (in pts)"""
2238         if self.subpaths:
2239             return self.subpaths[-1].end_pt()
2240         else:
2241             raise PathException("cannot return last point of empty path")
2242
2243     def end(self):
2244         """return coordinates of last point of last subpath in path"""
2245         x_pt, y_pt = self.end_pt()
2246         return x_pt * unit.t_pt, y_pt * unit.t_pt
2247
2248     def join(self, other):
2249         if not self.subpaths:
2250             raise PathException("cannot join to end of empty path")
2251         if self.subpaths[-1].closed:
2252             raise PathException("cannot join to end of closed sub path")
2253         other = normpath(other)
2254         if not other.subpaths:
2255             raise PathException("cannot join empty path")
2256
2257         self.subpaths[-1].normpathels += other.subpaths[0].normpathels
2258         self.subpaths += other.subpaths[1:]
2259
2260     def joined(self, other):
2261         result = normpath(self.subpaths)
2262         result.join(other)
2263         return result
2264
2265     def intersect(self, other):
2266         """intersect self with other path
2267
2268         returns a tuple of lists consisting of the parameter values
2269         of the intersection points of the corresponding normpath
2270
2271         """
2272         if not isinstance(other, normpath):
2273             other = normpath(other)
2274
2275         # here we build up the result
2276         intersections = ([], [])
2277
2278         # Intersect all subpaths of self with the subpaths of
2279         # other. Here, st_a, st_b are the parameter values
2280         # corresponding to the first point of the subpaths sp_a and
2281         # sp_b, respectively.
2282         st_a = 0
2283         for sp_a in self.subpaths:
2284             st_b =0
2285             for sp_b in other.subpaths:
2286                 for intersection in zip(*sp_a.intersect(sp_b)):
2287                     intersections[0].append(intersection[0]+st_a)
2288                     intersections[1].append(intersection[1]+st_b)
2289                 st_b += sp_b.range()
2290             st_a += sp_a.range()
2291         return intersections
2292
2293     def range(self):
2294         """return maximal value for parameter value param"""
2295         return sum([sp.range() for sp in self.subpaths])
2296
2297     def reverse(self):
2298         """reverse path"""
2299         self.subpaths.reverse()
2300         for sp in self.subpaths:
2301             sp.reverse()
2302
2303     def reversed(self):
2304         """return reversed path"""
2305         nnormpath = normpath()
2306         for i in range(len(self.subpaths)):
2307             nnormpath.subpaths.append(self.subpaths[-(i+1)].reversed())
2308         return nnormpath
2309
2310     def split(self, params):
2311         """split path at parameter values params
2312
2313         Note that the parameter list has to be sorted.
2314
2315         """
2316
2317         # check whether parameter list is really sorted
2318         sortedparams = list(params)
2319         sortedparams.sort()
2320         if sortedparams!=list(params):
2321             raise ValueError("split parameter list params has to be sorted")
2322
2323         # we construct this list of normpaths
2324         result = []
2325
2326         # the currently built up normpath
2327         np = normpath()
2328
2329         t0 = 0
2330         for subpath in self.subpaths:
2331             tf = t0+subpath.range()
2332             if params and tf>=params[0]:
2333                 # split this subpath
2334                 # determine the relevant splitting params
2335                 for i in range(len(params)):
2336                     if params[i]>tf: break
2337                 else:
2338                     i = len(params)
2339
2340                 splitsubpaths = subpath.split([x-t0 for x in params[:i]])
2341                 # handle first element, which may be None, separately
2342                 if splitsubpaths[0] is None:
2343                     if not np.subpaths:
2344                         result.append(None)
2345                     else:
2346                         result.append(np)
2347                         np = normpath()
2348                     splitsubpaths.pop(0)
2349
2350                 for sp in splitsubpaths[:-1]:
2351                     np.subpaths.append(sp)
2352                     result.append(np)
2353                     np = normpath()
2354
2355                 # handle last element which may be None, separately
2356                 if splitsubpaths:
2357                     if splitsubpaths[-1] is None:
2358                         if np.subpaths:
2359                             result.append(np)
2360                             np = normpath()
2361                     else:
2362                         np.subpaths.append(splitsubpaths[-1])
2363
2364                 params = params[i:]
2365             else:
2366                 # append whole subpath to current normpath
2367                 np.subpaths.append(subpath)
2368             t0 = tf
2369
2370         if np.subpaths:
2371             result.append(np)
2372         else:
2373             # mark split at the end of the normsubpath
2374             result.append(None)
2375
2376         return result
2377
2378     def tangent(self, param=None, arclen=None, length=None):
2379         """return tangent vector of path at either parameter value param
2380         or arc length arclen.
2381
2382         At discontinuities in the path, the limit from below is returned.
2383         If length is not None, the tangent vector will be scaled to
2384         the desired length.
2385         """
2386         sp, param = self._findsubpath(param, arclen)
2387         return sp.tangent(param, length)
2388
2389     def transform(self, trafo):
2390         """transform path according to trafo"""
2391         for sp in self.subpaths:
2392             sp.transform(trafo)
2393
2394     def transformed(self, trafo):
2395         """return path transformed according to trafo"""
2396         return normpath([sp.transformed(trafo) for sp in self.subpaths])
2397
2398     def trafo(self, param=None, arclen=None):
2399         """return transformation at either parameter value param or arc length arclen"""
2400         sp, param = self._findsubpath(param, arclen)
2401         return sp.trafo(param)
2402
2403     def outputPS(self, file):
2404         for sp in self.subpaths:
2405             sp.outputPS(file)
2406
2407     def outputPDF(self, file):
2408         for sp in self.subpaths:
2409             sp.outputPDF(file)