Computing polygon normal with Newell's method

[solub]

Hi guys,

I’m having difficulties computing the normal of a polygon (more than 3 vertices) based on the Newell’s method.

                 𝑛𝑥=∑𝑖=0𝑛(𝑦𝑖−𝑦𝑖+1)(𝑧𝑖+𝑧𝑖+1)

                 𝑛𝑦=∑𝑖=0𝑛(𝑧𝑖−𝑧𝑖+1)(𝑥𝑖+𝑥𝑖+1)

                 𝑛𝑧=∑𝑖=0𝑛(𝑥𝑖−𝑥𝑖+1)(𝑦𝑖+𝑦𝑖+1)

It seems pretty straightforward and I’ve seen various implementations (in Java (line 123), in C++ or in Python) and yet mine gives me very strange results:

def computePolyNormal(face):
    
    #Calculate normal with Newell's method
    nml = PVector()
    
    for i in range(face.getVertexCount() - 1):
        
        p0 = face.getVertex(i) #current vertex
        p1 = face.getVertex(i+1) #next vertex
        
        nml.x += (p0.y - p1.y) * (p0.z + p1.z)
        nml.y += (p0.z - p1.z) * (p0.x + p1.x)
        nml.z += (p0.x - p1.x) * (p0.y + p1.y)
        
    nml.normalize()

(Normals pointing towards strange directions)

Can someone help me find what’s wrong with this code ?

Please note that for the sake of clarity and simplicity I’m showing triangle faces but ideally the Newell’s formula should work for more complex polygons.

add_library('peasycam')


def setup():
    global tetrahedron, centroid
    hint(DISABLE_DEPTH_TEST)
    size(1000, 600, P3D)
    fill(55, 180, 220)
    smooth(8)
    
    cam = PeasyCam(this, 100)
    

    tetrahedron = createShape(GROUP)
    
    #First face
    f1 = createShape()
    f1.beginShape(TRIANGLES)
    f1.vertex(10, 10, 0)
    f1.vertex(40, 30, 0)
    f1.vertex(2, 35, 0)
    f1.endShape()
      
    #Second face  
    f2 = createShape()
    f2.beginShape(TRIANGLES)
    f2.vertex(10, 10, 0)
    f2.vertex(40, 30, 0)
    f2.vertex(17, 25, 40)
    f2.endShape()
    
    #Third face  
    f3 = createShape()
    f3.beginShape(TRIANGLES)
    f3.vertex(40, 30, 0)
    f3.vertex(2, 35, 0)
    f3.vertex(17, 25, 40)
    f3.endShape()
    
    #Fourth face  
    f4 = createShape()
    f4.beginShape(TRIANGLES)
    f4.vertex(10, 10, 0)
    f4.vertex(2, 35, 0)
    f4.vertex(17, 25, 40)
    f4.endShape()
    
    tetrahedron.addChild(f1)
    tetrahedron.addChild(f2)
    tetrahedron.addChild(f3)
    tetrahedron.addChild(f4)
    
    
    #Compute centroid
    centroid = PVector()
    count = 0
    
    for i in range(tetrahedron.getChildCount()):
        face = tetrahedron.getChild(i)
        for v in range(face.getVertexCount()):
            vert = face.getVertex(v)
            
            centroid.x += vert.x
            centroid.y += vert.y
            centroid.z += vert.z

            count += 1
            
    centroid.div(count)
    
    
    
def draw():
    background(255)
    strokeWeight(2)
    
    #Rotate
    rotateX(HALF_PI)
    rotate(frameCount*.004)
    
    #Tanslate to the center of the screen
    pushMatrix()
    translate(-centroid.x, -centroid.y, -centroid.z)
    
    #Display Tetrahedron
    shape(tetrahedron)
    
    
    #Display normals
    for i in range(tetrahedron.getChildCount()):
        face = tetrahedron.getChild(i)
        computePolyNormal(face)
        #computeTriangleNormal(face) #Work but only for triangle faces

    popMatrix()
    

    
    
def computePolyNormal(face):
    
    #Calculate normal with Newell's method
    nml = PVector()
    
    for i in range(face.getVertexCount() - 1):
        
        p0 = face.getVertex(i) #current vertex
        p1 = face.getVertex(i+1) #next vertex
        
        nml.x += (p0.y - p1.y) * (p0.z + p1.z)
        nml.y += (p0.z - p1.z) * (p0.x + p1.x)
        nml.z += (p0.x - p1.x) * (p0.y + p1.y)
        
    #Normalizing then scaling-up to see a line
    nml.normalize().mult(10)
    
    
    #Center of the polygon
    center = PVector()
    
    for i in range(face.getVertexCount()):
        center.x += face.getVertex(i).x
        center.y += face.getVertex(i).y
        center.z += face.getVertex(i).z
            
    center.div(face.getVertexCount())
    
        
    #line
    stroke(30, 255, 100)
    line(center.x, center.y, center.z, center.x + nml.x, center.y + nml.y, center.z + nml.z)
    
        

def computeTriangleNormal(face):
    
    #Vertices of the face
    v1 = face.getVertex(0)
    v2 = face.getVertex(1)
    v3 = face.getVertex(2)
    
    #Center of the face
    center = PVector((v1.x + v2.x + v3.x) / 3, (v1.y + v2.y + v3.y) / 3, (v1.z + v2.z + v3.z) / 3)
    
    #Normal
    a = v2.sub(v1)
    b = v3.sub(v1)
    nml = a.cross(b).normalize().mult(10)
    
    #Check if normal is heading inwards or outwards
    d = PVector.dot(center - centroid, nml)
    
    if d < 0:
        #If inwards -> reverse
        nml.mult(-1)
        stroke(30, 255, 100)
    else:
        stroke(255, 50, 100)
        
    #line
    line(center.x, center.y, center.z, center.x + nml.x, center.y + nml.y, center.z + nml.z)

[villares]

Hi @solub!

It’s beyond me why the Newell’s method is not working
But, if you have a planar polygon, pick any 3 points and use the triangle method that is working to get the angle!

[solub]

Thanks for the reply @villares,

I just realized the for loop needs to go full circle (p0-p1, p1-p2, p2-p0)
It seems to work now.

def computePolyNormal(face):
    
    #Calculate normal with Newell's method
    nml = PVector()
    
    for i in range(face.getVertexCount()):
                
        p0 = face.getVertex(i) #current vertex
        p1 = face.getVertex((i+1)%face.getVertexCount()) #next vertex
        
        nml.x += (p0.y - p1.y) * (p0.z + p1.z)
        nml.y += (p0.z - p1.z) * (p0.x + p1.x)
        nml.z += (p0.x - p1.x) * (p0.y + p1.y)

        
    #Normalizing then scaling-up to see a line
    nml.normalize().mult(10)