The intersection of geometric primitives is a fundamental construct in many computer graphics and modeling applications ([Foley et al, 1996], [O'Rourke, 1998]). Here we look at the algorithms for the simplest 2D and 3D linear primitives: lines, segments and planes.
Line and Segment Intersections
For computing intersections of lines and segments in 2D and 3D, it is best to use the parametric equation representation for lines. Other representations are discussed in Algorithm 2 about the Distance of a Point to a Line. It is shown there how to convert from other representations to the parametric one.
In any dimension, the parametric equation of a line defined by two points P_{0 }and P_{1} can be represented as:
,
where the parameter s is a real number and is a line direction vector. Using this representation , , and when , P(s) is a point on the finite segment P_{0}P_{1} where s is the fraction of P(s)'s distance along the segment. That is, s = d(P_{0},P(s)) / d(P_{0},P_{1}). Further, if s < 0 then P(s) is outside the segment on the P_{0} side, and if s > 1 then P(s) is outside the segment on the P_{1} side.
Let two lines be given by: and , either or both of which could be infinite, a ray, or a finite segment.
Parallel Lines
These lines are parallel when and only when their directions are collinear, namely when the two vectors and are linearly related as u = av for some real number a. For and , this means that all ratios have the value a, or that _{ }for all i_{}. This is equivalent to the conditions that all . In 2D, with and , this is the perp product [Hill, 1994] condition that where , the perp operator, is perpendicular to u. This condition says that two vectors in the Euclidean plane are parallel when they are both perpendicular to the same direction vector . When true, the two associated lines are either coincident or do not intersect at all.
Coincidence is easily checked by testing if a point on one line, say P_{0}, also lies on the other line Q(t). That is, there exists a number t_{0} such that: P_{0 }= Q(t_{0}) = Q_{0 }+ t_{0}v. If w = (w_{i}) = P_{0} – Q_{0}, then this is equivalent to the condition that w = t_{0}v for some t_{0} which is the same as for all i. In 2D, this is another perp product condition: . If this condition holds, one has , and the infinite lines are coincident. And if one line (but not the other) is a finite segment, then it is the coincident intersection. However, if both lines are finite segments, then they may (or may not) overlap. In this case, solve for t_{0} and t_{1} such that P_{0 }= Q(t_{0}) and P_{1 }= Q(t_{1}). If the segment intervals [t_{0},t_{1}] and [0,1] are disjoint, there is no intersection. Otherwise, intersect the intervals (using max and min operations) to get . Then the intersection segment is . This works in any dimension.
NonParallel Lines
When the two lines or segments are not parallel, they might intersect in a unique point. In 2D Euclidean space, infinite lines always intersect. In higher dimensions they usually miss each other and do not intersect. But if they intersect, then their linear projections onto a 2D plane will also intersect. So, one can simply restrict to two coordinates, for which u and v are not parallel, compute the 2D intersection point I at P(s_{I}) and Q(t_{I}) for those two coordinates, and then test if P(s_{I}) = Q(t_{I}) for all coordinates. To compute the 2D intersection point, consider the two lines and the associated vectors in the diagram:
To determine s_{I}, we have the vector equality ^{ }where . At the intersection, the vector is perpendicular to , and this is equivalent to the perp product condition that . Solving this equation, we get:
Note that the denominator only when the lines are parallel as previously discussed. Similarly, solving for Q(t_{I}), we get:
The denominators are the same up to sign, since , and should only be computed once if we want to know both s_{I} and t_{I}. However, knowing either is enough to get the intersection point I = P(s_{I}) = Q(t_{I}).
Further, if one of the two lines is a finite segment (or a ray), say P_{0}P_{1}, then the intersect point is in the segment only when (or for a ray). If both lines are segments, then both solution parameters, s_{I} and t_{I}, must be in the [0,1] interval for the segments to intersect. Although this sounds simple enough, the code for the intersection of two segments is a bit delicate since many special cases need to be checked (see our implementation intersect2D_2Segments()).
Plane Intersections
Planes are represented as described in Algorithm 4, see Planes.
LinePlane Intersection
In 3D, a line L is either parallel to a plane P or intersects it in a single point. Let L be given by the parametric equation: , and the plane P be given by a point V_{0} on it and a normal vector . We first check if L is parallel to P by testing if , which means that the line direction vector u is perpendicular to the plane normal n. If this is true, then L and P are parallel and either never intersect or else L lies totally in the plane P . Disjointness or coincidence can be determined by testing whether any one specific point of L, say P_{0}, is contained in P , that is whether it satisfies the implicit line equation: .
If the line and plane are not parallel, then L and P intersect in a unique point P(s_{I}) which is computed using a method similar to the one for the intersection of two lines in 2D. Consider the diagram:
At the intersect point, the vector is perpendicular to n,^{ }where . This is equivalent to the dot product condition: . Solving we get:
If the line L is a finite segment from P_{0} to P_{1}, then one just has to check that to verify that there is an intersection between the segment and the plane. For a positive ray, there is an intersection with the plane when .
Intersection of 2 Planes
In 3D, two planes P_{1} and P_{2} are either parallel or they intersect in a single straight line L. Let P _{i} (i = 1,2) be given by a point V_{i} and a normal vector n_{i}, and have an implicit equation: n_{i }· P + d_{i }= 0, where P = (x, y, z). The planes P_{1} and P_{2} are parallel whenever their normal vectors n_{1} and n_{2} are parallel, and this is equivalent to the condition that: n_{1} x n_{2 }= 0. In software, one would test if where division by would cause overflow, and uses this as the robust condition for n_{1} and n_{2} to be parallel in practice. When not parallel, is a direction vector for the intersection line L since u is perpendicular to both n_{1} and n_{2}, and thus is parallel to both planes as shown in the following diagram. If u is small, it is best to scale it so that u = 1, making u a unit direction vector.
After computing (3 adds + 6 multiplies), to fully determine the intersection line, we still need to find a specific point on it. That is, we need to find a point that lies in both planes. We can do this by finding a common solution of the implicit equations for P_{1} and P_{2}. But there are only two equations in the 3 unknowns since the point P_{0} can lie anywhere on the 1dimensional line L. So we need another constraint to solve for a specific P_{0}. There are a number of ways this could be done:
(A) Direct Linear Equation. One could set one coordinate to zero, say z_{0 }= 0, and then solve for the other two. But this will only work when L intersects the plane z_{0 }= 0. This is will be true when the zcoordinate u_{z} of is nonzero. So, one must first select a nonzero coordinate of u, and then set the corresponding coordinate of P_{0} to 0. Further, one should choose the coordinate with the largest absolute value, as this will give the most robust computations. Suppose that , then lies on L. Solving the two equations: and , one gets:
and the parametric equation of L is:
The denominator here is equal to the nonzero 3rd coordinate of u. So, ignoring the test for a large nonzero coordinate, and counting division as a multiplication, the total number of operations for this solution = 5 adds + 13 multiplies.
(B) Line Intersect Point. If one knows a specific line in one plane (for example, two points in the plane), and this line intersects the other plane, then its point of intersection, I, will lie in both planes. Thus, it is on the line of intersection for the two planes, and the parametric equation of L is: P(s) = I + s (n_{1} x n_{2}). To compute and the intersection point (given the line), the total number of operations = 11 adds + 19 multiplies.
One way of constructing a line in one plane that must intersect the other plane is to project one plane's normal vector onto the other plane. This gives a line that must always be orthogonal to the line of the planes' intersection. So, the projection of n_{2} on P_{1} defines a line that intersects P_{2} in the sought for point P_{0} on L. More specifically, project the two points 0 = (0,0,0) and n_{2 }= (nx_{2}, ny_{2}, nz_{2}) to P_{1}(0) and P_{1}(n_{2}) respectively. Then the projected line in P_{1} is L_{1}: Q(t) = P_{1}(0) + t (P_{1}(n_{2}) – P_{1}(0)), and intersection of it with P_{2} can be computed. In the most efficient case, where both n_{1} and n_{2} are unit normal vectors and the constant P_{1}(0) is prestored, the total operations = 17 adds + 22 multiplies.
(C) 3 Plane Intersect Point. Another method selects a third plane P_{3} with an implicit equation n_{3 }· P = 0 where n_{3 }= n_{1} x n_{2} and d_{3 }= 0 (meaning it passes through the origin). This always works since: (1) L is perpendicular to P_{3} and thus intersects it, and (2) the vectors n_{1}, n_{2}, and n_{3} are linearly independent. Thus the planes P_{1}, P_{2} and P_{3} intersect in a unique point P_{0} which must be on L. Using the formula for the intersection of 3 planes (see the next section), where d_{3 }= 0 for P_{3}, we get:
and the parametric equation of L is:
The number of operations for this solution = 11 adds + 23 multiplies.
The bottom line is that the most efficient method is the direct solution (A) that uses only 5 adds + 13 multiplies to compute the equation of the intersection line.
Intersection of 3 Planes
In 3D, three planes P_{1}, P_{2} and P_{3} can intersect (or not) in the following ways:
Geometric Relation

Intersection

Algebraic Condition

All 3 planes are parallel


for all

They coincide

A plane

_{}

They are disjoint

None

_{}




Only two planes are parallel, and the 3rd plane cuts each in a line [Note: the 2 parallel planes may coincide]

2 parallel lines
[planes coincide => 1 line]

Only one for




No two planes are parallel, so pairwise they intersect in 3 lines


for all

The intersect lines are parallel



They all coincide

1 line

Test a point of one line with another line

They are disjoint

3 parallel lines

Same test, but fails => do not coincide

No intersect lines are parallel

A unique point



As shown in the illustrations:
One should first test for the most frequent case of a unique intersect point, namely that , since this excludes all the other cases. When the intersection is a unique point, it is given by the formula:
which can verified by showing that this P_{0} satisfies the parametric equations for all planes P_{1}, P_{2} and P_{3}.
However, there can be a problem with the robustness of this computation when the denominator is very small. In that case, it would be best to get a robust line of intersection for two of the planes, and then compute the point where this line intersects the third plane.
Implementations
Here are some sample "C++" implementations of these algorithms.
// 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, z;} // operators for: // == to test equality // != to test inequality // Point = Point ± Vector // Vector = Point  Point // Vector = Scalar * Vector (scalar product) // Vector = Vector * Vector (3D cross product) // Line and Ray and Segment with defining points {Point P0, P1;} // (a Line is infinite, Rays and Segments start at P0) // (a Ray extends beyond P1, but a Segment ends at P1) // Plane with a point and a normal {Point V0; Vector n;} //===================================================================
#define SMALL_NUM 0.00000001 // anything that avoids division overflow // dot product (3D) which allows vector operations in arguments #define dot(u,v) ((u).x * (v).x + (u).y * (v).y + (u).z * (v).z) #define perp(u,v) ((u).x * (v).y  (u).y * (v).x) // perp product (2D)
// intersect2D_2Segments(): find the 2D intersection of 2 finite segments // Input: two finite segments S1 and S2 // Output: *I0 = intersect point (when it exists) // *I1 = endpoint of intersect segment [I0,I1] (when it exists) // Return: 0=disjoint (no intersect) // 1=intersect in unique point I0 // 2=overlap in segment from I0 to I1 int intersect2D_2Segments( Segment S1, Segment S2, Point* I0, Point* I1 ) { Vector u = S1.P1  S1.P0; Vector v = S2.P1  S2.P0; Vector w = S1.P0  S2.P0; float D = perp(u,v);
// test if they are parallel (includes either being a point) if (fabs(D) < SMALL_NUM) { // S1 and S2 are parallel if (perp(u,w) != 0  perp(v,w) != 0) { return 0; // they are NOT collinear } // they are collinear or degenerate // check if they are degenerate points float du = dot(u,u); float dv = dot(v,v); if (du==0 && dv==0) { // both segments are points if (S1.P0 != S2.P0) // they are distinct points return 0; *I0 = S1.P0; // they are the same point return 1; } if (du==0) { // S1 is a single point if (inSegment(S1.P0, S2) == 0) // but is not in S2 return 0; *I0 = S1.P0; return 1; } if (dv==0) { // S2 a single point if (inSegment(S2.P0, S1) == 0) // but is not in S1 return 0; *I0 = S2.P0; return 1; } // they are collinear segments  get overlap (or not) float t0, t1; // endpoints of S1 in eqn for S2 Vector w2 = S1.P1  S2.P0; if (v.x != 0) { t0 = w.x / v.x; t1 = w2.x / v.x; } else { t0 = w.y / v.y; t1 = w2.y / v.y; } if (t0 > t1) { // must have t0 smaller than t1 float t=t0; t0=t1; t1=t; // swap if not } if (t0 > 1  t1 < 0) { return 0; // NO overlap } t0 = t0<0? 0 : t0; // clip to min 0 t1 = t1>1? 1 : t1; // clip to max 1 if (t0 == t1) { // intersect is a point *I0 = S2.P0 + t0 * v; return 1; }
// they overlap in a valid subsegment *I0 = S2.P0 + t0 * v; *I1 = S2.P0 + t1 * v; return 2; }
// the segments are skew and may intersect in a point // get the intersect parameter for S1 float sI = perp(v,w) / D; if (sI < 0  sI > 1) // no intersect with S1 return 0;
// get the intersect parameter for S2 float tI = perp(u,w) / D; if (tI < 0  tI > 1) // no intersect with S2 return 0;
*I0 = S1.P0 + sI * u; // compute S1 intersect point return 1; } //===================================================================
// inSegment(): determine if a point is inside a segment // Input: a point P, and a collinear segment S // Return: 1 = P is inside S // 0 = P is not inside S int inSegment( Point P, Segment S) { if (S.P0.x != S.P1.x) { // S is not vertical if (S.P0.x <= P.x && P.x <= S.P1.x) return 1; if (S.P0.x >= P.x && P.x >= S.P1.x) return 1; } else { // S is vertical, so test y coordinate if (S.P0.y <= P.y && P.y <= S.P1.y) return 1; if (S.P0.y >= P.y && P.y >= S.P1.y) return 1; } return 0; } //===================================================================
// intersect3D_SegmentPlane(): find the 3D intersection of a segment and a plane // Input: S = a segment, and Pn = a plane = {Point V0; Vector n;} // Output: *I0 = the intersect point (when it exists) // Return: 0 = disjoint (no intersection) // 1 = intersection in the unique point *I0 // 2 = the segment lies in the plane int intersect3D_SegmentPlane( Segment S, Plane Pn, Point* I ) { Vector u = S.P1  S.P0; Vector w = S.P0  Pn.V0;
float D = dot(Pn.n, u); float N = dot(Pn.n, w);
if (fabs(D) < SMALL_NUM) { // segment is parallel to plane if (N == 0) // segment lies in plane return 2; else return 0; // no intersection } // they are not parallel // compute intersect param float sI = N / D; if (sI < 0  sI > 1) return 0; // no intersection
*I = S.P0 + sI * u; // compute segment intersect point return 1; } //===================================================================
// intersect3D_2Planes(): find the 3D intersection of two planes // Input: two planes Pn1 and Pn2 // Output: *L = the intersection line (when it exists) // Return: 0 = disjoint (no intersection) // 1 = the two planes coincide // 2 = intersection in the unique line *L int intersect3D_2Planes( Plane Pn1, Plane Pn2, Line* L ) { Vector u = Pn1.n * Pn2.n; // cross product float ax = (u.x >= 0 ? u.x : u.x); float ay = (u.y >= 0 ? u.y : u.y); float az = (u.z >= 0 ? u.z : u.z);
// test if the two planes are parallel if ((ax+ay+az) < SMALL_NUM) { // Pn1 and Pn2 are near parallel // test if disjoint or coincide Vector v = Pn2.V0  Pn1.V0; if (dot(Pn1.n, v) == 0) // Pn2.V0 lies in Pn1 return 1; // Pn1 and Pn2 coincide else return 0; // Pn1 and Pn2 are disjoint }
// Pn1 and Pn2 intersect in a line // first determine max abs coordinate of cross product int maxc; // max coordinate if (ax > ay) { if (ax > az) maxc = 1; else maxc = 3; } else { if (ay > az) maxc = 2; else maxc = 3; }
// next, to get a point on the intersect line // zero the max coord, and solve for the other two Point iP; // intersect point float d1, d2; // the constants in the 2 plane equations d1 = dot(Pn1.n, Pn1.V0); // note: could be prestored with plane d2 = dot(Pn2.n, Pn2.V0); // ditto
switch (maxc) { // select max coordinate case 1: // intersect with x=0 iP.x = 0; iP.y = (d2*Pn1.n.z  d1*Pn2.n.z) / u.x; iP.z = (d1*Pn2.n.y  d2*Pn1.n.y) / u.x; break; case 2: // intersect with y=0 iP.x = (d1*Pn2.n.z  d2*Pn1.n.z) / u.y; iP.y = 0; iP.z = (d2*Pn1.n.x  d1*Pn2.n.x) / u.y; break; case 3: // intersect with z=0 iP.x = (d2*Pn1.n.y  d1*Pn2.n.y) / u.z; iP.y = (d1*Pn2.n.x  d2*Pn1.n.x) / u.z; iP.z = 0; } L>P0 = iP; L>P1 = iP + u; return 2; } //===================================================================
References
James Foley, Andries van Dam, Steven Feiner & John Hughes, "Clipping Lines" in Computer Graphics (3rd Edition) (2013)
Joseph O'Rourke, "Search and Intersection" in Computational Geometry in C (2nd Edition) (1998)
