3D Line-Plane Intersection
If given a line (represented by either a vector or two points on the line) how do I find the point at which the line intersects a plane? I've found loads of resources on this but I can't understand the equations there (they don't seem to be standard algebraic). I would like an equation (no matter how long) that can be interpreted by a standard programming language (I'm using Java).
Solution 1:
Here is a Python example which finds the intersection of a line and a plane.
Where the plane can be either a point and a normal, or a 4d vector (normal form), In the examples below (code for both is provided).
Also note that this function calculates a value representing where the point is on the line, (called fac in the code below).
You may want to return this too, because values from 0 to 1 intersect the line segment - which may be useful for the caller.
Other details noted in the code-comments.
Note: This example uses pure functions, without any dependencies - to make it easy to move to other languages.
With a Vector data type and operator overloading, it can be more concise (included in example below).
# intersection function
def isect_line_plane_v3(p0, p1, p_co, p_no, epsilon=1e-6):
"""
p0, p1: Define the line.
p_co, p_no: define the plane:
p_co Is a point on the plane (plane coordinate).
p_no Is a normal vector defining the plane direction;
(does not need to be normalized).
Return a Vector or None (when the intersection can't be found).
"""
u = sub_v3v3(p1, p0)
dot = dot_v3v3(p_no, u)
if abs(dot) > epsilon:
# The factor of the point between p0 -> p1 (0 - 1)
# if 'fac' is between (0 - 1) the point intersects with the segment.
# Otherwise:
# < 0.0: behind p0.
# > 1.0: infront of p1.
w = sub_v3v3(p0, p_co)
fac = -dot_v3v3(p_no, w) / dot
u = mul_v3_fl(u, fac)
return add_v3v3(p0, u)
# The segment is parallel to plane.
return None
# ----------------------
# generic math functions
def add_v3v3(v0, v1):
return (
v0[0] + v1[0],
v0[1] + v1[1],
v0[2] + v1[2],
)
def sub_v3v3(v0, v1):
return (
v0[0] - v1[0],
v0[1] - v1[1],
v0[2] - v1[2],
)
def dot_v3v3(v0, v1):
return (
(v0[0] * v1[0]) +
(v0[1] * v1[1]) +
(v0[2] * v1[2])
)
def len_squared_v3(v0):
return dot_v3v3(v0, v0)
def mul_v3_fl(v0, f):
return (
v0[0] * f,
v0[1] * f,
v0[2] * f,
)
If the plane is defined as a 4d vector (normal form), we need to find a point on the plane, then calculate the intersection as before (see p_co assignment).
def isect_line_plane_v3_4d(p0, p1, plane, epsilon=1e-6):
u = sub_v3v3(p1, p0)
dot = dot_v3v3(plane, u)
if abs(dot) > epsilon:
# Calculate a point on the plane
# (divide can be omitted for unit hessian-normal form).
p_co = mul_v3_fl(plane, -plane[3] / len_squared_v3(plane))
w = sub_v3v3(p0, p_co)
fac = -dot_v3v3(plane, w) / dot
u = mul_v3_fl(u, fac)
return add_v3v3(p0, u)
return None
For further reference, this was taken from Blender and adapted to Python.
isect_line_plane_v3() in math_geom.c
For clarity, here are versions using the mathutils API (which can be modified for other math libraries with operator overloading).
# point-normal plane
def isect_line_plane_v3(p0, p1, p_co, p_no, epsilon=1e-6):
u = p1 - p0
dot = p_no * u
if abs(dot) > epsilon:
w = p0 - p_co
fac = -(plane * w) / dot
return p0 + (u * fac)
return None
# normal-form plane
def isect_line_plane_v3_4d(p0, p1, plane, epsilon=1e-6):
u = p1 - p0
dot = plane.xyz * u
if abs(dot) > epsilon:
p_co = plane.xyz * (-plane[3] / plane.xyz.length_squared)
w = p0 - p_co
fac = -(plane * w) / dot
return p0 + (u * fac)
return None
Solution 2:
Here is a method in Java that finds the intersection between a line and a plane. There are vector methods that aren't included but their functions are pretty self explanatory.
/**
* Determines the point of intersection between a plane defined by a point and a normal vector and a line defined by a point and a direction vector.
*
* @param planePoint A point on the plane.
* @param planeNormal The normal vector of the plane.
* @param linePoint A point on the line.
* @param lineDirection The direction vector of the line.
* @return The point of intersection between the line and the plane, null if the line is parallel to the plane.
*/
public static Vector lineIntersection(Vector planePoint, Vector planeNormal, Vector linePoint, Vector lineDirection) {
if (planeNormal.dot(lineDirection.normalize()) == 0) {
return null;
}
double t = (planeNormal.dot(planePoint) - planeNormal.dot(linePoint)) / planeNormal.dot(lineDirection.normalize());
return linePoint.plus(lineDirection.normalize().scale(t));
}
Solution 3:
You'll need to consider three cases:
- Plane is parallel to line, and line does not lie in plane (no intersection)
- Plane is not parallel to line (one point of intersection)
- Plane contains the line (line intersects at every point on it)
You can express the line in paramaterized form, like here:
http://answers.yahoo.com/question/index?qid=20080830195656AA3aEBr
The first few pages of this lecture do the same for the plane:
http://math.mit.edu/classes/18.02/notes/lecture5compl-09.pdf
If the normal to the plane is perpendicular to the direction along the line, then you have an edge case and need to see whether it intersects at all, or lies within the plane.
Otherwise, you have one point of intersection, and can solve for it.
I know this isn't code but to get a robust solution you'll probably want to put this in the context of your application.
EDIT: Here's an example for which there's exactly one point of intersection. Say you start with the parameterized equations in the first link:
x = 5 - 13t
y = 5 - 11t
z = 5 - 8t
The parameter t can be anything. The (infinite) set of all (x, y, z) that satisfy these equations comprise the line. Then, if you have the equation for a plane, say:
x + 2y + 2z = 5
(taken from here) you can substitute the equations for x, y, and z above into the equation for the plane, which is now in only the parameter t. Solve for t. This is the particular value of t for that line that lies in the plane. Then you can solve for x, y, and z by going back up to the line equations and substituting t back in.