Logo_iSurf_org

Subscribe
to eAlerts

for Geom Site
Updates

CLICK HERE

 

Intersect a Segment & a Convex Polygon

by Dan Sunday

 

 

One of the most common computer graphics computations is the clipping of a line segment with a convex polygonal object. There is an efficient algorithm for determining this that works in all dimensions, and we will describe its implementation for 2D and 3D convex polygons and polyhedrons. It is also an efficient method for determining if a line segment or ray do not intersect the convex object, in applications that need to determine visibility or collision avoidance/detection.

For this, it is very useful to have convex bounding containers for more complicated objects, as discussed in Algorithm 8. The convex bounding container will have a smaller number of facets (2D edges or 3D faces) than a complicated object, which may have hundreds or thousands of them. The axis-oriented box (AOB) container has only 2n facets in n dimensional space. Also, the convex hull is the smallest convex container that can closely approximate an object. Efficient algorithms for constructing 2D convex hulls were discussed in Algorithm 10, Algorithm 11, and Algorithm 12. The hull may have many facets, but there are usually far fewer than the complicated object it bounds. For efficency in applications, it can be useful to have a hierarchy of containers, both the AOB and the hull.

Several efficient algorithms have been developed for computing the intersection (or not) of an oriented line segment with with a convex polytope. One that is frequently used in implementations, referred to as "parametric line-clipping" [Foley et al, 1996], was originally developed by [Cyrus & Beck, 1978] and later improved by [Liang & Barsky, 1984]. Their algorithm is usually described as a method for clipping a segment with a 2D rectangular window. However, it easily generalizes to 2D convex polygons and 3D convex polyhedrons, and is a very efficient algorithm in both cases. The root of this efficiency is because one only has to determine the intersections with the hyperplanes (ie: 2D lines or 3D planes) for the facets of the convex container. That is, one does not have to determine if the intersection is interior to the object’s facets. This greatly simplifies and speeds up the algorithm.

Throughout, let S be a line segment between endpoints P0 and P1. The extended line through P0 and P1 is given by the parametric equation: P(t)=P0+t(P1-P0)=P0+tv with v=P1-P0 the line direction vector. Then, the segment S contains those points P(t) with t.ge-0.le-1. This representation is valid for both 2D and 3D lines. [Note: a more detailed discussion about line equations is in Algorithm 2 about the Distance of a Point to a Line.]


Intersect a Segment and a Convex Polygon (2D)

Let a convex polygon OMEGA be given by n vertices V0,V1,..,Vn-1 going counterclockwise (ccw) around the polygon, and let Vn=V0. Also let ei be the i-th edge (line segment) ViVi+1 for i=0,n-1; and evi=(ei1,ei2)=Vi+1-Vi be the edge vector. Then, an outward-pointing normal vector for ei is given by: ni=-evi=(ei2,-ei1), where "perp_prod_symbol" is the 2D perp-operator described in Math Vector Products.

Using the technique from Algorithm 5 for Line and Segment Intersections, we first compute the intersection of the (extended) line P(t) with the extended line for a single edge ei. Consider the following diagram:

edge_clip

As indicated, intersection occurs when (P(t)-Vi)_dot_ni=0, since any vector parallel to the edge ei is perpendicular to the edge normal vector. Substituting for P(t) and solving for t, we get:

Eqn_seg-edge1

at the intersection point Ii=P(ti). The denominator here is nonzero unless the segment and the edge are parallel, and we treat that as a special case. In particular, when S is outside the edge ei, then the segment cannot intersect the convex polygon at all, and we are done. This condition can be checked by testing if the vector from Vi to P0 points to the outside of the edge. This is true when (P0-Vi)_dot_ni.gt.0, in which case processing can stop since the segment is completely outside the polygon. Otherwise, the segment S is inside the parallel edge, so we can simply ignore this edge, and continue processing the other edges.

But when they are not parallel and a unique Ii exists, then for ei, classify the associated ti as entering or leaving by the criteria:

Eqn_seg-edge2

as illustrated in the following diagram.

polygon_clip

Observe that the line P(t) must have crossed from the outside to the inside of all the edges where it is entering before it can be inside the whole polygon. This happens when it reaches the maximum value of all the entering ti values. Conversely, once P(t) has crossed from the inside to the outside across any edge where it is leaving, then it can never re-enter the polygon [Note: this is only true for convex polygons!]. So, once t gets larger than the minimum of all the leaving ti values, then it is outside the polygon for good. Put:

Eqn_seg-edge3

Then, the subsegment of the line inside the polygon will start at P(tE) and end at P(tL) for increasing t. But this will be non-empty only when 0.le.tE.le.tL.le.1. If tL.lt.tE, the segment essentially leaves the polygon before it fully enters it, and is never inside all faces at the same time.

Pseudo-Code: Intersect Segment with Polygon

The following pseudo-code efficiently encapsulates the logic of this algorithm.

    Input: a 2D segment S from point P0 to point P1
          
a 2D convex polygon OMEGA with n vertices V0,...,Vn-1,Vn=V0

    if (P0 == P1) then S is a single point, so {
        test for point inclusion of P0 in OMEGA; and
        return the test result (TRUE or FALSE);
    }

    Initialize:
        tE = 0 for the maximum entering segment parameter;
        tL = 1 for the minimum leaving segment parameter;
       dS = P1 - P0 is the segment direction vector;

    for each (edge ei = ViVi+1 of OMEGA; i=0,n-1)
    {
        Let ni = an outward normal of the edge ei;
        N = - dot product of (P0-Vi) and ni;
        D = dot product of dS and ni;
        if (D == 0) then S is parallel to the edge ei {
            if (N < 0) then P0 is outside the edge ei
                 return FALSE since S cannot intersect OMEGA;
            else S cannot enter or leave OMEGA across edge ei {
                 ignore edge ei and
                 continue to process the next edge;
            }
        }

        Put t = N / D;
        if (D < 0) then segment S is entering OMEGA across edge ei {
            New tE = max of current tE and this t
            if (tE > tL) then segment S enters OMEGA after leaving
                 return FALSE since S cannot intersect OMEGA
        }
        else (D > 0) then segment S is leaving OMEGA across edge ei {
            New tL = min of current tL and this t
            if (tL < tE) then segment S leaves OMEGA before entering
                 return FALSE since S cannot intersect OMEGA
        }
    }

    Output: [Note: to get here, one must have tE <= tL]
    there is a valid intersection of S with OMEGA
        from the entering point: P(tE) = P0 + tE * dS
        to the leaving point:    P(tL) = P0 + tL * dS
    return TRUE

 

A C++ implementation intersect2D_SegPoly() is given below.


Intersect a Segment and a Convex Polyhedron (3D)

The 2D computations and algorithm extend to 3D with very few changes. The main difference concerns specifying the polyhedron data structure. We let a polyhedron OMEGA be given by a list of n (planar) polygon faces Fi for i=0,n-1. We do not assume any relationship between succesive faces Fi and Fi+1, such as having a common edge or vertex (like adjacent edges of a polygon). It is just an unstructured list of all the distinct faces. However, to make the algorithm work, we need some extra information about each face, namely:

  1. Each face Fi lies completely in a plane.
  2. Vi = a point in the plane of face Fi (e.g.: a vertex of Fi).
  3. ni = an outward-pointing normal vector to Fi ; that is, a normal vector pointing to the side of Fi that is exterior to the polyhedron OMEGA.

This information can be precomputed from any decent data structure for a polyhedron.

Again, the 3D line segment S = P0P1 is given by a parametric equation P(t). For the intersection of the extended line segment with the plane of a specific face Fi, consider the following diagram.

face_clip

This is exactly the same condition as in the 2D case; namely, the intersection point occurs when (P(t)-Vi)_dot_ni=0 as shown. Again, solving this gives:

Eqn_seg-edge1

When the denominator is zero, then the segment S is parallel to the face Fi which is treated exactly the same as before. Nevertheless, when ti exists, there is a unique intersection point Ii=P(ti).

Further, since ni is an outward-pointing normal vector to the face plane, we have exactly the same criteria for classifying each ti as entering or leaving; namely:

Eqn_seg-face2

And again we compute:

Eqn_seg-edge3

Finally, the part of the segment inside the convex polyhedron is the subsegment from P(tE) to P(tL) when 0.le.tE.le.tL.le.1. Again, if tL.lt.tE, then the segment is completely outside the polyhedron.

Pseudo-Code: Intersect Segment with Polyhedron

The 3D pseudo-code is almost the same as for 2D polygons, but with the obvious modifications.

    Input: a 3D segment S from point P0 to point P1
          
a 3D convex polyhedron OMEGA with n faces F0,...,Fn-1 and
                Vi = a vertex for each face Fi
                ni = an outward normal vector for each face Fi

    if (P0 == P1) then S is a single point, so {
        test for point inclusion of P0 in OMEGA; and
        return the test result (TRUE or FALSE);
    }

    Initialize:
        tE = 0 for the maximum entering segment parameter;
        tL = 1 for the minimum leaving segment parameter;
        dS = P1 - P0 is the segment direction vector;

    for (each face Fi of OMEGA; i=0,n-1)
    {
        N = - dot product of (P0 - Vi) and ni;
        D = dot product of dS and ni;
        if (D == 0) then S is parallel to the face Fi {
            if (N < 0) then P0 is outside the face Fi
                 return FALSE since S cannot intersect OMEGA;
            else S cannot enter or leave OMEGA across face Fi {
                 ignore face Fi and
                 continue to process the next face;
            }
        }

        Put t = N / D;
        if (D < 0) then segment S is entering OMEGA across face Fi {
            New tE = max of current tE and this t
            if (tE > tL) then segment S enters OMEGA after leaving
                 return FALSE since S cannot intersect OMEGA
        }
        else (D > 0) then segment S is leaving OMEGA across face Fi {
            New tL = min of current tL and this t
            if (tL < tE) then segment S leaves OMEGA before entering
                 return FALSE since S cannot intersect OMEGA
        }
    }

    Output: [Note: to get here, one must have tE <= tL]
    there is a valid intersection of S with OMEGA
        from the entering point: P(tE) = P0 + tE * dS
        to the leaving point:    P(tL) = P0 + tL * dS
    return TRUE


Beyond 3D, this algorithm can be easily generalized to n-dimensional space for the intersection of a parametric segment with the interior of a convex polytope. The details are basically the same.


Implementation

Here is a sample "C++" implementation of this algorithm.

// Copyright 2001 softSurfer, 2012 Dan Sunday
// This code may be freely used and modified for any purpose
// providing that this copyright notice is included with it.
// SoftSurfer makes no warranty for this code, and cannot be held
// liable for any real or imagined damage resulting from its use.
// Users of this code must verify correctness for their application.
 

// Assume that classes are already given for the objects:
//    Point and Vector with
//        coordinates {float x, y;}
//        operators for:
//            == to test equality
//            != to test inequality
//            Point  = Point Vector
//            Vector = Point - Point
//            Vector = Vector Vector
//            Vector = Scalar * Vector    (scalar product)
//    Segment with defining endpoints {Point P0, P1;}
//===================================================================
 

#define TRUE       1
#define FALSE      0
#define SMALL_NUM  0.00000001 // anything that avoids division overflow

#define dot(u,v)   ((u).x * (v).x + (u).y * (v).y)   // 2D dot product
#define perp(u,v)  ((u).x * (v).y - (u).y * (v).x)   // 2D perp product
 


// intersect2D_SegPoly(): intersect a 2D segment with a convex polygon
//    Input:  S = 2D segment to intersect with the convex polygon V[]
//            n = number of 2D points in the polygon
//            V[] = array of n+1 vertex points with V[n] = V[0]
//      Note: The polygon MUST be convex and
//                have vertices oriented counterclockwise (ccw).
//            This code does not check for and verify these conditions.
//    Output: *IS = the intersection segment (when it exists)
//    Return: FALSE = no intersection
//            TRUE  = a valid intersection segment exists
int
intersect2D_SegPoly( Segment S, Point* V, int n, Segment* IS)
{
    if (S.P0 == S.P1) {         // the segment S is a single point
        // test for inclusion of S.P0 in the polygon
        *IS = S;                // same point if inside polygon
        return cn_PnPoly( S.P0, V, n );
    }

    float  tE = 0;              // the maximum entering segment parameter
    float  tL = 1;              // the minimum leaving segment parameter
    float  t, N, D;             // intersect parameter t = N / D
    Vector dS = S.P1- S.P0;     // the  segment direction vector
    Vector e;                   // edge vector
    // Vector ne;               // edge outward normal (not explicit in code)

    for (int i=0; i < n; i++)   // process polygon edge V[i]V[i+1]
    {
        e = V[i+1] - V[i];
        N = perp(e, S.P0-V[i]); // = -dot(ne, S.P0 - V[i])
        D = -perp(e, dS);       // = dot(ne, dS)
        if (fabs(D) < SMALL_NUM) {  // S is nearly parallel to this edge
            if (N < 0)              // P0 is outside this edge, so
                 return FALSE;      // S is outside the polygon
            else                    // S cannot cross this edge, so
                 continue;          // ignore this edge
        }

        t = N / D;
        if (D < 0) {            // segment S is entering across this edge
            if (t > tE) {       // new max tE
                 tE = t;
                 if (tE > tL)   // S enters after leaving polygon
                     return FALSE;
            }
        }
        else {                  // segment S is leaving across this edge
            if (t < tL) {       // new min tL
                 tL = t;
                 if (tL < tE)   // S leaves before entering polygon
                     return FALSE;
            }
        }
    }

    // tE <= tL implies that there is a valid intersection subsegment
    IS->P0 = S.P0 + tE * dS;   // = P(tE) = point where S enters polygon
    IS->P1 = S.P0 + tL * dS;   // = P(tL) = point where S leaves polygon
    return TRUE;
}
 


References

M. Cyrus & J. Beck, "Generalized Two- and Three-Dimensional Clipping", Computers and Graphics 3(1), 23-28 (1978)

James Foley, Andries van Dam, Steven Feiner & John Hughes, "A Parametric Line-Clipping Algorithm" in Computer Graphics (3rd Edition) (2013)

Y-D. Liang & B. Barsky, "A New Concept and Method for Line Clipping", ACM TOG 3(1), 1-22 (1984)

 

© Copyright 2012 Dan Sunday, 2001 softSurfer