Projecting to a 2D Plane for Barycentric Calculations

Question

I have three vertices which make up a plane/polygon in 3D Space, v0, v1 & v2.

To calculate barycentric co-ordinates for a 3D point upon this plane I must first project both the plane and point into 2D space.

After trawling the web I have a good understanding of how to calculate barycentric co-ordinates in 2D space, but I am stuck at finding the best way to project my 3D points into a suitable 2D plane.

It was suggested to me that the best way to achieve this was to "drop the axis with the smallest projection". Without testing the area of the polygon formed when projected on each world axis (xy, yz, xz) how can I determine which projection is best (has the largest area), and therefore is most suitable for calculating the most accurate barycentric co-ordinate?

Answer 1

Example of computation of barycentric coordinates in 3D space as requested by the OP. Given:

• 3D points v0, v1, v2 that define the triangle
• 3D point p that lies on the plane defined by v0, v1 and v2 and inside the triangle spanned by the same points.

"x" denotes the cross product between two 3D vectors.
"len" denotes the length of a 3D vector.
"u", "v", "w" are the barycentric coordinates belonging to v0, v1 and v2 respectively.

triArea =   len((v1 - v0) x (v2 - v0)) * 0.5
u =       ( len((v1 - p ) x (v2 - p )) * 0.5 ) / triArea
v =       ( len((v0 - p ) x (v2 - p )) * 0.5 ) / triArea
w =       ( len((v0 - p ) x (v1 - p )) * 0.5 ) / triArea

=> p == u * v0 + v * v1 + w * v2

The cross product is defined like this:

v0 x v1 := { v0.y * v1.z - v0.z * v1.y,
             v0.z * v1.x - v0.x * v1.z,
             v0.x * v1.y - v0.y * v1.x }

Answer 2

WARNING - Almost every thing I know about using barycentric coordinates, and using matrices to solve linear equations, was learned last night because I found this question so interesting. So the following may be wrong, wrong, wrong - but some test values I have put in do seem to work. Guys and girls, please feel free to rip this apart if I screwed up completely - but here goes.

Finding barycentric coords in 3D space (with a little help from Wikipedia)

Given:

v0 = (x0, y0, z0)
v1 = (x1, y1, z1)
v2 = (x2, y2, z2)

p = (xp, yp, zp)

Find the barycentric coordinates: b0, b1, b2 of point p relative to the triangle defined by v0, v1 and v2

Knowing that:

xp = b0*x0 + b1*x1 + b2*x2
yp = b0*y0 + b1*y1 + b2*y2
zp = b0*z0 + b1*z1 + b2*z2

Which can be written as

[xp]      [x0]      [x1]      [x2]
[yp] = b0*[y0] + b1*[y1] + b2*[y2]
[zp]      [z0]      [z1]      [z2]

or

[xp]   [x0  x1  x2]   [b0]
[yp] = [y0  y1  y2] . [b1]
[zp]   [z0  z1  z2]   [b2]

re-arranged as

                   -1
[b0]   [x0  x1  x2]     [xp]
[b1] = [y0  y1  y2]   . [yp]
[b2]   [z0  z1  z2]     [zp]

the determinant of the 3x3 matrix is:

det = x0(y1*z2 - y2*z1) + x1(y2*z0 - z2*y0) + x2(y0*z1 - y1*z0)

its adjoint is

[y1*z2-y2*z1  x2*z1-x1*z2  x1*y2-x2*y1]
[y2*z0-y0*z2  x0*z2-x2*z0  x2*y0-x0*y2]
[y0*z1-y1*z0  x1*z0-x0*z1  x0*y1-x1*y0]

giving:

[b0]     [y1*z2-y2*z1  x2*z1-x1*z2  x1*y2-x2*y1]   [xp]
[b1] = ( [y2*z0-y0*z2  x0*z2-x2*z0  x2*y0-x0*y2] . [yp] ) / det
[b2]     [y0*z1-y1*z0  x1*z0-x0*z1  x0*y1-x1*y0]   [zp]

If you need to test a number of points against the triangle, stop here. Calculate the above 3x3 matrix once for the triangle (dividing it by the determinant as well), and then dot product that result to each point to get the barycentric coords for each point.

If you are only doing it once per triangle, then here is the above multiplied out (courtesy of Maxima):

b0 = ((x1*y2-x2*y1)*zp+xp*(y1*z2-y2*z1)+yp*(x2*z1-x1*z2)) / det
b1 = ((x2*y0-x0*y2)*zp+xp*(y2*z0-y0*z2)+yp*(x0*z2-x2*z0)) / det
b2 = ((x0*y1-x1*y0)*zp+xp*(y0*z1-y1*z0)+yp*(x1*z0-x0*z1)) / det

That's quite a few additions, subtractions and multiplications - three divisions - but no sqrts or trig functions. It obviously does take longer than the pure 2D calcs, but depending on the complexity of your projection heuristics and calcs, this might end up being the fastest route.

As I mentioned - I have no idea what I'm talking about - but maybe this will work, or maybe someone else can come along and correct it.

Answer 3

Update: Disregard, this approach does not work in all cases

I think I have found a valid solution to this problem.

NB: I require a projection to 2D space rather than working with 3D Barycentric co-ordinates as I am challenged to make the most efficient algorithm possible. The additional overhead incurred by finding a suitable projection plane should still be smaller than the overhead incurred when using more complex operations such as sqrt or sin() cos() functions (I guess I could use lookup tables for sin/cos but this would increase the memory footprint and defeats the purpose of this assignment).

My first attempts found the delta between the min/max values on each axis of the polygon, then eliminated the axis with the smallest delta. However, as suggested by @PeterTaylor there are cases where dropping the axis with the smallest delta, can yeild a straight line rather than a triangle when projected into 2D space. THIS IS BAD.

Therefore my revised solution is as follows...

1. Find each sub delta on each axis for the polygon { abs(v1.x-v0.x), abs(v2.x-v1.x), abs(v0.x-v2.x) }, this results in 3 scalar values per axis.
2. Next, multiply these scaler values to compute a score. Repeat this, calculating a score for each axis. (This way any 0 deltas force the score to 0, automatically eliminating this axis, avoiding triangle degeneration)
3. Eliminate the axis with the lowest score to form the projection, e.g. If the lowest score is in the x-axis, project onto the y-z plane.

I have not had time to unit test this approach but after preliminary tests it seems to work rather well. I would be eager to know if this is in-fact the best approach?

Answer 4

After much discussion there is actually a pretty simple way to solve the original problem of knowing which axis to drop when projecting to 2D space. The answer is described in 3D Math Primer for Graphics and Game Development as follows...

"A solution to this dilemma is to choose the plane of projection so as to maximize the area of the projected triangle. This can be done by examining the plane normal; the coordinate that has the largest absolute value is the coordinate that we will discard. For example, if the normal is [–1, 0, 0], then we would discard the x values of the vertices and p, projecting onto the yz plane."

My original solution which involved computing a score per axis (using sub deltas) is flawed as it is possible to generate a zero score for all three axis, in which case the axis to drop remains undetermined.

Using the normal of the collision plane (which can be precomputed for efficiency) to determine which axis to drop when projecting into 2D is therefore the best approach.

Answer 5

To project a point p onto the plane defined by the vertices v0, v1 & v2 you must calculate a rotation matrix. Let us call the projected point pd

e1 = v1-v0
e2 = v2-v0

r = normalise(e1)
n = normalise(cross(e1,e2))
u = normalise(n X r)

temp = p-v0

pd.x = dot(temp, r)
pd.y = dot(temp, u)
pd.z = dot(temp, n)

Now pd can be projected onto the plane by setting pd.z=0 Also pd.z is the distance between the point and the plane defined by the 3 triangles. i.e. if the projected point lies within the triangle, pd.z is the distance to the triangle.

Another point to note above is that after rotation and projection onto this plane, the vertex v0 lies is at the origin and v1 lies along the x axis.

HTH

Answer 6

I'm not sure that the suggestion is actually the best one. It's not too hard to project to the plane containing the triangle. I assume here that p is actually in that plane.

Let d1 = sqrt((v1-v0).(v1-v0)) - i.e. the distance v0-v1. Similarly let d2 = sqrt((v2-v0).(v2-v0))

v0 -> (0,0)
v1 -> (d1, 0)

What about v2? Well, you know the distance v0-v2 = d2. All you need is the angle v1-v0-v2. (v1-v0).(v2-v0) = d1 d2 cos(theta). Wlog you can take v2 as having positive y.

Then apply a similar process to p, with one exception: you can't necessarily take it as having positive y. Instead you can check whether it has the same sign of y as v2 by taking the sign of (v1-v0)x(v2-v0) . (v1-v0)x(p-v0).

---

As an alternative solution, you could use a linear algebra solver on the matrix equation for the tetrahedral case, taking as the fourth vertex of the tetrahedron v0 + (v1-v0)x(v2-v0) and normalising if necessary.

Answer 7

You shouldn't need to determine the optimal area to find a decent projection.

It's not strictly necessary to find the "best" projection at all, just one that's good enough, and that doesn't degenerate to a line when projected into 2D.

EDIT - algorithm deleted due to degenerate case I hadn't thought of