有没有一种简单的方法来确定一个点是否在三角形内?是2D的,不是3D的。


当前回答

我在JavaScript中改编的高性能代码(文章如下):

function pointInTriangle (p, p0, p1, p2) {
  return (((p1.y - p0.y) * (p.x - p0.x) - (p1.x - p0.x) * (p.y - p0.y)) | ((p2.y - p1.y) * (p.x - p1.x) - (p2.x - p1.x) * (p.y - p1.y)) | ((p0.y - p2.y) * (p.x - p2.x) - (p0.x - p2.x) * (p.y - p2.y))) >= 0;
}

pointInTriangle(p, p0, p1, p2) -用于逆时针方向的三角形 pointInTriangle(p, p0, p1, p2) -用于顺时针三角形

在jsFiddle(包括性能测试)中,在一个单独的函数中也有缠绕检查。或按下面的“运行代码片段”

var ctx = $("canvas")[0].getContext("2d"); var W = 500; var H = 500; var point = { x: W / 2, y: H / 2 }; var triangle = randomTriangle(); $("canvas").click(function(evt) { point.x = evt.pageX - $(this).offset().left; point.y = evt.pageY - $(this).offset().top; test(); }); $("canvas").dblclick(function(evt) { triangle = randomTriangle(); test(); }); document.querySelector('#performance').addEventListener('click', _testPerformance); test(); function test() { var result = checkClockwise(triangle.a, triangle.b, triangle.c) ? pointInTriangle(point, triangle.a, triangle.c, triangle.b) : pointInTriangle(point, triangle.a, triangle.b, triangle.c); var info = "point = (" + point.x + "," + point.y + ")\n"; info += "triangle.a = (" + triangle.a.x + "," + triangle.a.y + ")\n"; info += "triangle.b = (" + triangle.b.x + "," + triangle.b.y + ")\n"; info += "triangle.c = (" + triangle.c.x + "," + triangle.c.y + ")\n"; info += "result = " + (result ? "true" : "false"); $("#result").text(info); render(); } function _testPerformance () { var px = [], py = [], p0x = [], p0y = [], p1x = [], p1y = [], p2x = [], p2y = [], p = [], p0 = [], p1 = [], p2 = []; for(var i = 0; i < 1000000; i++) { p[i] = {x: Math.random() * 100, y: Math.random() * 100}; p0[i] = {x: Math.random() * 100, y: Math.random() * 100}; p1[i] = {x: Math.random() * 100, y: Math.random() * 100}; p2[i] = {x: Math.random() * 100, y: Math.random() * 100}; } console.time('optimal: pointInTriangle'); for(var i = 0; i < 1000000; i++) { pointInTriangle(p[i], p0[i], p1[i], p2[i]); } console.timeEnd('optimal: pointInTriangle'); console.time('original: ptInTriangle'); for(var i = 0; i < 1000000; i++) { ptInTriangle(p[i], p0[i], p1[i], p2[i]); } console.timeEnd('original: ptInTriangle'); } function pointInTriangle (p, p0, p1, p2) { return (((p1.y - p0.y) * (p.x - p0.x) - (p1.x - p0.x) * (p.y - p0.y)) | ((p2.y - p1.y) * (p.x - p1.x) - (p2.x - p1.x) * (p.y - p1.y)) | ((p0.y - p2.y) * (p.x - p2.x) - (p0.x - p2.x) * (p.y - p2.y))) >= 0; } function ptInTriangle(p, p0, p1, p2) { var s = (p0.y * p2.x - p0.x * p2.y + (p2.y - p0.y) * p.x + (p0.x - p2.x) * p.y); var t = (p0.x * p1.y - p0.y * p1.x + (p0.y - p1.y) * p.x + (p1.x - p0.x) * p.y); if (s <= 0 || t <= 0) return false; var A = (-p1.y * p2.x + p0.y * (-p1.x + p2.x) + p0.x * (p1.y - p2.y) + p1.x * p2.y); return (s + t) < A; } function render() { ctx.fillStyle = "#CCC"; ctx.fillRect(0, 0, 500, 500); drawTriangle(triangle.a, triangle.b, triangle.c); drawPoint(point); } function checkClockwise(p0, p1, p2) { var A = (-p1.y * p2.x + p0.y * (-p1.x + p2.x) + p0.x * (p1.y - p2.y) + p1.x * p2.y); return A > 0; } function drawTriangle(p0, p1, p2) { ctx.fillStyle = "#999"; ctx.beginPath(); ctx.moveTo(p0.x, p0.y); ctx.lineTo(p1.x, p1.y); ctx.lineTo(p2.x, p2.y); ctx.closePath(); ctx.fill(); ctx.fillStyle = "#000"; ctx.font = "12px monospace"; ctx.fillText("1", p0.x, p0.y); ctx.fillText("2", p1.x, p1.y); ctx.fillText("3", p2.x, p2.y); } function drawPoint(p) { ctx.fillStyle = "#F00"; ctx.beginPath(); ctx.arc(p.x, p.y, 5, 0, 2 * Math.PI); ctx.fill(); } function rand(min, max) { return Math.floor(Math.random() * (max - min + 1)) + min; } function randomTriangle() { return { a: { x: rand(0, W), y: rand(0, H) }, b: { x: rand(0, W), y: rand(0, H) }, c: { x: rand(0, W), y: rand(0, H) } }; } <script src="https://cdnjs.cloudflare.com/ajax/libs/jquery/1.9.1/jquery.min.js"></script> <button id="performance">Run performance test (open console)</button> <pre>Click: place the point. Double click: random triangle.</pre> <pre id="result"></pre> <canvas width="500" height="500"></canvas>

受此启发: http://www.phatcode.net/articles.php?id=459

其他回答

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]。

下面是一个高效的Python实现:

def PointInsideTriangle2(pt,tri):
    '''checks if point pt(2) is inside triangle tri(3x2). @Developer'''
    a = 1/(-tri[1,1]*tri[2,0]+tri[0,1]*(-tri[1,0]+tri[2,0])+ \
        tri[0,0]*(tri[1,1]-tri[2,1])+tri[1,0]*tri[2,1])
    s = a*(tri[2,0]*tri[0,1]-tri[0,0]*tri[2,1]+(tri[2,1]-tri[0,1])*pt[0]+ \
        (tri[0,0]-tri[2,0])*pt[1])
    if s<0: return False
    else: t = a*(tri[0,0]*tri[1,1]-tri[1,0]*tri[0,1]+(tri[0,1]-tri[1,1])*pt[0]+ \
              (tri[1,0]-tri[0,0])*pt[1])
    return ((t>0) and (1-s-t>0))

和一个示例输出:

bool isInside( float x, float y, float x1, float y1, float x2, float y2, float x3, float y3 ) {
  float l1 = (x-x1)*(y3-y1) - (x3-x1)*(y-y1), 
    l2 = (x-x2)*(y1-y2) - (x1-x2)*(y-y2), 
    l3 = (x-x3)*(y2-y3) - (x2-x3)*(y-y3);
  return (l1>0 && l2>0  && l3>0) || (l1<0 && l2<0 && l3<0);
}

没有比这更有效率的了!三角形的每边都可以有独立的位置和方向,因此需要进行l1、l2和l3三个计算,每个计算需要进行2次乘法。一旦l1, l2和l3是已知的,结果只是一些基本的比较和布尔运算。

由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,而不管缠绕方向如何。