使用 Newell 方法在 Python 中计算表面法线

Question

我正在尝试根据来自here 的以下伪代码，在 Python 中实现 Newell 方法来计算表面法线向量。

Begin Function CalculateSurfaceNormal (Input Polygon) Returns Vector

   Set Vertex Normal to (0, 0, 0)

   Begin Cycle for Index in [0, Polygon.vertexNumber)

      Set Vertex Current to Polygon.verts[Index]
      Set Vertex Next    to Polygon.verts[(Index plus 1) mod Polygon.vertexNumber]

      Set Normal.x to Sum of Normal.x and (multiply (Current.y minus Next.y) by (Current.z plus Next.z))
      Set Normal.y to Sum of Normal.y and (multiply (Current.z minus Next.z) by (Current.x plus Next.x))
      Set Normal.z to Sum of Normal.z and (multiply (Current.x minus Next.x) by (Current.y plus Next.y))

   End Cycle

   Returning Normalize(Normal)

End Function

这是我的代码：

Point3D = collections.namedtuple('Point3D', 'x y z')

def surface_normal(poly):
    n = [0.0, 0.0, 0.0]

    for i, v_curr in enumerate(poly):
        v_next = poly[(i+1) % len(poly)]
        n[0] += (v_curr.y - v_next.y) * (v_curr.z - v_next.z)
        n[1] += (v_curr.z - v_next.z) * (v_curr.x - v_next.x)
        n[2] += (v_curr.x - v_next.x) * (v_curr.y - v_next.y)

    normalised = [i/sum(n) for i in n]

    return normalised

def test_surface_normal():
    poly = [Point3D(0.0, 0.0, 0.0),
            Point3D(0.0, 1.0, 0.0),
            Point3D(1.0, 1.0, 0.0),
            Point3D(1.0, 0.0, 0.0)]

    assert surface_normal(poly) == [0.0, 0.0, 1.0]

这在标准化步骤中失败，因为此时的n 是[0.0, 0.0, 0.0]。如果我理解正确，应该是[0.0, 0.0, 1.0]（Wolfram Alpha 的confirmed）。

我在这里做错了什么？有没有更好的方法在 python 中计算表面法线？我的多边形将始终是平面的，因此如果有其他方法，纽厄尔方法并不是绝对必要的。

Answer 1

好吧，这个问题实际上是一个愚蠢的问题。

像这样的行：

n[0] += (v_curr.y - v_next.y) * (v_curr.z - v_next.z)

应该是：

n[0] += (v_curr.y - v_next.y) * (v_curr.z + v_next.z) 

第二组括号中的值应该相加，而不是相减。

Answer 2

如果您想要 Newell 方法的替代方法，您可以使用 2 个非平行向量的叉积。这应该适用于您提供的任何平面形状。我知道理论说它适用于凸多边形，但我们在 Wolfram Alpha 上看到的示例会为偶数凹多边形返回适当的表面法线（例如 bowtie 多边形）。

Answer 3

我认为有必要澄清一些关于上面 Newell 方法的信息实现的几点。在第一次阅读https://stackoverflow.com/users/1706564/jamie-bull 更新时，我认为符号的变化仅涉及第一个等式。但事实上它们都是。其次，提供一个比较叉积的方法。最后，处理除以零的问题 - 但可以只除以 epsilon。我也使用 numpy 代替。

import numpy as np

def surface_normal_newell(poly):

    n = np.array([0.0, 0.0, 0.0])

    for i, v_curr in enumerate(poly):
        v_next = poly[(i+1) % len(poly),:]
        n[0] += (v_curr[1] - v_next[1]) * (v_curr[2] + v_next[2]) 
        n[1] += (v_curr[2] - v_next[2]) * (v_curr[0] + v_next[0])
        n[2] += (v_curr[0] - v_next[0]) * (v_curr[1] + v_next[1])

    norm = np.linalg.norm(n)
    if norm==0:
        raise ValueError('zero norm')
    else:
        normalised = n/norm

    return normalised


def surface_normal_cross(poly):

    n = np.cross(poly[1,:]-poly[0,:],poly[2,:]-poly[0,:])

    norm = np.linalg.norm(n)
    if norm==0:
        raise ValueError('zero norm')
    else:
        normalised = n/norm

    return normalised


def test_surface_normal3():
    """ should return:

        Newell:
        Traceback (most recent call last):
          File "demo_newells_surface_normals.py", line 96, in <module>
            test_surface_normal3()
          File "demo_newells_surface_normals.py", line 58, in test_surface_normal3
            print "Newell:", surface_normal_newell(poly) 
          File "demo_newells_surface_normals.py", line 24, in surface_normal_newell
            raise ValueError('zero norm')
        ValueError: zero norm
    """
    poly = np.array([[1.0,0.0,0.0],
                     [1.0,0.0,0.0],
                     [1.0,0.0,0.0]])
    print "Newell:", surface_normal_newell(poly) 


def test_surface_normal2():
    """ should return:

        Newell: [ 0.08466675 -0.97366764 -0.21166688]
        Cross : [ 0.08466675 -0.97366764 -0.21166688]
    """
    poly = np.array([[6.0,1.0,4.0],
                     [7.0,0.0,9.0],
                     [1.0,1.0,2.0]])
    print "Newell:", surface_normal_newell(poly)
    print "Cross :", surface_normal_cross(poly)


def test_surface_normal1():
    """ should return:

        Newell: [ 0.  0. -1.]
        Cross : [ 0.  0. -1.]
    """
    poly = np.array([[0.0,1.0,0.0],
                     [1.0,1.0,0.0],
                     [1.0,0.0,0.0]])
    print "Newell:", surface_normal_newell(poly) 
    print "Cross :", surface_normal_cross(poly)


print "Test 1:"
test_surface_normal1()
print "\n"

print "Test 2:"
test_surface_normal2()
print "\n"

print "Test 3:"
test_surface_normal3()