Updated demo and removed call to 'mathutils.geometry'.

.gitignore

@@ -0,0 +1,2 @@
+.DS_Store
+**/__pycache__

bpypolyskel/bpyeuclid.py

@@ -1,109 +1,146 @@
 import mathutils
 
-# -------------------------------------------------------------------------
-# from https://bryceboe.com/2006/10/23/line-segment-intersection-algorithm/,
-# works fine with mathutils.Vector
-def ccw(A,B,C):
-    return (C.y-A.y) * (B.x-A.x) > (B.y-A.y) * (C.x-A.x)
-# Return true if line segments AB and CD intersect
-def intersect(A,B,C,D):
-    return ccw(A,C,D) != ccw(B,C,D) and ccw(A,B,C) != ccw(A,B,D)
-# -------------------------------------------------------------------------
+# From https://bryceboe.com/2006/10/23/line-segment-intersection-algorithm/.
+def ccw(A, B, C):
+    return (C.y - A.y) * (B.x - A.x) > (B.y - A.y) * (C.x - A.x)
+
+# Return true if line segments AB and CD intersect.
+def intersect(A, B, C, D):
+    return ccw(A, C, D) != ccw(B, C, D) and ccw(A, B, C) != ccw(A, B, D)
 
 def _intersect_line2_line2(A, B):
     d = B.v.y * A.v.x - B.v.x * A.v.y
+
     if d == 0:
         return None
 
     dy = A.p.y - B.p.y
     dx = A.p.x - B.p.x
     ua = (B.v.x * dy - B.v.y * dx) / d
+
     if not A.intsecttest(ua):
         return None
+
     ub = (A.v.x * dy - A.v.y * dx) / d
+
     if not B.intsecttest(ub):
         return None
 
     return mathutils.Vector((A.p.x + ua * A.v.x, A.p.y + ua * A.v.y))
 
+def _intersect_point_line(C, A, B):
+    # Drop a perpendicular from C to line AB and return
+    # the coordinates of the intersection point.
+    x1, y1 = A.x, A.y
+    x2, y2 = B.x, B.y
+    x3, y3 = C.x, C.y
+
+    if x1 == x2:
+        x, y = x1, y3
+    elif y1 == y2:
+        x, y = x3, y1
+    else:
+        # Calculate the equation of the line AB (y = mx + n).
+        m_AB = (y2 - y1) / (x2 - x1)
+        n_AB = y1 - m_AB * x1
+
+        # Calculate the equation of the perpendicular line through point C.
+        m_perp = -1 / m_AB
+        n_perp = y3 - m_perp * x3
+
+        # Solve the system of equations to find the intersection point (x, y).
+        x = (n_perp - n_AB) / (m_AB - m_perp)
+        y = m_AB * x + n_AB
+
+    return mathutils.Vector((x, y))
+
+
 class Edge2:
-    def __init__(self, p1, p2, norm=None, verts=None, center=mathutils.Vector((0,0))):  # evtl. conversion from 3D to 2D
+    def __init__(self, p1, p2, norm=None, verts=None, center=mathutils.Vector((0, 0))):
+        # Eventual conversion from 3D to 2D.
         """
         Args:
             p1 (mathutils.Vector | int): Index of the first edge vertex if <verts> is given or
-                the vector of the first edge vertex otherwise
+                the vector of the first edge vertex otherwise.
             p2 (mathutils.Vector | int): Index of the second edge vertex if <verts> is given or
-                the vector of the seconf edge vertex otherwise
+                the vector of the seconf edge vertex otherwise.
             norm (mathutils.Vector): Normalized edge vector
-            verts (list): Python list of vertices
+            verts (list): Python list of vertices.
         """
         if verts:
             self.i1 = p1
             self.i2 = p2
-            p1 = verts[p1]-center
-            p2 = verts[p2]-center
+            p1 = verts[p1] - center
+            p2 = verts[p2] - center
+
         self.p1 = mathutils.Vector((p1[0], p1[1]))
         self.p2 = mathutils.Vector((p2[0], p2[1]))
+
         if norm:
             self.norm = mathutils.Vector((norm[0], norm[1]))
         else:
-            norm = self.p2-self.p1
+            norm = self.p2 - self.p1
             norm.normalize()
             self.norm = norm
 
     def length_squared(self):
-        return (self.p2-self.p1).length_squared
+        return (self.p2 - self.p1).length_squared
 
 class Ray2:
     def __init__(self, _p, _v):
         self.p = _p
         self.p1 = _p
-        self.p2 = _p+_v
+        self.p2 = _p + _v
         self.v = _v
 
-    def intsecttest(self,u):
-        return u>=0.0
+    def intsecttest(self, u):
+        return u >= 0.0
 
-    def intersect(self,other):
-        return _intersect_line2_line2(self,other)
+    def intersect(self, other):
+        return _intersect_line2_line2(self, other)
 
 class Line2:
     def __init__(self, p1, p2=None, ptype=None):
-        if p2 is None: # p1 is a LineSegment2, a Line2 or a Ray2
+        # Note that 'p1' is a Line2 or Ray2 object.
+        if p2 is None:
             self.p = p1.p1.copy()
-            self.v = (p1.p2-p1.p1).copy()
+            self.v = (p1.p2 - p1.p1).copy()
         elif ptype == 'pp':
             self.p = p1.p.copy()
-            self.v = p2-p1
+            self.v = p2 - p1
         elif ptype == 'pv':
             self.p = p1.copy()
             self.v = p2.copy()
+
         self.p1 = self.p
-        self.p2 = self.p+self.v
+        self.p2 = self.p + self.v
 
-    def intsecttest(self,u):
+    def intsecttest(self, u):
         return True
 
-    def intersect(self,other):
-        intsect = _intersect_line2_line2(self,other)
-        return intsect
+    def intersect(self, other):
+        return _intersect_line2_line2(self, other)
 
-    def distance(self,other): # other is a vector
-        nearest = mathutils.geometry.intersect_point_line(other, self.p, self.p+self.v)[0]
-        dist = (other-nearest).length
-        return dist
+    def distance(self, other):
+        # Note that 'other' is a vector.
+        nearest = _intersect_point_line(other, self.p1, self.p2)
+        return (other - nearest).length
 
 def fitCircle3Points(points):
+    # Circle through three points using complex math, see answer in
+    # https://stackoverflow.com/questions/28910718/give-3-points-and-a-plot-circle.
     N = len(points)
-    # circle through three points usingg complex math, see answer in 
-    # https://stackoverflow.com/questions/28910718/give-3-points-and-a-plot-circle
+
     x = complex(points[0].x, points[0].y)
-    y = complex(points[N//2].x, points[N//2].y)
+    y = complex(points[N // 2].x, points[N // 2].y)
     z = complex(points[-1].x, points[-1].y)
-    w = z-x
-    w /= y-x
-    c = (x-y)*(w-abs(w)**2)/2j/w.imag-x
+
+    w = z - x
+    w /= y - x
+    c = (x - y) * (w - abs(w) ** 2) / 2j / w.imag - x
+
     x0 = -c.real
     y0 = -c.imag
-    R = abs(c+x)
-    return mathutils.Vector((x0,y0)), R
+    R = abs(c + x)
+
+    return mathutils.Vector((x0, y0)), R

demo.py

@@ -9,38 +9,78 @@
 import math
 import mathutils
 import matplotlib.pyplot as plt
+import sys
+import time
 
 from bpypolyskel import bpypolyskel
 from matplotlib.widgets import Slider
 
 
-# Define vertices of a polygon and a hole.
-_COORDS = [
-    # Polygon contour in counterclockwise order, seen from top.
-    # The polygon is on the left of this contour.
-    [0, 0, 0],
-    [10, 0, 0],
-    [10, 5, 0],
-    [45, 5, 0],
-    [45, 20, 0],
-    [10, 20, 0],
-    [10, 25, 0],
-    [0, 25, 0],
-    # Hole contour in clockwise order, seen from top.
-    # The polygon is on the left of this contour.
-    [5, 16, 0],
-    [35, 16, 0],
-    [35, 9, 0],
-    [5, 9, 0],
+# Define the outer loop of the floor plan (CCW order).
+_OUTER_LOOP = [
+    [1482.83, 495.548, 0],
+    [1623.76, 593.700, 0],
+    [1734.48, 724.983, 0],
+    [1807.45, 880.449, 0],
+    [1837.70, 1049.50, 0],
+    [1823.16, 1220.63, 0],
+    [1764.82, 1382.16, 0],
+    [1666.67, 1523.08, 0],
+    [1535.39, 1633.81, 0],
+    [1379.92, 1706.78, 0],
+    [1210.87, 1737.02, 0],
+    [1039.75, 1722.49, 0],
+    [878.216, 1664.15, 0],
+    [737.288, 1566.00, 0],
+    [626.566, 1434.72, 0],
+    [553.596, 1279.25, 0],
+    [523.349, 1110.20, 0],
+    [537.887, 939.072, 0],
+    [596.221, 777.543, 0],
+    [694.373, 636.616, 0],
+    [825.656, 525.894, 0],
+    [981.122, 452.923, 0],
+    [1150.18, 422.676, 0],
+    [1321.30, 437.215, 0],
 ]
 
+# Define the inner loop of the floor plan (CW order).
+_INNER_LOOP = [
+    [1351.82, 1444.80, 0],
+    [1440.44, 1388.03, 0],
+    [1511.34, 1310.26, 0],
+    [1559.70, 1216.78, 0],
+    [1582.23, 1113.98, 0],
+    [1577.37, 1008.85, 0],
+    [1545.47, 908.555, 0],
+    [1488.70, 819.936, 0],
+    [1410.93, 749.030, 0],
+    [1317.46, 700.668, 0],
+    [1214.65, 678.148, 0],
+    [1109.52, 683.002, 0],
+    [1009.23, 714.901, 0],
+    [920.609, 771.671, 0],
+    [849.702, 849.442, 0],
+    [801.341, 942.916, 0],
+    [778.820, 1045.72, 0],
+    [783.675, 1150.85, 0],
+    [815.574, 1251.14, 0],
+    [872.344, 1339.76, 0],
+    [950.115, 1410.67, 0],
+    [1043.59, 1459.03, 0],
+    [1146.39, 1481.55, 0],
+    [1251.53, 1476.70, 0],
+]
+
+_COORDS = _OUTER_LOOP + _INNER_LOOP
+
 # Convert the coordinates to 'mathutils.Vector' objects.
 VERTS = [mathutils.Vector(coords) for coords in _COORDS]
 
 # Define indices and lengths of polygon and hole.
 FIRST_VERTEX_INDEX = 0
-NUM_VERTS = 8
-HOLES_INFO = [(8, 4)]
+NUM_VERTS = len(_OUTER_LOOP)
+HOLES_INFO = [(NUM_VERTS, len(_INNER_LOOP))]
 
 # We let polygonize() compute the unit vectors and have no faces yet.
 FACES, UNIT_VECTORS = None, None
@@ -51,10 +91,8 @@
 DEFAULT_PITCH = 30.0
 
 
-def create_roof(pitch, ax):
-    ax.clear()
-
-    faces = bpypolyskel.polygonize(
+def get_roof_faces(pitch):
+    return bpypolyskel.polygonize(
         VERTS,
         FIRST_VERTEX_INDEX,
         NUM_VERTS,
@@ -65,14 +103,20 @@ def create_roof(pitch, ax):
         UNIT_VECTORS,
     )
 
+def create_roof(pitch, ax):
+    ax.clear()
+
+    # Create the roof faces.
+    faces = get_roof_faces(pitch)
+
     # Plot the hipped roof in 3D.
     for face in faces:
         for edge in zip(face, face[1:] + face[:1]):
             p1, p2 = VERTS[edge[0]], VERTS[edge[1]]
             ax.plot([p1.x, p2.x], [p1.y, p2.y], [p1.z, p2.z], 'k')
 
     ax.axis('equal')
-    ax.set_zlim(0, 30)
+    ax.set_zlim(0, 1200)
 
 
 def run_demo():
@@ -85,18 +129,33 @@ def run_demo():
     # Add a vertical slider to control the pitch.
     pitch_slider = Slider(
         ax=fig.add_axes([0.9, 0.2, 0.02, 0.6]),
-        label="Roof pitch\n(degrees)",
+        label='Roof pitch\n(degrees)',
         valmin=MIN_PITCH,
         valmax=MAX_PITCH,
         valinit=DEFAULT_PITCH,
         valstep=1.0,
-        orientation="vertical"
+        orientation='vertical'
     )
     pitch_slider.on_changed(lambda pitch: create_roof(pitch, ax))
 
     plt.show()
 
 
-if __name__ == "__main__":
+def run_benchmark():
+    start = time.perf_counter()
+    pitch = 0
+
+    while pitch < 90:
+        _ = get_roof_faces(pitch)
+        pitch += 0.1
+
+    duration = time.perf_counter() - start
+    print(f'Execution time: {round(duration, 4)} seconds.')
+
+
+if __name__ == '__main__':
     # Execute only if run as a script.
-    run_demo()
+    if len(sys.argv) == 2 and sys.argv[1] == '--benchmark':
+        run_benchmark()
+    else:
+        run_demo()