3d point-in-triangle algorithm

Started by
25 comments, last by EpicArtist 16 years, 2 months ago
I am in need of a fast 3d point-in-triangle algorithm.
Advertisement
There is this new thing called Google. All the cool kids are doing it.
First result of the Google.
"When you die, if you get a choice between going to regular heaven or pie heaven, choose pie heaven. It might be a trick, but if it's not, mmmmmmm, boy."
How to Ask Questions the Smart Way.
Google.

http://www.google.ca/search?hl=en&q=3d+point+in+triangle+algorithm&btnG=Google+Search&meta=

I don't know the context in which you're using this algorithm so I couldn't begin to give you hints on optimization. Give us some context and maybe you'll get more help.

EDIT: CodeMunkie, you beat me to it. :P

I used this function for ray collisions with an md2 mesh:

bool point_inside_triangle(g_vector contact_point, triangle * P)	{		double accumilator = 0;	double dot1 = 0;	double dot2 = 0;	double dot3 = 0;		g_vector V1,V2,V3,line1,line2;		V1=P->vertex_array[0];	V2=P->vertex_array[1];	V3=P->vertex_array[2];		//first angle	line1=(V1-contact_point);line1.normalize();	line2=(V2-contact_point);line2.normalize();	dot1=dot_product(line1,line2);	accumilator += acos(dot1);		//second angle	line1=(V2-contact_point);line1.normalize();	line2=(V3-contact_point);line2.normalize();	dot2=dot_product(line1,line2);	accumilator += acos(dot1);		//third angle	line1=(V3-contact_point);line1.normalize();	line2=(V1-contact_point);line2.normalize();	dot3=dot_product(line1,line2);		accumilator = acos(dot1) + acos (dot2) + acos(dot3);		if(accumilator < (359.9*pirad))		return false;	else		return true;};


One advantage of this is that you can pre-calculate the lines of the triangle allowing you to simply do three dot products, three multiplies, and two additions.
Don't thank me, thank the moon's gravitation pull! Post in My Journal and help me to not procrastinate!
A couple of different approaches I've taken:

After finding the line intersection point on the triangle plane you can:

1) Collapse on an axis to create a 2d problem. Take a line from the point to an aritrary point outside the triangle on one of the remaining axis'. If this crosses an odd number of edges your intersection point is inside, if not its outside. I believe this works for non convex planar polygons too and though I haven't used it in years remember it as being good for non vectorized hardware.

2) Store edge normals (normals pointing consistently into or out of the triangle lying orthogonal to the each edge and the plane normal). Use these for reasonably (very with sse5 / ps3) fast rejections of points outside the triangle. You can precalculate compress these normals to 64bit values use multi media extractions to 128 bit values and still get pretty accurate results. One of my ex-colleages did some excellent compression to create a packed quad representation (we collapsed an edge for triangles and used the same math!) consisting of two points and four edge normals. The remaining two points were rarely needed and could be sythesised quicker than the cache misses caused by storing them and bloating data.

This second approach boils down to (for a triangle) three subractions, three dot products and optional early outs.
I am fully aware of BlackPawn's website, his method only works for 2 dimensions, although I haven't implemented the first I'm pretty sure it would too.

SpeciesUnkown, could you post a version with a more psuedo syntax? I'm not quite sure I known what everything means (for example,V1=P->vertex_array[0];.)

Tom, I am not entirely sure your second method would work in 3d without performing a plane intersection test.

Quote:Original post by EpicArtist
I am fully aware of BlackPawn's website, his method only works for 2 dimensions, although I haven't implemented the first I'm pretty sure it would too.


I'm not sure what you mean by this - determining if a point is inside a polygon is conceptually a 2D problem. Both methods given work perfectly well for 3D if the point is on the plane of the polygon.
Quote:Original post by EpicArtist
I am fully aware of BlackPawn's website, his method only works for 2 dimensions, although I haven't implemented the first I'm pretty sure it would too.


In fact, triangles are exclusively a 2d primitive and, despite the number of dimensions, the point-in-polygon test algorithm is always in 2d.

Quote:
SpeciesUnkown, could you post a version with a more psuedo syntax? I'm not quite sure I known what everything means (for example,V1=P->vertex_array[0];.)

Tom, I am not entirely sure your second method would work in 3d without performing a plane intersection test.


OK, I just sort of dropped my code in as it was. each of my triangles has a pointer to 3 vertices, hence P->vertex_array[0] for the first point, P->vertex_array[1]for the second etc.

Here is the psudocode version, unoptimised (mine is designed to reuse variables):


BOOL point_in_triangle (VECTOR3 c, TRIANGLE p)  VECTOR3 L1,L2,L3  FLOAT angle  FLOAT dot1,dot2,dot3  L1 = [vector from c to p.point(0)]  L2 = [vector from c to p.point(1)]  L3 = [vector from c to p.point(2)]  normalise(L1)  normalise(L2)  normalise(L3)  dot1 = dot_product(L1,L2)  dot2 = dot_product(L2,L3)  dot3 = dot_product(L3,L1)  angle = acos(dot1) + acos(dot2) + acos(dot3)  IF angle CLOSE TO 2*pi    RETURN true  ELSE    RETURN false  END IFRETURN


This works on the principle that all the internal angles of a shape should add up to 2PI. You could also do it within a for loop and therefore it will work on any convex polygon. (although it fails miserably for concave ones). This can tell you if any point is on the polygon surface, even if it is not on the polygon's plane. This test fails on concave polygons, but it works on polygons in 3d space, since it simply relies on the functions dot_product and "vector from" both being in n dimensions.

n.b. can somebody tell me what the type property of the source bbcode tag is for psudocode? My psudocode always gets changed to c++. THX
[/source]
Don't thank me, thank the moon's gravitation pull! Post in My Journal and help me to not procrastinate!
Fastest way I know of is using the barycentric technique described in the link in the first response. You can find a pretty fast version of it in the book "Real Time Rendering" by Moller and Haines.
Quote:Original post by Leo_E_49
Fastest way I know of is using the barycentric technique described in the link in the first response. You can find a pretty fast version of it in the book "Real Time Rendering" by Moller and Haines.


Unfortunately, the OP asked for one that works in 3d, as the sum of angles method does. Otherwise, I would have suggested the barycentric technique. If you know that the point is on the plane of the polygon (e.g. it came from a line-vs-plane test) use the barycentric method (its faster). If you are not sure, use the sum of angles method.

You could do a precursory point-on-plane test followed by the barycentric method, thus adding a third dimension.
Don't thank me, thank the moon's gravitation pull! Post in My Journal and help me to not procrastinate!

This topic is closed to new replies.

Advertisement