# Ray-triangle intersection, triangle edges being missed?

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

## Recommended Posts

I'm currently implementing a system to fire many rays against a triangle mesh, calculating intersections. Everything works great in general, except for the occasional case of rays firing between adjoining triangles without detecting a hit on either. I'm only using floats for all the math so I know my precision isn't as high as it could be, but the triangles are all using shared vertices, etc so they should be nicely locked together without 'real' gaps.

Can anyone suggest a nice way to handle this and detect this case where my rays would be striking the 'edge' of the triangle? Performance is a a big consideration. Cheers!

##### Share on other sites
I guess you should post the math for testing the intersection so someone can suggest a better algorithm or add something to it to make it work properly.

##### Share on other sites
Good call, Waterlimon, I was hoping for a general 'this algorithm handles edges nicely', but if anyone wants to suggest an improvement to my existing algorithm, thar she blows:

 //return values: // -1 triangle is degenerate, no processing done // 0 no intersection // 1 valid intersection // 2 ray lies on triangle plane //rayOrigin = x,y,z floats //rayDir = normalized x,y,z floats //srcTriangleVert* = x,y,z floats //triNormal = x,y,z triangle normal //destIntersectData = x,y,z intersect + distance to intersect. (< 0 if ray was headed in wrong direction) int intersect_RayTriangleCPU_Native( float *rayOrigin, float *rayDir, float *srcTriangleVertA, float *srcTriangleVertB, float *srcTriangleVertC, float *triNormal, float *destIntersectData ) { float u[3], v[3]; // triangle vectors float w0[3], w[3]; // ray vectors float r, a, b; // params to calc ray-plane intersect float uu, uv, vv, wu, wv, D; float s, t; // get triangle edge vectors and plane normal u[0] = srcTriangleVertB[0] - srcTriangleVertA[0]; u[1] = srcTriangleVertB[1] - srcTriangleVertA[1]; u[2] = srcTriangleVertB[2] - srcTriangleVertA[2]; v[0] = srcTriangleVertC[0] - srcTriangleVertA[0]; v[1] = srcTriangleVertC[1] - srcTriangleVertA[1]; v[2] = srcTriangleVertC[2] - srcTriangleVertA[2]; w0[0] = rayOrigin[0] - srcTriangleVertA[0]; w0[1] = rayOrigin[1] - srcTriangleVertA[1]; w0[2] = rayOrigin[2] - srcTriangleVertA[2]; //manually compute dot product and negate result a = -((triNormal[0] * w0[0]) + (triNormal[1] * w0[1]) + (triNormal[2] * w0[2])); //manually compute dot product b = ((triNormal[0] * rayDir[0]) + (triNormal[1] * rayDir[1]) + (triNormal[2] * rayDir[2])); if (fabs(b) < SMALL_NUM) // ray is parallel to triangle plane { if (a == 0.0) // ray lies in triangle plane { return 2; } else { return 0; // ray disjoint from plane } } // get intersect point of ray with triangle plane r = a / b; if (r < 0.0f) // ray goes away from triangle { return 0; // => no intersect } // for a segment, also test if (r > 1.0) => no intersect (we are firing raw rays, no endpoint) destIntersectData[0] = rayOrigin[0] + (rayDir[0] * r); // intersect point of ray and plane destIntersectData[1] = rayOrigin[1] + (rayDir[1] * r); // intersect point of ray and plane destIntersectData[2] = rayOrigin[2] + (rayDir[2] * r); // intersect point of ray and plane //hijack the 4th element to return distance to intersection destIntersectData[3] = r; // is I inside srcTriangle? //manually compute dot products uu = (u[0] * u[0]) + (u[1] * u[1]) + (u[2] * u[2]); uv = (u[0] * v[0]) + (u[1] * v[1]) + (u[2] * v[2]); vv = (v[0] * v[0]) + (v[1] * v[1]) + (v[2] * v[2]); w[0] = destIntersectData[0] - srcTriangleVertA[0]; w[1] = destIntersectData[1] - srcTriangleVertA[1]; w[2] = destIntersectData[2] - srcTriangleVertA[2]; wu = (w[0] * u[0]) + (w[1] * u[1]) + (w[2] * u[2]); wv = (w[0] * v[0]) + (w[1] * v[1]) + (w[2] * v[2]); D = uv * uv - uu * vv; // get and test parametric coords s = (uv * wv - vv * wu) / D; if (s < 0.0 || s > 1.0) // I is outside T { return 0; } t = (uv * wu - uu * wv) / D; if ((t < 0.0) || ((s + t) > 1.0)) // I is outside T { return 0; } return 1; // I is in T } 

##### Share on other sites
Just so I can understand correctly, Basically you have an object lets say its a model, and for the sake of argument lets say its a door.

This door is made up of a bunch of triangles, and your issue is that your code works correctly when it comes to picking seperate triangles, but sometimes when moving your mouse over the door, there are moments when the picking is not detected, and what you want is that whenever your mouse is over any part of that door, the door infact is recognized as being "highlighted" over the mouse.

Is this correct?

##### Share on other sites

Just so I can understand correctly, Basically you have an object lets say its a model, and for the sake of argument lets say its a door.

This door is made up of a bunch of triangles, and your issue is that your code works correctly when it comes to picking seperate triangles, but sometimes when moving your mouse over the door, there are moments when the picking is not detected, and what you want is that whenever your mouse is over any part of that door, the door infact is recognized as being "highlighted" over the mouse.

Is this correct?

Yeah exactly. In my case I have a triangle mesh such as a door, and I'm doing something very similar to picking in order to check where a ray intersects the mesh. Each triangle is tested individually. Normally it works perfectly, but as my picking ray crosses over two triangles of the door, it sometimes manages to pass between the two triangles, generating no hit with either. The triangles share vertices, so they don't really have gaps beyond the margin of error from float precision.

I've noticed that barycentric coord based solutions can technically detect edge cases such as this, but I wonder if in reality they have the same issues with float precision? I'm really hoping to avoid 'fudging' the numbers to make my triangles bigger than they are as a solution.

##### Share on other sites
I honestly haven't taken a thorough look at your code, but when you get the intersection point, how are you calculating if the point falls inside the triangle?

For example the method I use is that i calculat the angles that the point makes with respect to the 3 points of the triangle, if the sum of the angles it makes equals 2*pi or 360 degrees, its inside the triangle.

Is this what your using? or another method?

##### Share on other sites
I did something very similar in my ray tracer for ray-triangle intersection. For what it's worth, I got the ray-plane intersection from Wikipedia and the triangle intersection from SoftSurfer. While I used pretty much the same method that Grumple posted, I'll post my code just to give another perspective/implementation that will hopefully help it be a little more clear:

struct Triangle { Vector3f a, b, c; Vector3f normal; // normalized bool intersects(const Ray& ray, float& t, Vector3f& n) const { // from http://softsurfer.com/Archive/algorithm_0105/algorithm_0105.htm (and Plane.intersect()) n = normal; float num = normal.dot(a - ray.start); float den = normal.dot(ray.direction); if (num != 0 && den == 0) { return false; // parallel, but never intersect } else if (num == 0 && den == 0) { t = 0; // parallel and in same plane return true; } else { t = num / den; } if (t < 0) { return false; // intersection point is behind the ray's start } Vector3f u = b - a; Vector3f v = c - a; Vector3f w = ray.start + ray.direction * t - a; float uv = u.dot(v); float wv = w.dot(v); float vv = v.squaredLength(); float wu = w.dot(u); float uu = u.squaredLength(); float s1 = (uv * wv - vv * wu) / (uv * uv - uu * vv); float t1 = (uv * wu - uu * wv) / (uv * uv - uu * vv); return s1 >= 0 && t1 >= 0 && (s1 + t1) <= 1; } };

It's important to be aware that this degenerates when the ray lies in the triangle's plane, because I simply set [font=courier new,courier,monospace]t = 0[/font], when really [font=courier new,courier,monospace]t = 0[/font] could be outside the triangle. You could optionally return false if [font=courier new,courier,monospace]t = 0[/font] and since a triangle is infinitely thin as well. It just depends on how you want to handle this situation.

##### Share on other sites
Thanks for posting that, Cornstalks. Your implementation is essentially identical (although admittedly more readable =P) than mine. Referring to your code, in situations where I should be detecting intersection with the triangle, I'm getting values of s1 along the lines of -0.0008, or (s1 + t1) of 1.0002.

I guess I could accept fudging my boundaries a bit if I could easily relate s1 and t1 to the units of measurement the overall world adheres to. For example if I can guarantee I am adding a set number of centimeters to my triangle edges as opposed to a fraction of the unknown triangle size.

Have you run into this issue with your implementation?

##### Share on other sites
That's a strange one, and I'm not sure how this edge case is happening. Can you characterise the nature of the bug, e.g. find an angle between triangles or camera angle where this happens repeatably? You could perhaps create a test harness to try a bunch of options with known outcomes (e.g. two triangles meeting at a known point, loop through different angles and test points) until you find a situation where it occurs.

One slightly unrelated question, do you actually *need* to ray test against a raw mesh? Often games use simpler models for collision etc for speed reasons. A simpler object may yield faster code, e.g. ray testing against transformed cubes, spheres, cylinders, etc is very easy.

##### Share on other sites
With really small errors like that, it could very easily be floating point rounding errors. Try replacing the intersection code with double instead of float and see if the error is smaller (this is just to make sure it is indeed floating point accuracy errors).

I'm going to go ahead and assume it is indeed floating point rounding issues, and to be honest there isn't a fool-proof way around this. One commonly used solution would be to use some epsilon to allow for a controllable degree of error (which is similar to what you mentioned). For example:

float epsilon = 0.001f; return (s1 >= 0.0f - epsilon) && (t1 >= 0.0f - epsilon) && ((s1 + t1) <= 1.0f + epsilon); 

Another thing would be to look at your compiler's optimization settings. There are usually some flags that can be set in the compiler that will trade accuracy for speed (usually I've found it to be accuracy is favored over speed by default, and you have to explicitly say to favor speed over accuracy). Double check them and make sure you aren't sacrificing accuracy for speed.

I've never had precision errors this small cause any problems for me though, so I can't give any more solutions. It also depends on your data. If you're mixing large numbers and small numbers in calculations, the errors will be much more noticeable. Try to make sure your floats don't differ by too many orders of magnitude if you want more precision.

1. 1
2. 2
3. 3
4. 4
5. 5
Rutin
15

• 14
• 9
• 9
• 9
• 10
• ### Forum Statistics

• Total Topics
632912
• Total Posts
3009199
• ### Who's Online (See full list)

There are no registered users currently online

×