Jump to content
  • Advertisement
Sign in to follow this  
_Sauce_

Ray-triangle intersection

This topic is 3066 days old which is more than the 365 day threshold we allow for new replies. Please post a new topic.

If you intended to correct an error in the post then please contact us.

Recommended Posts

I'm working on my raytracer and in an effort to speed it up a bit, I'm ditching the algebraic method of ray-triangle intersection, and am replacing it with a parametric method, as described here. I'm gonna start talking some math stuff here so I made this diagram to help explain the method I'm using. Ro = origin of the ray. Rd = direction of the ray. P = point of intersection on the plane constructed from the three vertices of the triangle v1, v2, v3 = vertices of the triangle, in clockwise winding order. now according to the slides in the Princeton lecture, in order for P to be inside the triangle, a and b must both be in the range 0..1. My 3D math isn't the best but I recognised that a is simply the scalar projection of Z onto X, and b is the scalar projection of Z onto Y, where Z = P - v3, X = v2 - v3, and Y = v1 - v3. Correct? So here's my implementation...
float RTTriangle::intersect(const Rayf &ray)
{
	float distance = mPlane.intersect(ray);
	if(distance == std::numeric_limits<float>::infinity())
	{
		return distance;
	}

	Vector3f P = ray.mOrigin + ray.mDir * res.distance;

	Vector3f Z(P - mVerts[2]); //P - v3
	Vector3f X(mVerts[1] - mVerts[2]); //v2 - v3

	X.normaliseSelf();
	float a = Z.dot(X);
	if(a < 0 || a > 1)
	{
		return std::numeric_limits<float>::infinity();
	}

	Vector3f Y(mVerts[0] - mVerts[2]); //v1 - v3

	Y.normaliseSelf();
	float b = Z.dot(Y);
	if(b < 0 || b > 1)
	{
		return std::numeric_limits<float>::infinity();
	}

	return res;
}

And the output... Doesn't quite look like the famous Utah Teapot, now does it? Where exactly did I go wrong here?

Share this post


Link to post
Share on other sites
Advertisement
in case of CCW winding (I suggest you to use CCW)

if [(v3-p) cross (v1-p)] dot normal > 0
and [(v1-p) cross (v2-p)] dot normal > 0
and [(v2-p) cross (v3-p)] dot normal > 0

then we have a hit.

I use that for convex polygons (with any number of sides). Notice that you don't need to normalize that way at all, even the normal can be non-normalized.

There is a much simpler way too (google for methods, I can't recall the name), that will use some u,v parameters for checking.
It's much faster in case of triangles, and u,v coordinates can be used to calculate texture lookup coordinates with a very simple formula.
So I suggest you to use that.

[Edited by - szecs on February 21, 2010 6:44:49 AM]

Share this post


Link to post
Share on other sites
By the way, your logic is wrong too: You would have to project with the direction of the other side of the triangle (and that way you only checked that the ray is inside a parallelogram), but you are projecting in a perpendicular direction to the side, so you check two sides and a some curve, I can't decide what kind at the moment (maybe ellipsis-arc, but more likely not).

Share this post


Link to post
Share on other sites
Oh well, either I totally misunderstand something, or the slide that is cited by the OP is totally wrong!?

One can express p as
p = v3 + a * y + b * x
so that the condition for a triangle would be
a,b >=0 and a+b <= 1
The less restrict condition
0 <= a,b <= 1
would make up a parallelogram instead.

Now, looking at the slide, the image shows a and b being vectors but the text denotes them being scalars. Well ... okay.

Then, with a and b being scalars, the shown formula
P = a ( T2-T1 ) + b ( T3-T1 )
would indicate the same procedure as I've shown above: Adding 2 vectors that are along T2T1 and T3T1, resp., but then the difference vectors to P would not neccessarily be perpendicular to the sides.

Further, assume a triangle in 2D where
T1 = (0,0) ; T2 = (1,0) ; T3 = (0,1)
and a point of interest
P = ( 1,1 )
Doing the calculations gives me a=b=1, so that according to the slide, P is inside the triangle. But obviously, it isn't! Here the condition is wrong (see above).


Can anybody tell me what the slide is about?

Share this post


Link to post
Share on other sites
But that seems plain wrong. Do I miss something?

image
See: the point successes the test, but is sooo outside.

All points succeed that are inside the blue quad:
image

Share this post


Link to post
Share on other sites
Quote:
Original post by szecs
The slide is the same as the OP's drawing. But that seems plain wrong. Do I miss something?
Well, I came to the same conclusion. Still waiting for somebody to enlighten me ;)

Share this post


Link to post
Share on other sites
I guess you have to do it in all three directions, but in that case, my method seems a bit faster.

Share this post


Link to post
Share on other sites
Quote:
Original post by szecs
in case of CCW winding (I suggest you to use CCW)

if [(v3-p) cross (v1-p)] dot normal > 0
and [(v1-p) cross (v2-p)] dot normal > 0
and [(v2-p) cross (v3-p)] dot normal > 0

then we have a hit.

I use that for convex polygons (with any number of sides). Notice that you don't need to normalize that way at all, even the normal can be non-normalized.

There is a much simpler way too (google for methods, I can't recall the name), that will use some u,v parameters for checking.
It's much faster in case of triangles, and u,v coordinates can be used to calculate texture lookup coordinates with a very simple formula.
So I suggest you to use that.


A quick implementation of your method, doesn't seem to work. It looks less chaotic than my other attempts but it's still far from correct.

float RTTriangle::intersect(const Rayf &ray)
{
float distance = mPlane.intersect(ray);

if(distance != std::numeric_limits<float>::infinity())
{
Vector3f intersection(ray.mOrigin + ray.mDir * res.distance),
a(mVerts[0] - intersection),
b(mVerts[1] - intersection),
c(mVerts[2] - intersection),
normal = mPlane.getNormal();

if((c.cross(a)).dot(normal) > 0)
{
if((a.cross(b)).dot(normal) > 0)
{
if((b.cross(c)).dot(normal) > 0)
{
return distance;
}
}
}
}

return std::numeric_limits<float>::infinity();
}





[Edited by - _Sauce_ on February 22, 2010 7:26:14 AM]

Share this post


Link to post
Share on other sites
Quote:
Original post by haegarr
Quote:
Original post by szecs
The slide is the same as the OP's drawing. But that seems plain wrong. Do I miss something?
Well, I came to the same conclusion. Still waiting for somebody to enlighten me ;)
Sorry, I totally misunderstood your first sentence (haven't slept too much).

Share this post


Link to post
Share on other sites
bool Poin_On_Poly(int sides,float *point,float *poly_pts, float *norm)
{
float v1[3], v2[3], n[3];
float skal;

for( int i = 0; i < sides; i++ )
{ v1[0] = poly_pts[ i *3 + 0] - point[0];
v1[1] = poly_pts[ i *3 + 1] - point[1];
v1[2] = poly_pts[ i *3 + 2] - point[2];

v2[0] = poly_pts[((i+1)%sides)*3 + 0] - point[0];
v2[1] = poly_pts[((i+1)%sides)*3 + 1] - point[1];
v2[2] = poly_pts[((i+1)%sides)*3 + 2] - point[2];

vekt_szorz(v1,v2,n);
skal = skal_szorz(n,norm);

if( skal > 0 ) //I have no idea why not skal < 0
return false;
}

return true;
}

Share this post


Link to post
Share on other sites
Sign in to follow this  

  • Advertisement
×

Important Information

By using GameDev.net, you agree to our community Guidelines, Terms of Use, and Privacy Policy.

We are the game development community.

Whether you are an indie, hobbyist, AAA developer, or just trying to learn, GameDev.net is the place for you to learn, share, and connect with the games industry. Learn more About Us or sign up!

Sign me up!