有没有一种简单的方法来确定一个点是否在三角形内?是2D的,不是3D的。
当前回答
这是确定一个点是在三角形的内、外还是在三角形的臂上的最简单的概念。
用行列式确定三角形内的点:
最简单的工作代码:
#-*- coding: utf-8 -*-
import numpy as np
tri_points = [(1,1),(2,3),(3,1)]
def pisinTri(point,tri_points):
Dx , Dy = point
A,B,C = tri_points
Ax, Ay = A
Bx, By = B
Cx, Cy = C
M1 = np.array([ [Dx - Bx, Dy - By, 0],
[Ax - Bx, Ay - By, 0],
[1 , 1 , 1]
])
M2 = np.array([ [Dx - Ax, Dy - Ay, 0],
[Cx - Ax, Cy - Ay, 0],
[1 , 1 , 1]
])
M3 = np.array([ [Dx - Cx, Dy - Cy, 0],
[Bx - Cx, By - Cy, 0],
[1 , 1 , 1]
])
M1 = np.linalg.det(M1)
M2 = np.linalg.det(M2)
M3 = np.linalg.det(M3)
print(M1,M2,M3)
if(M1 == 0 or M2 == 0 or M3 ==0):
print("Point: ",point," lies on the arms of Triangle")
elif((M1 > 0 and M2 > 0 and M3 > 0)or(M1 < 0 and M2 < 0 and M3 < 0)):
#if products is non 0 check if all of their sign is same
print("Point: ",point," lies inside the Triangle")
else:
print("Point: ",point," lies outside the Triangle")
print("Vertices of Triangle: ",tri_points)
points = [(0,0),(1,1),(2,3),(3,1),(2,2),(4,4),(1,0),(0,4)]
for c in points:
pisinTri(c,tri_points)
其他回答
由andreasdr和Perro Azul发布的重心方法的c#版本。我添加了一个检查,当s和t有相反的符号(而且都不为零)时,放弃面积计算,因为潜在地避免三分之一的乘法成本似乎是合理的。
public static bool PointInTriangle(Point p, Point p0, Point p1, Point p2)
{
var s = (p0.X - p2.X) * (p.Y - p2.Y) - (p0.Y - p2.Y) * (p.X - p2.X);
var t = (p1.X - p0.X) * (p.Y - p0.Y) - (p1.Y - p0.Y) * (p.X - p0.X);
if ((s < 0) != (t < 0) && s != 0 && t != 0)
return false;
var d = (p2.X - p1.X) * (p.Y - p1.Y) - (p2.Y - p1.Y) * (p.X - p1.X);
return d == 0 || (d < 0) == (s + t <= 0);
}
2021年更新:这个版本正确处理任意一个缠绕方向(顺时针和逆时针)指定的三角形。请注意,对于恰好位于三角形边缘上的点,本页上的一些其他答案会给出不一致的结果,这取决于三角形三个点的排列顺序。这些点被认为是“在”三角形中,这段代码正确地返回true,而不管缠绕方向如何。
一个简单的方法是:
找出连接 分别指向三角形的三个点 顶点和夹角之和 这些向量。如果它们的和 角度是2*那么点是 在三角形里面。
两个解释替代方案的好网站是:
黑卒和沃尔夫勒姆
python中的其他函数,比Developer的方法更快(至少对我来说),并受到Cédric Dufour解决方案的启发:
def ptInTriang(p_test, p0, p1, p2):
dX = p_test[0] - p0[0]
dY = p_test[1] - p0[1]
dX20 = p2[0] - p0[0]
dY20 = p2[1] - p0[1]
dX10 = p1[0] - p0[0]
dY10 = p1[1] - p0[1]
s_p = (dY20*dX) - (dX20*dY)
t_p = (dX10*dY) - (dY10*dX)
D = (dX10*dY20) - (dY10*dX20)
if D > 0:
return ( (s_p >= 0) and (t_p >= 0) and (s_p + t_p) <= D )
else:
return ( (s_p <= 0) and (t_p <= 0) and (s_p + t_p) >= D )
你可以用:
X_size = 64
Y_size = 64
ax_x = np.arange(X_size).astype(np.float32)
ax_y = np.arange(Y_size).astype(np.float32)
coords=np.meshgrid(ax_x,ax_y)
points_unif = (coords[0].reshape(X_size*Y_size,),coords[1].reshape(X_size*Y_size,))
p_test = np.array([0 , 0])
p0 = np.array([22 , 8])
p1 = np.array([12 , 55])
p2 = np.array([7 , 19])
fig = plt.figure(dpi=300)
for i in range(0,X_size*Y_size):
p_test[0] = points_unif[0][i]
p_test[1] = points_unif[1][i]
if ptInTriang(p_test, p0, p1, p2):
plt.plot(p_test[0], p_test[1], '.g')
else:
plt.plot(p_test[0], p_test[1], '.r')
绘制网格需要花费很多时间,但是该网格在0.0195319652557秒内测试,而开发人员代码为0.0844349861145秒。
最后是代码注释:
# Using barycentric coordintes, any point inside can be described as:
# X = p0.x * r + p1.x * s + p2.x * t
# Y = p0.y * r + p1.y * s + p2.y * t
# with:
# r + s + t = 1 and 0 < r,s,t < 1
# then: r = 1 - s - t
# and then:
# X = p0.x * (1 - s - t) + p1.x * s + p2.x * t
# Y = p0.y * (1 - s - t) + p1.y * s + p2.y * t
#
# X = p0.x + (p1.x-p0.x) * s + (p2.x-p0.x) * t
# Y = p0.y + (p1.y-p0.y) * s + (p2.y-p0.y) * t
#
# X - p0.x = (p1.x-p0.x) * s + (p2.x-p0.x) * t
# Y - p0.y = (p1.y-p0.y) * s + (p2.y-p0.y) * t
#
# we have to solve:
#
# [ X - p0.x ] = [(p1.x-p0.x) (p2.x-p0.x)] * [ s ]
# [ Y - p0.Y ] [(p1.y-p0.y) (p2.y-p0.y)] [ t ]
#
# ---> b = A*x ; ---> x = A^-1 * b
#
# [ s ] = A^-1 * [ X - p0.x ]
# [ t ] [ Y - p0.Y ]
#
# A^-1 = 1/D * adj(A)
#
# The adjugate of A:
#
# adj(A) = [(p2.y-p0.y) -(p2.x-p0.x)]
# [-(p1.y-p0.y) (p1.x-p0.x)]
#
# The determinant of A:
#
# D = (p1.x-p0.x)*(p2.y-p0.y) - (p1.y-p0.y)*(p2.x-p0.x)
#
# Then:
#
# s_p = { (p2.y-p0.y)*(X - p0.x) - (p2.x-p0.x)*(Y - p0.Y) }
# t_p = { (p1.x-p0.x)*(Y - p0.Y) - (p1.y-p0.y)*(X - p0.x) }
#
# s = s_p / D
# t = t_p / D
#
# Recovering r:
#
# r = 1 - (s_p + t_p)/D
#
# Since we only want to know if it is insidem not the barycentric coordinate:
#
# 0 < 1 - (s_p + t_p)/D < 1
# 0 < (s_p + t_p)/D < 1
# 0 < (s_p + t_p) < D
#
# The condition is:
# if D > 0:
# s_p > 0 and t_p > 0 and (s_p + t_p) < D
# else:
# s_p < 0 and t_p < 0 and (s_p + t_p) > D
#
# s_p = { dY20*dX - dX20*dY }
# t_p = { dX10*dY - dY10*dX }
# D = dX10*dY20 - dY10*dX20
我只是想用一些简单的向量数学来解释安德里亚斯给出的重心坐标解,它会更容易理解。
区域A定义为s * v02 + t * v01给出的任意向量,条件s >= 0, t >= 0。如果三角形v0 v1 v2内的任意一点,它一定在区域A内。
如果进一步限制s, t属于[0,1]。得到包含s * v02 + t * v01的所有向量的区域B,条件s, t属于[0,1]。值得注意的是,区域B的下部是三角形v0, v1, v2的镜像。问题来了,我们是否可以给定一定的s和t条件,来进一步排除区域B的低部分。
假设我们给出一个值s, t在[0,1]内变化。在下图中,点p位于v1v2的边缘。s * v02 + t * v01的所有向量沿着虚线通过简单向量和得到。在v1v2和虚线交点p处,我们有:
(1-S)|V0v2|/ |v0v2|= tp|v0v1|/ |v0v1|
得到1 - s = tp,然后1 = s + tp。如果任意t > tp,即1 < s + t where在双虚线上,则该向量在三角形外,任意t <= tp,即1 >= s + t where在单虚线上,则该向量在三角形内。
如果我们给出[0,1]中的任意s,对应的t必须满足1 >= s + t,对于三角形内的向量。
最后我们得到v = s * v02 +t * v01, v在三角形内,条件s, t, s+t属于[0,1]。然后翻译到点,我们有
P - p0 = s * (p1 - p0) + t * (p2 - p0), and s, t, s + t in [0,1]
和Andreas解方程组的解是一样的 P = p0 + s * (p1 - p0) + t * (p2 - p0),带s, t, s + t属于[0,1]。
其中一个最简单的方法来检查是否由三角形的顶点组成的面积 (x1,y1) (x2,y2) (x3,y3)是否为正。
面积可由公式计算:
1/2 [x1(y2–y3) + x2(y3–y1) + x3(y1–y2)]
或者python代码可以写成:
def triangleornot(p1,p2,p3):
return (1/ 2) [p1[0](p2[1]–p3[1]) + p2[0] (p3[1]–p1[1]) + p3[0] (p1[0]–p2[0])]