@@ -1561,14 +1561,12 @@ def draw(self, renderer):
1561
1561
calculation much easier than doing rotated ellipse
1562
1562
intersection directly).
1563
1563
1564
- This uses the "line intersecting a circle" algorithm
1565
- from:
1564
+ This uses the "line intersecting a circle" algorithm from:
1566
1565
1567
1566
Vince, John. *Geometry for Computer Graphics: Formulae,
1568
1567
Examples & Proofs.* London: Springer-Verlag, 2005.
1569
1568
1570
- 2. The angles of each of the intersection points are
1571
- calculated.
1569
+ 2. The angles of each of the intersection points are calculated.
1572
1570
1573
1571
3. Proceeding counterclockwise starting in the positive
1574
1572
x-direction, each of the visible arc-segments between the
@@ -1601,7 +1599,7 @@ def theta_stretch(theta, scale):
1601
1599
self ._path = Path .arc (theta1 , theta2 )
1602
1600
return Patch .draw (self , renderer )
1603
1601
1604
- def iter_circle_intersect_on_line (x0 , y0 , x1 , y1 ):
1602
+ def line_circle_intersect (x0 , y0 , x1 , y1 ):
1605
1603
dx = x1 - x0
1606
1604
dy = y1 - y0
1607
1605
dr2 = dx * dx + dy * dy
@@ -1613,7 +1611,7 @@ def iter_circle_intersect_on_line(x0, y0, x1, y1):
1613
1611
if discrim == 0.0 :
1614
1612
x = (D * dy ) / dr2
1615
1613
y = (- D * dx ) / dr2
1616
- yield x , y
1614
+ return np . array ([[ x , y ]])
1617
1615
elif discrim > 0.0 :
1618
1616
# The definition of "sign" here is different from
1619
1617
# np.sign: we never want to get 0.0
@@ -1622,12 +1620,15 @@ def iter_circle_intersect_on_line(x0, y0, x1, y1):
1622
1620
else :
1623
1621
sign_dy = 1.0
1624
1622
sqrt_discrim = np .sqrt (discrim )
1625
- for sign in (1. , - 1. ):
1626
- x = (D * dy + sign * sign_dy * dx * sqrt_discrim ) / dr2
1627
- y = (- D * dx + sign * np .abs (dy ) * sqrt_discrim ) / dr2
1628
- yield x , y
1623
+ return np .array (
1624
+ [[(D * dy + sign_dy * dx * sqrt_discrim ) / dr2 ,
1625
+ (- D * dx + abs (dy ) * sqrt_discrim ) / dr2 ],
1626
+ [(D * dy - sign_dy * dx * sqrt_discrim ) / dr2 ,
1627
+ (- D * dx - abs (dy ) * sqrt_discrim ) / dr2 ]])
1628
+ else :
1629
+ return np .empty ((0 , 2 ))
1629
1630
1630
- def iter_circle_intersect_on_line_seg (x0 , y0 , x1 , y1 ):
1631
+ def segment_circle_intersect (x0 , y0 , x1 , y1 ):
1631
1632
epsilon = 1e-9
1632
1633
if x1 < x0 :
1633
1634
x0e , x1e = x1 , x0
@@ -1637,17 +1638,13 @@ def iter_circle_intersect_on_line_seg(x0, y0, x1, y1):
1637
1638
y0e , y1e = y1 , y0
1638
1639
else :
1639
1640
y0e , y1e = y0 , y1
1640
- x0e -= epsilon
1641
- y0e -= epsilon
1642
- x1e += epsilon
1643
- y1e += epsilon
1644
- for x , y in iter_circle_intersect_on_line (x0 , y0 , x1 , y1 ):
1645
- if x0e <= x <= x1e and y0e <= y <= y1e :
1646
- yield x , y
1641
+ xys = line_circle_intersect (x0 , y0 , x1 , y1 )
1642
+ xs , ys = xys .T
1643
+ return xys [(x0e - epsilon < xs ) & (xs < x1e + epsilon )
1644
+ & (y0e - epsilon < ys ) & (ys < y1e + epsilon )]
1647
1645
1648
1646
# Transforms the axes box_path so that it is relative to the unit
1649
- # circle in the same way that it is relative to the desired
1650
- # ellipse.
1647
+ # circle in the same way that it is relative to the desired ellipse.
1651
1648
box_path = Path .unit_rectangle ()
1652
1649
box_path_transform = transforms .BboxTransformTo (self .axes .bbox ) + \
1653
1650
self .get_transform ().inverted ()
@@ -1656,16 +1653,10 @@ def iter_circle_intersect_on_line_seg(x0, y0, x1, y1):
1656
1653
thetas = set ()
1657
1654
# For each of the point pairs, there is a line segment
1658
1655
for p0 , p1 in zip (box_path .vertices [:- 1 ], box_path .vertices [1 :]):
1659
- x0 , y0 = p0
1660
- x1 , y1 = p1
1661
- for x , y in iter_circle_intersect_on_line_seg (x0 , y0 , x1 , y1 ):
1662
- theta = np .arccos (x )
1663
- if y < 0 :
1664
- theta = 2 * np .pi - theta
1665
- # Convert radians to angles
1666
- theta = np .rad2deg (theta )
1667
- if theta1 < theta < theta2 :
1668
- thetas .add (theta )
1656
+ xy = segment_circle_intersect (* p0 , * p1 )
1657
+ x , y = xy .T
1658
+ theta = np .rad2deg (np .arctan2 (y , x ))
1659
+ thetas .update (theta [(theta1 < theta ) & (theta < theta2 )])
1669
1660
thetas = sorted (thetas ) + [theta2 ]
1670
1661
1671
1662
last_theta = theta1
0 commit comments