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         cumlengths = self.seglengths(1, epsilon)
1565         l = len(cumlengths)
1566         parlengths = [cumlengths[i][1] for i in range(l)]
1567         cumlengths[0] = cumlengths[0][0]
1568         for i in range(1,l):
1569             cumlengths[i] = cumlengths[i][0] + cumlengths[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(cumlengths, length)
1577             except: # workaround for python 2.0
1578                 lindex = bisect.bisect(cumlengths, length)
1579                 while lindex and (lindex >= len(cumlengths) or
1580                                   cumlengths[lindex] >= length):
1581                     lindex -= 1
1582             if lindex == 0:
1583                 param = length * 1.0 / cumlengths[0]
1584                 param *= parlengths[0]
1585             elif lindex >= l-2:
1586                 param = 1
1587             else:
1588                 param = (length - cumlengths[lindex]) * 1.0 / (cumlengths[lindex+1] - cumlengths[lindex])
1589                 param *= parlengths[lindex+1]
1590                 for i in range(lindex+1):
1591                     param += parlengths[i]
1592             param = max(min(param,1),0)
1593             params.append(param)
1594         return (params, cumlengths[-1])
1595
1596     def arclen_pt(self, epsilon=1e-5):
1597         """computes arclen of bpathel in pts using successive midpoint split"""
1598         if self.isstraight(epsilon):
1599             return math.hypot(self.x3-self.x0, self.y3-self.y0)
1600         else:
1601             (a, b) = self.midpointsplit()
1602             return a.arclen_pt(epsilon) + b.arclen_pt(epsilon)
1603
1604
1605     def at_pt(self, t):
1606         xt = (  (-self.x0+3*self.x1-3*self.x2+self.x3)*t*t*t +
1607                (3*self.x0-6*self.x1+3*self.x2        )*t*t +
1608               (-3*self.x0+3*self.x1                  )*t +
1609               self.x0)
1610         yt = (  (-self.y0+3*self.y1-3*self.y2+self.y3)*t*t*t +
1611                (3*self.y0-6*self.y1+3*self.y2        )*t*t +
1612               (-3*self.y0+3*self.y1                  )*t +
1613               self.y0)
1614         return (xt, yt)
1615
1616     def bbox(self):
1617         return bbox._bbox(min(self.x0, self.x1, self.x2, self.x3),
1618                           min(self.y0, self.y1, self.y2, self.y3),
1619                           max(self.x0, self.x1, self.x2, self.x3),
1620                           max(self.y0, self.y1, self.y2, self.y3))
1621
1622     def begin_pt(self):
1623         return self.x0, self.y0
1624
1625     def curvradius_pt(self, param):
1626         xdot = 3 * (1-param)*(1-param) * (-self.x0 + self.x1) \
1627              + 6 * (1-param)*param * (-self.x1 + self.x2) \
1628              + 3 * param*param * (-self.x2 + self.x3)
1629         ydot = 3 * (1-param)*(1-param) * (-self.y0 + self.y1) \
1630              + 6 * (1-param)*param * (-self.y1 + self.y2) \
1631              + 3 * param*param * (-self.y2 + self.y3)
1632         xddot = 6 * (1-param) * (self.x0 - 2*self.x1 + self.x2) \
1633               + 6 * param * (self.x1 - 2*self.x2 + self.x3)
1634         yddot = 6 * (1-param) * (self.y0 - 2*self.y1 + self.y2) \
1635               + 6 * param * (self.y1 - 2*self.y2 + self.y3)
1636         return (xdot**2 + ydot**2)**1.5 / (xdot*yddot - ydot*xddot)
1637
1638     def end_pt(self):
1639         return self.x3, self.y3
1640
1641     def intersect(self, other, epsilon=1e-5):
1642         if isinstance(other, normline):
1643             return  _intersectnormcurves(self, 0, 1, other._normcurve(), 0, 1, epsilon)
1644         else:
1645             return  _intersectnormcurves(self, 0, 1, other, 0, 1, epsilon)
1646
1647     def isstraight(self, epsilon=1e-5):
1648         """check wheter the normcurve is approximately straight"""
1649
1650         # just check, whether the modulus of the difference between
1651         # the length of the control polygon
1652         # (i.e. |P1-P0|+|P2-P1|+|P3-P2|) and the length of the
1653         # straight line between starting and ending point of the
1654         # normcurve (i.e. |P3-P1|) is smaller the epsilon
1655         return abs(math.hypot(self.x1-self.x0, self.y1-self.y0)+
1656                    math.hypot(self.x2-self.x1, self.y2-self.y1)+
1657                    math.hypot(self.x3-self.x2, self.y3-self.y2)-
1658                    math.hypot(self.x3-self.x0, self.y3-self.y0))<epsilon
1659
1660     def midpointsplit(self):
1661         """splits bpathel at midpoint returning bpath with two bpathels"""
1662
1663         # for efficiency reason, we do not use self.split(0.5)!
1664
1665         # first, we have to calculate the  midpoints between adjacent
1666         # control points
1667         x01 = 0.5*(self.x0+self.x1)
1668         y01 = 0.5*(self.y0+self.y1)
1669         x12 = 0.5*(self.x1+self.x2)
1670         y12 = 0.5*(self.y1+self.y2)
1671         x23 = 0.5*(self.x2+self.x3)
1672         y23 = 0.5*(self.y2+self.y3)
1673
1674         # In the next iterative step, we need the midpoints between 01 and 12
1675         # and between 12 and 23
1676         x01_12 = 0.5*(x01+x12)
1677         y01_12 = 0.5*(y01+y12)
1678         x12_23 = 0.5*(x12+x23)
1679         y12_23 = 0.5*(y12+y23)
1680
1681         # Finally the midpoint is given by
1682         xmidpoint = 0.5*(x01_12+x12_23)
1683         ymidpoint = 0.5*(y01_12+y12_23)
1684
1685         return (normcurve(self.x0, self.y0,
1686                           x01, y01,
1687                           x01_12, y01_12,
1688                           xmidpoint, ymidpoint),
1689                 normcurve(xmidpoint, ymidpoint,
1690                           x12_23, y12_23,
1691                           x23, y23,
1692                           self.x3, self.y3))
1693
1694     def reverse(self):
1695         self.x0, self.y0, self.x1, self.y1, self.x2, self.y2, self.x3, self.y3 = \
1696         self.x3, self.y3, self.x2, self.y2, self.x1, self.y1, self.x0, self.y0
1697
1698     def reversed(self):
1699         return normcurve(self.x3, self.y3, self.x2, self.y2, self.x1, self.y1, self.x0, self.y0)
1700
1701     def seglengths(self, paraminterval, epsilon=1e-5):
1702         """returns the list of segment line lengths (in pts) of the bpathel
1703            together with the length of the parameterinterval"""
1704
1705         # lower and upper bounds for the arclen
1706         lowerlen = math.hypot(self.x3-self.x0, self.y3-self.y0)
1707         upperlen = ( math.hypot(self.x1-self.x0, self.y1-self.y0) +
1708                      math.hypot(self.x2-self.x1, self.y2-self.y1) +
1709                      math.hypot(self.x3-self.x2, self.y3-self.y2) )
1710
1711         # instead of isstraight method:
1712         if abs(upperlen-lowerlen)<epsilon:
1713             return [( 0.5*(upperlen+lowerlen), paraminterval )]
1714         else:
1715             (a, b) = self.midpointsplit()
1716             return a.seglengths(0.5*paraminterval, epsilon) + b.seglengths(0.5*paraminterval, epsilon)
1717
1718     def _split(self, parameters):
1719         """return list of normcurve corresponding to split at parameters"""
1720
1721         # first, we calculate the coefficients corresponding to our
1722         # original bezier curve. These represent a useful starting
1723         # point for the following change of the polynomial parameter
1724         a0x = self.x0
1725         a0y = self.y0
1726         a1x = 3*(-self.x0+self.x1)
1727         a1y = 3*(-self.y0+self.y1)
1728         a2x = 3*(self.x0-2*self.x1+self.x2)
1729         a2y = 3*(self.y0-2*self.y1+self.y2)
1730         a3x = -self.x0+3*(self.x1-self.x2)+self.x3
1731         a3y = -self.y0+3*(self.y1-self.y2)+self.y3
1732
1733         if parameters[0]!=0:
1734             parameters = [0] + parameters
1735         if parameters[-1]!=1:
1736             parameters = parameters + [1]
1737
1738         result = []
1739
1740         for i in range(len(parameters)-1):
1741             t1 = parameters[i]
1742             dt = parameters[i+1]-t1
1743
1744             # [t1,t2] part
1745             #
1746             # the new coefficients of the [t1,t1+dt] part of the bezier curve
1747             # are then given by expanding
1748             #  a0 + a1*(t1+dt*u) + a2*(t1+dt*u)**2 +
1749             #  a3*(t1+dt*u)**3 in u, yielding
1750             #
1751             #   a0 + a1*t1 + a2*t1**2 + a3*t1**3        +
1752             #   ( a1 + 2*a2 + 3*a3*t1**2 )*dt    * u    +
1753             #   ( a2 + 3*a3*t1 )*dt**2           * u**2 +
1754             #   a3*dt**3                         * u**3
1755             #
1756             # from this values we obtain the new control points by inversion
1757             #
1758             # XXX: we could do this more efficiently by reusing for
1759             # (x0, y0) the control point (x3, y3) from the previous
1760             # Bezier curve
1761
1762             x0 = a0x + a1x*t1 + a2x*t1*t1 + a3x*t1*t1*t1
1763             y0 = a0y + a1y*t1 + a2y*t1*t1 + a3y*t1*t1*t1
1764             x1 = (a1x+2*a2x*t1+3*a3x*t1*t1)*dt/3.0 + x0
1765             y1 = (a1y+2*a2y*t1+3*a3y*t1*t1)*dt/3.0 + y0
1766             x2 = (a2x+3*a3x*t1)*dt*dt/3.0 - x0 + 2*x1
1767             y2 = (a2y+3*a3y*t1)*dt*dt/3.0 - y0 + 2*y1
1768             x3 = a3x*dt*dt*dt + x0 - 3*x1 + 3*x2
1769             y3 = a3y*dt*dt*dt + y0 - 3*y1 + 3*y2
1770
1771             result.append(normcurve(x0, y0, x1, y1, x2, y2, x3, y3))
1772
1773         return result
1774
1775     def split(self, parameters):
1776         if parameters:
1777             # we need to split
1778             bps = self._split(list(parameters))
1779
1780             if parameters[0]==0:
1781                 result = [None]
1782             else:
1783                 bp0 = bps[0]
1784                 result = [normcurve(self.x0, self.y0, bp0.x1, bp0.y1, bp0.x2, bp0.y2, bp0.x3, bp0.y3)]
1785                 bps = bps[1:]
1786
1787             for bp in bps:
1788                 result.append(normcurve(bp.x0, bp.y0, bp.x1, bp.y1, bp.x2, bp.y2, bp.x3, bp.y3))
1789
1790             if parameters[-1]==1:
1791                 result.append(None)
1792         else:
1793             result = []
1794         return result
1795
1796     def tangentvector_pt(self, t):
1797         tvectx = (3*(  -self.x0+3*self.x1-3*self.x2+self.x3)*t*t +
1798                   2*( 3*self.x0-6*self.x1+3*self.x2        )*t +
1799                     (-3*self.x0+3*self.x1                  ))
1800         tvecty = (3*(  -self.y0+3*self.y1-3*self.y2+self.y3)*t*t +
1801                   2*( 3*self.y0-6*self.y1+3*self.y2        )*t +
1802                     (-3*self.y0+3*self.y1                  ))
1803         return (tvectx, tvecty)
1804
1805     def transform(self, trafo):
1806         self.x0, self.y0 = trafo._apply(self.x0, self.y0)
1807         self.x1, self.y1 = trafo._apply(self.x1, self.y1)
1808         self.x2, self.y2 = trafo._apply(self.x2, self.y2)
1809         self.x3, self.y3 = trafo._apply(self.x3, self.y3)
1810
1811     def transformed(self, trafo):
1812         return normcurve(*(trafo._apply(self.x0, self.y0)+
1813                            trafo._apply(self.x1, self.y1)+
1814                            trafo._apply(self.x2, self.y2)+
1815                            trafo._apply(self.x3, self.y3)))
1816
1817     def outputPS(self, file):
1818         file.write("%g %g %g %g %g %g curveto\n" % (self.x1, self.y1, self.x2, self.y2, self.x3, self.y3))
1819
1820     def outputPDF(self, file):
1821         file.write("%f %f %f %f %f %f c\n" % (self.x1, self.y1, self.x2, self.y2, self.x3, self.y3))
1822
1823 #
1824 # normpaths are made up of normsubpaths, which represent connected line segments
1825 #
1826
1827 class normsubpath:
1828
1829     """sub path of a normalized path
1830
1831     A subpath consists of a list of normpathels, i.e., lines and bcurves
1832     and can either be closed or not.
1833
1834     Some invariants, which have to be obeyed:
1835     - All normpathels have to be longer than epsilon pts.
1836     - The last point of a normpathel and the first point of the next
1837     element have to be equal.
1838     - When the path is closed, the last normpathel has to be a
1839     normline and the last point of this normline has to be equal
1840     to the first point of the first normpathel, except when
1841     this normline would be too short.
1842     """
1843
1844     __slots__ = "normpathels", "closed", "epsilon"
1845
1846     def __init__(self, normpathels, closed, epsilon=1e-5):
1847         self.normpathels = [npel for npel in normpathels if not npel.isstraight(epsilon) or npel.arclen_pt(epsilon)>epsilon]
1848         self.closed = closed
1849         self.epsilon = epsilon
1850
1851     def __str__(self):
1852         return "subpath(%s, [%s])" % (self.closed and "closed" or "open",
1853                                     ", ".join(map(str, self.normpathels)))
1854
1855     def arclen_pt(self):
1856         """returns total arc length of normsubpath in pts with accuracy epsilon"""
1857         return sum([npel.arclen_pt(self.epsilon) for npel in self.normpathels])
1858
1859     def _arclentoparam_pt(self, lengths):
1860         """returns [t, l] where t are parameter value(s) matching given length(s)
1861         and l is the total length of the normsubpath
1862         The parameters are with respect to the normsubpath: t in [0, self.range()]
1863         lengths that are < 0 give parameter 0"""
1864
1865         allarclen = 0
1866         allparams = [0] * len(lengths)
1867         rests = copy.copy(lengths)
1868
1869         for pel in self.normpathels:
1870             params, arclen = pel._arclentoparam_pt(rests, self.epsilon)
1871             allarclen += arclen
1872             for i in range(len(rests)):
1873                 if rests[i] >= 0:
1874                     rests[i] -= arclen
1875                     allparams[i] += params[i]
1876
1877         return (allparams, allarclen)
1878
1879     def at_pt(self, param):
1880         """return coordinates in pts of sub path at parameter value param
1881
1882         The parameter param must be smaller or equal to the number of
1883         segments in the normpath, otherwise None is returned.
1884         """
1885         try:
1886             return self.normpathels[int(param-self.epsilon)].at_pt(param-int(param-self.epsilon))
1887         except:
1888             raise PathException("parameter value param out of range")
1889
1890     def bbox(self):
1891         if self.normpathels:
1892             abbox = self.normpathels[0].bbox()
1893             for anormpathel in self.normpathels[1:]:
1894                 abbox += anormpathel.bbox()
1895             return abbox
1896         else:
1897             return None
1898
1899     def begin_pt(self):
1900         return self.normpathels[0].begin_pt()
1901
1902     def curvradius_pt(self, param):
1903         try:
1904             return self.normpathels[int(param-self.epsilon)].curvradius_pt(param-int(param-self.epsilon))
1905         except:
1906             raise PathException("parameter value param out of range")
1907
1908     def end_pt(self):
1909         return self.normpathels[-1].end_pt()
1910
1911     def intersect(self, other):
1912         """intersect self with other normsubpath
1913
1914         returns a tuple of lists consisting of the parameter values
1915         of the intersection points of the corresponding normsubpath
1916
1917         """
1918         intersections = ([], [])
1919         epsilon = min(self.epsilon, other.epsilon)
1920         # Intersect all subpaths of self with the subpaths of other
1921         for t_a, pel_a  in enumerate(self.normpathels):
1922             for t_b, pel_b in enumerate(other.normpathels):
1923                 for intersection in pel_a.intersect(pel_b, epsilon):
1924                     # check whether an intersection occurs at the end
1925                     # of a closed subpath. If yes, we don't include it
1926                     # in the list of intersections to prevent a
1927                     # duplication of intersection points
1928                     if not ((self.closed and self.range()-intersection[0]-t_a<epsilon) or
1929                             (other.closed and other.range()-intersection[1]-t_b<epsilon)):
1930                         intersections[0].append(intersection[0]+t_a)
1931                         intersections[1].append(intersection[1]+t_b)
1932         return intersections
1933
1934     def range(self):
1935         """return maximal parameter value, i.e. number of line/curve segments"""
1936         return len(self.normpathels)
1937
1938     def reverse(self):
1939         self.normpathels.reverse()
1940         for npel in self.normpathels:
1941             npel.reverse()
1942
1943     def reversed(self):
1944         nnormpathels = []
1945         for i in range(len(self.normpathels)):
1946             nnormpathels.append(self.normpathels[-(i+1)].reversed())
1947         return normsubpath(nnormpathels, self.closed)
1948
1949     def split(self, params):
1950         """split normsubpath at list of parameter values params and return list
1951         of normsubpaths
1952
1953         The parameter list params has to be sorted. Note that each element of
1954         the resulting list is an open normsubpath.
1955         """
1956
1957         if min(params) < -self.epsilon or max(params) > self.range()+self.epsilon:
1958             raise PathException("parameter for split of subpath out of range")
1959
1960         result = []
1961         npels = None
1962         for t, pel in enumerate(self.normpathels):
1963             # determine list of splitting parameters relevant for pel
1964             nparams = []
1965             for nt in params:
1966                 if t+1 >= nt:
1967                     nparams.append(nt-t)
1968                     params = params[1:]
1969
1970             # now we split the path at the filtered parameter values
1971             # This yields a list of normpathels and possibly empty
1972             # segments marked by None
1973             splitresult = pel.split(nparams)
1974             if splitresult:
1975                 # first split?
1976                 if npels is None:
1977                     if splitresult[0] is None:
1978                         # mark split at the beginning of the normsubpath
1979                         result = [None]
1980                     else:
1981                         result.append(normsubpath([splitresult[0]], 0))
1982                 else:
1983                     npels.append(splitresult[0])
1984                     result.append(normsubpath(npels, 0))
1985                 for npel in splitresult[1:-1]:
1986                     result.append(normsubpath([npel], 0))
1987                 if len(splitresult)>1 and splitresult[-1] is not None:
1988                     npels = [splitresult[-1]]
1989                 else:
1990                     npels = []
1991             else:
1992                 if npels is None:
1993                     npels = [pel]
1994                 else:
1995                     npels.append(pel)
1996
1997         if npels:
1998             result.append(normsubpath(npels, 0))
1999         else:
2000             # mark split at the end of the normsubpath
2001             result.append(None)
2002
2003         # glue last and first segment together if the normsubpath was originally closed
2004         if self.closed:
2005             if result[0] is None:
2006                 result = result[1:]
2007             elif result[-1] is None:
2008                 result = result[:-1]
2009             else:
2010                 result[-1].normpathels.extend(result[0].normpathels)
2011                 result = result[1:]
2012         return result
2013
2014     def tangent(self, param, length=None):
2015         tx, ty = self.at_pt(param)
2016         try:
2017             tdx, tdy = self.normpathels[int(param-self.epsilon)].tangentvector_pt(param-int(param-self.epsilon))
2018         except:
2019             raise PathException("parameter value param out of range")
2020         tlen = math.hypot(tdx, tdy)
2021         if not (length is None or tlen==0):
2022             sfactor = unit.topt(length)/tlen
2023             tdx *= sfactor
2024             tdy *= sfactor
2025         return line_pt(tx, ty, tx+tdx, ty+tdy)
2026
2027     def trafo(self, param):
2028         tx, ty = self.at_pt(param)
2029         try:
2030             tdx, tdy = self.normpathels[int(param-self.epsilon)].tangentvector_pt(param-int(param-self.epsilon))
2031         except:
2032             raise PathException("parameter value param out of range")
2033         return trafo.translate_pt(tx, ty)*trafo.rotate(degrees(math.atan2(tdy, tdx)))
2034
2035     def transform(self, trafo):
2036         """transform sub path according to trafo"""
2037         for pel in self.normpathels:
2038             pel.transform(trafo)
2039
2040     def transformed(self, trafo):
2041         """return sub path transformed according to trafo"""
2042         nnormpathels = []
2043         for pel in self.normpathels:
2044             nnormpathels.append(pel.transformed(trafo))
2045         return normsubpath(nnormpathels, self.closed)
2046
2047     def outputPS(self, file):
2048         # if the normsubpath is closed, we must not output a normline at
2049         # the end
2050         if not self.normpathels:
2051             return
2052         if self.closed and isinstance(self.normpathels[-1], normline):
2053             normpathels = self.normpathels[:-1]
2054         else:
2055             normpathels = self.normpathels
2056         if normpathels:
2057             file.write("%g %g moveto\n" % self.begin_pt())
2058             for anormpathel in normpathels:
2059                 anormpathel.outputPS(file)
2060         if self.closed:
2061             file.write("closepath\n")
2062
2063     def outputPDF(self, file):
2064         # if the normsubpath is closed, we must not output a normline at
2065         # the end
2066         if not self.normpathels:
2067             return
2068         if self.closed and isinstance(self.normpathels[-1], normline):
2069             normpathels = self.normpathels[:-1]
2070         else:
2071             normpathels = self.normpathels
2072         if normpathels:
2073             file.write("%f %f m\n" % self.begin_pt())
2074             for anormpathel in normpathels:
2075                 anormpathel.outputPDF(file)
2076         if self.closed:
2077             file.write("h\n")
2078
2079 #
2080 # the normpath class
2081 #
2082
2083 class normpath(path):
2084
2085     """normalized path
2086
2087     A normalized path consists of a list of normalized sub paths.
2088
2089     """
2090
2091     def __init__(self, arg=[], epsilon=1e-5):
2092         """ construct a normpath from another normpath passed as arg,
2093         a path or a list of normsubpaths. An accuracy of epsilon pts
2094         is used for numerical calculations.
2095         """
2096
2097         self.epsilon = epsilon
2098         if isinstance(arg, normpath):
2099             self.subpaths = copy.copy(arg.subpaths)
2100             return
2101         elif isinstance(arg, path):
2102             # split path in sub paths
2103             self.subpaths = []
2104             currentsubpathels = []
2105             context = _pathcontext()
2106             for pel in arg.path:
2107                 for npel in pel._normalized(context):
2108                     if isinstance(npel, moveto_pt):
2109                         if currentsubpathels:
2110                             # append open sub path
2111                             self.subpaths.append(normsubpath(currentsubpathels, 0, epsilon))
2112                         # start new sub path
2113                         currentsubpathels = []
2114                     elif isinstance(npel, closepath):
2115                         if currentsubpathels:
2116                             # append closed sub path
2117                             currentsubpathels.append(normline(context.currentpoint[0], context.currentpoint[1],
2118                                                               context.currentsubpath[0], context.currentsubpath[1]))
2119                         self.subpaths.append(normsubpath(currentsubpathels, 1, epsilon))
2120                         currentsubpathels = []
2121                     else:
2122                         currentsubpathels.append(npel)
2123                 pel._updatecontext(context)
2124
2125             if currentsubpathels:
2126                 # append open sub path
2127                 self.subpaths.append(normsubpath(currentsubpathels, 0, epsilon))
2128         else:
2129             # we expect a list of normsubpaths
2130             self.subpaths = list(arg)
2131
2132     def __add__(self, other):
2133         result = normpath(other)
2134         result.subpaths = self.subpaths + result.subpaths
2135         return result
2136
2137     def __iadd__(self, other):
2138         self.subpaths += normpath(other).subpaths
2139         return self
2140
2141     def __nonzero__(self):
2142         return len(self.subpaths)>0
2143
2144     def __str__(self):
2145         return "normpath(%s)" % ", ".join(map(str, self.subpaths))
2146
2147     def _findsubpath(self, param, arclen):
2148         """return a tuple (subpath, rparam), where subpath is the subpath
2149         containing the position specified by either param or arclen and rparam
2150         is the corresponding parameter value in this subpath.
2151         """
2152
2153         if param is not None and arclen is not None:
2154             raise PathException("either param or arclen has to be specified, but not both")
2155         elif arclen is not None:
2156             param = self.arclentoparam(arclen)
2157
2158         spt = 0
2159         for sp in self.subpaths:
2160             sprange = sp.range()
2161             if spt <= param <= sprange+spt+self.epsilon:
2162                 return sp, param-spt
2163             spt += sprange
2164         raise PathException("parameter value out of range")
2165
2166     def append(self, pathel):
2167         # XXX factor parts of this code out
2168         if self.subpaths[-1].closed:
2169             context = _pathcontext(self.end_pt(), None)
2170             currentsubpathels = []
2171         else:
2172             context = _pathcontext(self.end_pt(), self.subpaths[-1].begin_pt())
2173             currentsubpathels = self.subpaths[-1].normpathels
2174             self.subpaths = self.subpaths[:-1]
2175         for npel in pathel._normalized(context):
2176             if isinstance(npel, moveto_pt):
2177                 if currentsubpathels:
2178                     # append open sub path
2179                     self.subpaths.append(normsubpath(currentsubpathels, 0, self.epsilon))
2180                 # start new sub path
2181                 currentsubpathels = []
2182             elif isinstance(npel, closepath):
2183                 if currentsubpathels:
2184                     # append closed sub path
2185                     currentsubpathels.append(normline(context.currentpoint[0], context.currentpoint[1],
2186                                                       context.currentsubpath[0], context.currentsubpath[1]))
2187                     self.subpaths.append(normsubpath(currentsubpathels, 1, self.epsilon))
2188                 currentsubpathels = []
2189             else:
2190                 currentsubpathels.append(npel)
2191
2192         if currentsubpathels:
2193             # append open sub path
2194             self.subpaths.append(normsubpath(currentsubpathels, 0, self.epsilon))
2195
2196     def arclen_pt(self):
2197         """returns total arc length of normpath in pts"""
2198         return sum([sp.arclen_pt() for sp in self.subpaths])
2199
2200     def arclen(self):
2201         """returns total arc length of normpath"""
2202         return unit.t_pt(self.arclen_pt())
2203
2204     def arclentoparam(self, lengths):
2205         """returns the parameter value(s) matching the given length(s)
2206
2207         all given lengths must be positive.
2208         A length greater than the total arclength will give self.range()"""
2209
2210         rests = [unit.topt(length) for length in helper.ensuresequence(lengths)]
2211         allparams = [0] * len(helper.ensuresequence(lengths))
2212
2213         for sp in self.subpaths:
2214             # we need arclen for knowing when all the parameters are done
2215             # for lengths that are done: rests[i] is negative
2216             # sp._arclentoparam has to ignore such lengths
2217             params, arclen = sp._arclentoparam_pt(rests)
2218             finis = 0 # number of lengths that are done
2219             for i in range(len(rests)):
2220                 if rests[i] >= 0:
2221                   rests[i] -= arclen
2222                   allparams[i] += params[i]
2223                 else:
2224                     finis += 1
2225             if finis == len(rests): break
2226
2227         if not helper.issequence(lengths): allparams = allparams[0]
2228         return allparams
2229
2230     def at_pt(self, param=None, arclen=None):
2231         """return coordinates in pts of path at either parameter value param
2232         or arc length arclen.
2233
2234         At discontinuities in the path, the limit from below is returned.
2235         """
2236         sp, param = self._findsubpath(param, arclen)
2237         return sp.at_pt(param)
2238
2239     def at(self, param=None, arclen=None):
2240         """return coordinates of path at either parameter value param
2241         or arc length arclen.
2242
2243         At discontinuities in the path, the limit from below is returned
2244         """
2245         x, y = self.at_pt(param, arclen)
2246         return unit.t_pt(x), unit.t_pt(y)
2247
2248     def bbox(self):
2249         abbox = None
2250         for sp in self.subpaths:
2251             nbbox =  sp.bbox()
2252             if abbox is None:
2253                 abbox = nbbox
2254             elif nbbox:
2255                 abbox += nbbox
2256         return abbox
2257
2258     def begin_pt(self):
2259         """return coordinates of first point of first subpath in path (in pts)"""
2260         if self.subpaths:
2261             return self.subpaths[0].begin_pt()
2262         else:
2263             raise PathException("cannot return first point of empty path")
2264
2265     def begin(self):
2266         """return coordinates of first point of first subpath in path"""
2267         x, y = self.begin_pt()
2268         return unit.t_pt(x), unit.t_pt(y)
2269
2270     def curvradius_pt(self, param=None, arclen=None):
2271         """Returns the curvature radius in pts (or None if infinite)
2272         at parameter param or arc length arclen.  This is the inverse
2273         of the curvature at this parameter
2274
2275         Please note that this radius can be negative or positive,
2276         depending on the sign of the curvature"""
2277         sp, param = self._findsubpath(param, arclen)
2278         return sp.curvradius_pt(param)
2279
2280     def curvradius(self, param=None, arclen=None):
2281         """Returns the curvature radius (or None if infinite) at
2282         parameter param or arc length arclen.  This is the inverse of
2283         the curvature at this parameter
2284
2285         Please note that this radius can be negative or positive,
2286         depending on the sign of the curvature"""
2287         radius = self.curvradius_pt(param, arclen)
2288         if radius is not None:
2289             radius = unit.t_pt(radius)
2290         return radius
2291
2292     def end_pt(self):
2293         """return coordinates of last point of last subpath in path (in pts)"""
2294         if self.subpaths:
2295             return self.subpaths[-1].end_pt()
2296         else:
2297             raise PathException("cannot return last point of empty path")
2298
2299     def end(self):
2300         """return coordinates of last point of last subpath in path"""
2301         x, y = self.end_pt()
2302         return unit.t_pt(x), unit.t_pt(y)
2303
2304     def glue(self, other):
2305         if not self.subpaths:
2306             raise PathException("cannot glue to end of empty path")
2307         if self.subpaths[-1].closed:
2308             raise PathException("cannot glue to end of closed sub path")
2309         other = normpath(other)
2310         if not other.subpaths:
2311             raise PathException("cannot glue empty path")
2312
2313         self.subpaths[-1].normpathels += other.subpaths[0].normpathels
2314         self.subpaths += other.subpaths[1:]
2315         return self
2316
2317     def intersect(self, other):
2318         """intersect self with other path
2319
2320         returns a tuple of lists consisting of the parameter values
2321         of the intersection points of the corresponding normpath
2322
2323         """
2324         if not isinstance(other, normpath):
2325             other = normpath(other)
2326
2327         # here we build up the result
2328         intersections = ([], [])
2329
2330         # Intersect all subpaths of self with the subpaths of
2331         # other. Here, st_a, st_b are the parameter values
2332         # corresponding to the first point of the subpaths sp_a and
2333         # sp_b, respectively.
2334         st_a = 0
2335         for sp_a in self.subpaths:
2336             st_b =0
2337             for sp_b in other.subpaths:
2338                 for intersection in zip(*sp_a.intersect(sp_b)):
2339                     intersections[0].append(intersection[0]+st_a)
2340                     intersections[1].append(intersection[1]+st_b)
2341                 st_b += sp_b.range()
2342             st_a += sp_a.range()
2343         return intersections
2344
2345     def range(self):
2346         """return maximal value for parameter value param"""
2347         return sum([sp.range() for sp in self.subpaths])
2348
2349     def reverse(self):
2350         """reverse path"""
2351         self.subpaths.reverse()
2352         for sp in self.subpaths:
2353             sp.reverse()
2354
2355     def reversed(self):
2356         """return reversed path"""
2357         nnormpath = normpath()
2358         for i in range(len(self.subpaths)):
2359             nnormpath.subpaths.append(self.subpaths[-(i+1)].reversed())
2360         return nnormpath
2361
2362     def split(self, params):
2363         """split path at parameter values params
2364
2365         Note that the parameter list has to be sorted.
2366
2367         """
2368
2369         # check whether parameter list is really sorted
2370         sortedparams = list(params)
2371         sortedparams.sort()
2372         if sortedparams!=list(params):
2373             raise ValueError("split parameter list params has to be sorted")
2374
2375         # we construct this list of normpaths
2376         result = []
2377
2378         # the currently built up normpath
2379         np = normpath()
2380
2381         t0 = 0
2382         for subpath in self.subpaths:
2383             tf = t0+subpath.range()
2384             if params and tf>=params[0]:
2385                 # split this subpath
2386                 # determine the relevant splitting params
2387                 for i in range(len(params)):
2388                     if params[i]>tf: break
2389                 else:
2390                     i = len(params)
2391
2392                 splitsubpaths = subpath.split([x-t0 for x in params[:i]])
2393                 # handle first element, which may be None, separately
2394                 if splitsubpaths[0] is None:
2395                     if not np.subpaths:
2396                         result.append(None)
2397                     else:
2398                         result.append(np)
2399                         np = normpath()
2400                     splitsubpaths.pop(0)
2401
2402                 for sp in splitsubpaths[:-1]:
2403                     np.subpaths.append(sp)
2404                     result.append(np)
2405                     np = normpath()
2406
2407                 # handle last element which may be None, separately
2408                 if splitsubpaths:
2409                     if splitsubpaths[-1] is None:
2410                         if np.subpaths:
2411                             result.append(np)
2412                             np = normpath()
2413                     else:
2414                         np.subpaths.append(splitsubpaths[-1])
2415
2416                 params = params[i:]
2417             else:
2418                 # append whole subpath to current normpath
2419                 np.subpaths.append(subpath)
2420             t0 = tf
2421
2422         if np.subpaths:
2423             result.append(np)
2424         else:
2425             # mark split at the end of the normsubpath
2426             result.append(None)
2427
2428         return result
2429
2430     def tangent(self, param=None, arclen=None, length=None):
2431         """return tangent vector of path at either parameter value param
2432         or arc length arclen.
2433
2434         At discontinuities in the path, the limit from below is returned.
2435         If length is not None, the tangent vector will be scaled to
2436         the desired length.
2437         """
2438         sp, param = self._findsubpath(param, arclen)
2439         return sp.tangent(param, length)
2440
2441     def transform(self, trafo):
2442         """transform path according to trafo"""
2443         for sp in self.subpaths:
2444             sp.transform(trafo)
2445
2446     def transformed(self, trafo):
2447         """return path transformed according to trafo"""
2448         return normpath([sp.transformed(trafo) for sp in self.subpaths])
2449
2450     def trafo(self, param=None, arclen=None):
2451         """return transformation at either parameter value param or arc length arclen"""
2452         sp, param = self._findsubpath(param, arclen)
2453         return sp.trafo(param)
2454
2455     def outputPS(self, file):
2456         for sp in self.subpaths:
2457             sp.outputPS(file)
2458
2459     def outputPDF(self, file):
2460         for sp in self.subpaths:
2461             sp.outputPDF(file)
2462