Search code examples
c#mathgeometryline-segment

Segment Segment intersection in 3D


Let's suppose I have a triangle defined by three coordinates A, B and C and a line segment defined by two coordinates D and E in 3D euclidean space. Let's further suppose that the segment is intersecting the triangle. Now, I need to find out whether the line segment intersects the triangle at its "border". So my idea was to take each edge of the triangle (AB, BC and CA) and test whether one of them intersects with the segment DE. The problem is that I haven't found a solution for segment/segment intersection. I have found this one but it doesn't work. Help very appreciated.

Edit Based on MBo's answer, I have implemented a C# method:

// return: 0 -> segment parallel to plane
// return: 1 -> segment intersects an edge of the face
// return: -1 -> segments intersects the interior or the triangle or not at all
public static int SegmentTriangleIntersection(Vector3D a, Vector3D b, Vector3D c, Vector3D d, Vector3D e)
{
    var ed = e - d;
    var ba = b - a;
    var ca = c - a;
    var ad = a - d;

    var det = (ed.X*-ba.Y*-ca.Z) + (-ba.X*-ca.Y*ed.Z) + (-ca.X*ed.Y*-ba.Z) - (ed.Z*-ba.Y*-ca.X) - (-ba.Z*-ca.Y*ed.X) -
                  (-ca.Z*ed.Y*-ba.X);

    if (Vector3D.IsNearlyEqual(det, 0)) // segment parallel to triangle
        return 0;

    var det_t = (ad.X * -ba.Y * -ca.Z) + (-ba.X * -ca.Y * ad.Z) + (-ca.X * ad.Y * -ba.Z) - (ad.Z * -ba.Y * -ca.X) - (-ba.Z * -ca.Y * ad.X) -
                    (-ca.Z * ad.Y * -ba.X);

    var t = det_t/det;
    if (t >= 0 & t <= 1) // segment intersects plane of triangle
    {
        var det_u = (ed.X*ad.Y*-ca.Z) + (ad.X*-ca.Y*ed.Z) + (-ca.X*ed.Y*ad.Z) - (ed.Z*ad.Y*-ca.X) - (ad.Z*-ca.Y*ed.X) -
                        (-ca.Z*ed.Y*ad.X);
        var u = det_u/det;
        var det_v = (ed.X*-ba.Y*ad.Z) + (-ba.X*ad.Y*ed.Z) + (ad.X*ed.Y*-ba.Z) - (ed.Z*-ba.Y*ad.X) -
                        (-ba.Z*ad.Y*ed.X)-(ad.Z*ed.Y*-ba.X);
        var v = det_v/det;

        if (Vector3D.IsNearlyEqual(u, 0) && v >= 0 && v <= 1)
            return 1;
        if (Vector3D.IsNearlyEqual(v, 0) && u >= 0 && u <= 1)
            return 1;
        if (Vector3D.IsNearlyEqual(v + u, 1) && u >= 0 && v >= 0)
            return 1;
    }
    return -1;
}

Solution

  • Let's use parametric equations for coordinates (vectors in bold):

    Segment DE:
    sDE(t) = D + t * (E - D) = D + t * ED
    //ED = E - D

    Segment AB:
    sAB(u) = A + u * (B - A) = A + u * BA

    Segment AC:
    sAC(v) = A + v * (C - A) = A + u * CA

    Plane of ABC triangle is described by bi-parametric equation:
    pABC(u,v) = A + u * BA + v * CA

    Point in the plane with coordinates (u,v) lies inside the triangle, if

    u in range [0..1]
    and 
    v in range [0..1]
    and
    u+v is in range [0..1]
    

    Point with coordinates (u,v) lies at the triangle edge, if [edge condition]

    u = 0 and v in range [0..1]    //at AC edge
    or
    v = 0 and u in range [0..1]    //at AC edge
    or
    v + u = 1 and u in range (0..1) and v in range (0..1)      //at BC edge
    

    Now let's write equation for intersection of DE segment and ABC plane

    sDE(t) = pABC(u,v)
    D + t * ED = A + u * BA + v * CA

    We can solve the last equation (system of 3 linear equations for every coordinate) and find parameters t,u,v

    D.X + t * ED.X = A.X + u * BA.X + v * CA.X
    D.Y + t * ED.Y = A.Y + u * BA.Y + v * CA.Y
    D.Z + t * ED.Z = A.Z + u * BA.Z + v * CA.Z
    
    t * ED.X - u * BA.X - v * CA.X = A.X - D.X
    t * ED.Y - u * BA.Y - v * CA.Y = A.Y - D.Y
    t * ED.Z - u * BA.Z - v * CA.Z = A.Z - D.Z
    

    At first find determinant of this system - if it is zero, segment is parallel to the plane.

    If not, check for parameter t. Segment intersects the plane if t is in range [0..1]

    If yes, check whether u, v parameters comply to edge condition above