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