3D Ligne d'Intersection des plans

Heres un exemple Python qui trouve à l'intersection d'une droite et d'un plan.

Où l'avion peut être soit un point et un normal, ou un vecteur 4d (forme normale),
Dans les exemples ci-dessous (code pour les deux sont fournis).

Également noter que cette fonction calcule une valeur représentant l'endroit où le point est sur la ligne, (appelé fac dans le code ci-dessous). Vous pouvez retourner ce trop, parce que les valeurs de 0 à 1 en intersection avec le segment de ligne - ce qui peut être utile pour l'appelant.

D'autres détails à noter dans le code-commentaires.

Remarque: Cet exemple utilise les fonctions pures, sans dépendances - pour le rendre facile à se déplacer vers d'autres langues. Avec un Vector type de données et la surcharge d'opérateur, il peut être plus concis (inclus dans l'exemple ci-dessous).

# 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)
else:
# 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,
)

Si l'avion est défini comme un vecteur 4d (forme normale), nous avons besoin de trouver un point sur le plan, puis de calculer l'intersection avant (voir p_co affectation).

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)
else:
return None

Pour référence ultérieure, ceci a été pris à partir de Blender et adapté à Python.
isect_line_plane_v3() dans math_geom.c

Pour plus de clarté, voici les versions utilisant le mathutils API (qui peut être modifié pour d'autres bibliothèques de mathématiques avec la surcharge d'opérateur).

# 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)
else:
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)
else:
return None

Utilisation de numpy et python:

#Based on http://geomalgorithms.com/a05-_intersect-1.html
from __future__ import print_function
import numpy as np
epsilon=1e-6
#Define plane
planeNormal = np.array([0, 0, 1])
planePoint = np.array([0, 0, 5]) #Any point on the plane
#Define ray
rayDirection = np.array([0, -1, -1])
rayPoint = np.array([0, 0, 10]) #Any point along the ray
ndotu = planeNormal.dot(rayDirection) 
if abs(ndotu) < epsilon:
print ("no intersection or line is within plane")
w = rayPoint - planePoint
si = -planeNormal.dot(w) /ndotu
Psi = w + si * rayDirection + planePoint
print ("intersection at", Psi)

Si vous avez deux points p et q qui définissent une ligne et un avion dans le grand cartésienne de la forme ax+by+cz+d = 0, vous pouvez utiliser la méthode paramétrique.

Si vous avez besoin de ce aux fins de codage, voici un extrait de code javascript:

JS:

/**
* findLinePlaneIntersectionCoords (to avoid requiring unnecessary instantiation)
* Given points p with px py pz and q that define a line, and the plane
* of formula ax+by+cz+d = 0, returns the intersection point or null if none.
*/
function findLinePlaneIntersectionCoords(px, py, pz, qx, qy, qz, a, b, c, d) {
var tDenom = a*(qx-px) + b*(qy-py) + c*(qz-pz);
if (tDenom == 0) return null;
var t = - ( a*px + b*py + c*pz + d ) / tDenom;
return {
x: (px+t*(qx-px)),
y: (py+t*(qy-py)),
z: (pz+t*(qz-pz))
};
}
//Example (plane at y = 10  and perpendicular line from the origin)
console.log(JSON.stringify(findLinePlaneIntersectionCoords(0,0,0,0,1,0,0,1,0,-10)));
//Example (no intersection, plane and line are parallel)
console.log(JSON.stringify(findLinePlaneIntersectionCoords(0,0,0,0,0,1,0,1,0,-10)));

Juste pour élargir sur ZGorlock de l' réponse, j'ai fait le produit scalaire, plus et scalse de la 3D Vecteurs. Les références de ces calculs sont Produit Scalaire, Ajouter deux vecteurs 3D et Mise à l'échelle. Remarque: Vec3D est juste une classe personnalisée qui a points: x, y et z.

/**
* 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 Vec3D lineIntersection(Vec3D planePoint, Vec3D planeNormal, Vec3D linePoint, Vec3D lineDirection) {
//ax × bx + ay × by
int dot = (int) (planeNormal.x * lineDirection.x + planeNormal.y * lineDirection.y);
if (dot == 0) {
return null;
}
//Ref for dot product calculation: https://www.mathsisfun.com/algebra/vectors-dot-product.html
int dot2 = (int) (planeNormal.x * planePoint.x + planeNormal.y * planePoint.y);
int dot3 = (int) (planeNormal.x * linePoint.x + planeNormal.y * linePoint.y);
int dot4 = (int) (planeNormal.x * lineDirection.x + planeNormal.y * lineDirection.y);
double t = (dot2 - dot3) /dot4;
float xs = (float) (lineDirection.x * t);
float ys = (float) (lineDirection.y * t);
float zs = (float) (lineDirection.z * t);
Vec3D lineDirectionScale = new Vec3D( xs, ys, zs);
float xa = (linePoint.x + lineDirectionScale.x);
float ya = (linePoint.y + lineDirectionScale.y);
float za = (linePoint.z + lineDirectionScale.z);
return new Vec3D(xa, ya, za);
}