bridge edges

bridge_edges.py

'''
in verts v
in edges s
out new_verts v
out new_edges s
'''

import numpy as np

def shortest_edge_between_edges(edge1, edge2):
    """
    Finds the shortest edge connecting two edges in 3D space.
    
    Parameters:
        edge1: Tuple of two points (P1, Q1), where each point is a 3D vector.
        edge2: Tuple of two points (P2, Q2), where each point is a 3D vector.

    Returns:
        Tuple containing:
            - Closest point on edge1
            - Closest point on edge2
            - Distance between the two closest points
    """
    def point_to_edge_closest_point(P, A, B):
        """Find the closest point on the edge AB to the point P."""
        AB = B - A
        t = np.dot(P - A, AB) / np.dot(AB, AB)
        t = np.clip(t, 0, 1)  # Clamp t to [0, 1]
        return A + t * AB
    
    P1, Q1 = np.array(edge1[0]), np.array(edge1[1])
    P2, Q2 = np.array(edge2[0]), np.array(edge2[1])

    # Vectors for edges
    U1 = Q1 - P1
    U2 = Q2 - P2

    # Direction vectors' cross product
    W0 = P1 - P2
    A = np.dot(U1, U1)
    B = np.dot(U1, U2)
    C = np.dot(U2, U2)
    D = np.dot(U1, W0)
    E = np.dot(U2, W0)

    denominator = A * C - B * B

    if np.isclose(denominator, 0):
        # Edges are nearly parallel, find the closest point to one edge from the other
        P1_on_edge2 = point_to_edge_closest_point(P1, P2, Q2)
        Q1_on_edge2 = point_to_edge_closest_point(Q1, P2, Q2)
        P2_on_edge1 = point_to_edge_closest_point(P2, P1, Q1)
        Q2_on_edge1 = point_to_edge_closest_point(Q2, P1, Q1)

        # Compute distances and choose the smallest
        points = [
            (P1, P1_on_edge2),
            (Q1, Q1_on_edge2),
            (P2_on_edge1, P2),
            (Q2_on_edge1, Q2)
        ]
        distances = [np.linalg.norm(p1 - p2) for p1, p2 in points]
        idx = np.argmin(distances)
        return points[idx][0], points[idx][1], distances[idx]

    # Compute parameters s and t
    s = (B * E - C * D) / denominator
    t = (A * E - B * D) / denominator
    s = np.clip(s, 0, 1)  # Clamp s to [0, 1]
    t = np.clip(t, 0, 1)  # Clamp t to [0, 1]

    # Closest points on each edge
    closest_point_edge1 = P1 + s * U1
    closest_point_edge2 = P2 + t * U2

    return closest_point_edge1, closest_point_edge2

inverts = verts[0]
inedges = edges[0]
v1, v2 = inverts[inedges[0][0]], inverts[inedges[0][1]]
v3, v4 = inverts[inedges[1][0]], inverts[inedges[1][1]]

p1, p2 = shortest_edge_between_edges((v1, v2), (v3, v4))

new_verts.append([p1.round(4).tolist(), p2.round(4).tolist()])
new_edges.append([[0, 1]])