How to efficiently determine the normal to a polygon in 3D space?

Question

I have a bunch of coplanar points defining a polygon in 3D space. These are always wound the same way (e.g. clockwise). I need to determine the signed normal to the plane containing this polygon, i.e., know which way is "up" for that polygon.

This seems easy enough at first: take two edges (vertex differences) and compute the cross product. But that fails if the edges happen to be colinear (you get a zero-magnitude cross product).

So then I tried walking the vertex list until I find a second edge that makes a reasonably large angle with the first edge. This works reliably on a convex polygon, but it can fail (point the opposite direction) on a non-convex polygon if the two edges I end up with don't define a triangle that's inside the polygon.

I know that if I triangulated the polygon first, then I could easily and reliably check the facing of any triangle... but the problem is that my triangulation library requires knowing the plane normal. So, egg must come before chicken.

How do I pick two edges (or three vertices) in an non-convex polygon that will reliably define which way the polygon is facing?

Hazelwood

Offtopic. Not a programming question. This is more geometry/math. Try math.stackexchange.com – Colonial

Asked here because mathematicians often give an answer that is not easy (or not clear how) to implement, e.g., "a vertex is an ear if a line connecting the neighboring points is entirely within the polygon"... a true statement, but not helpful for actually coding it. I'm asking for an algorithm that I can code. – Hazelwood

@JoeStrout I know you already have a satisfactory answer, so just wanted to draw your attention that I added a different (simpler?) approach and I wasn't sure that you'll "come looking" for another answer. – Nayarit

Answer 1

If I were you, I would have done it in the following way:

1. Choose any point C near the polygon (any vertex or mass center).
2. Sum cross products (P[i] - C) x (P[i+1] - C) for all i (including last and first points pair).
3. Normalize the sum vector.

Note that after step 2 you have a vector which has normal direction with proper orientation, and its magnitude is 2 S, where S is the area of your polygon. That's why it should work unless your polygon has zero or almost zero area.

By the way, point C is used here only to make calculation a bit more precise for small polygons located far from the origin. You can choose C = 0, efficiently removing it from calculations.

Decahedron

Tested and confirmed to work (in all my test cases, at least). Thanks very much! – Hazelwood

@YvesDaoust: I'm note sure that you read the question completely. Your solution does not produce correct orientation of normal. Moreover, it can be unstable (depending on implementation). – Decahedron

Answer 2

The sum of consecutive cross products provided by stgatilov is not robust. See Robust polygon normal calculation for example.

A robust solution is to find the largest cross product (P[i] - C) x (P[j] - C) for all i, j, (i < j) and normalize it. It will correspond to the largest inscribed triangle of the polygon.

Sheetfed

couldn't this end up pointing the wrong way for concave polygons? What's the choice of C? – Abaxial

Answer 3

1. Get 3 points (p1, p2, p3) from your polygon which aren't lying on a same line (you can use dot product check to avoid points which lie between colinear segments).

2. Make a plane from these points (defined by normal vector and origin distance: Ax + By + Cz + D = 0). Plane equation from 3 points in space is well known formula. For arbitrary polygon this plane might have normal which points in the opposite direction, but this doesn't matter at this step.

3. Project any segment from your polygon (two connected consecutive points) onto this plane. Now you have 2 points in plane's 2D space. This is also well known op - you just need to get plane origin (which is n * D) and x and y plane basis vectors which you can build out of plane normal and project segment points vectors onto basis vectors to get their 2D x and y coords in plane's space.

4. Take a 2D normal of this segment, in 2D space it's very easy - just rotate 2D vector by 90 (or -90, depending on your coord. system) - swap (x, y) and negate one of the coords.

5. Transform new rotated 2D vector back to 3D space - just use plane origin point and x, y basis vectors you already have from step 3. Let's say this is your binormal vector.

6. Now you have your original segment direction vector (just subtract segment points, let's say it's tangent) and your binormal vector from step 5 - cross product them and you'll get the arbitrary polygon normal which will always point in correct direction no matter which segment you use to produce it because your stated that points always wound the same direction.

Pericarditis