Jump to content
  • Advertisement
Sign in to follow this  
xcodedave

OpenGL problem with ray/box intersection in GLSL

This topic is 647 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

Hi all,
 
For part of an OpenGL ray-tracer I'm working on (as a hobby) I'm attempting to find the interval t in which a ray intersects an area-aligned-bounding-box (AABB).

I've created a GLSL function (left un-optimised for clarity) based on the slab method described here: https://tavianator.com/fast-branchless-raybounding-box-intersections/:

#version 330

const float INFINITY = 1.0 / 0.0;

struct AABB
{
    vec3 min;
    vec3 max;
};

struct Ray
{
    vec3 origin;
    vec3 direction;
};

float rayAABBIntersection (Ray ray, AABB aabb)
{
        vec3 inverseRayDirection = 1.0 / ray.direction;

	float t1;
	float t2;
	float tmin = -INFINITY;
	float tmax = INFINITY;

	if (ray.direction.x != 0.0) 
	{
	    t1 = (aabb.min.x - ray.origin.x) * inverseRayDirection.x;
	    t2 = (aabb.max.x - ray.origin.x) * inverseRayDirection.x;
	    tmin = max(tmin, min(t1, t2));
	    tmax = min(tmax, max(t1, t2));
	}
    
	if (ray.direction.y != 0.0) 
	{
	    t1 = (aabb.min.y - ray.origin.y) * inverseRayDirection.y;
	    t2 = (aabb.max.y - ray.origin.y) * inverseRayDirection.y;
	    tmin = max(tmin, min(t1, t2));
	    tmax = min(tmax, max(t1, t2));
	}
	
	if (ray.direction.z != 0.0) 
	{
	    t1 = (aabb.min.z - ray.origin.z) * inverseRayDirection.z;
	    t2 = (aabb.max.z - ray.origin.z) * inverseRayDirection.z;
	    tmin = max(tmin, min(t1, t2));
	    tmax = min(tmax, max(t1, t2));
	}
	
	if (tmin < 0)
	{
	    return tmax;
	}
	else
	{
	    return tmin;
	}
}

I've tested the function first in C by replacing `min` and `max` with `fmin` and `fmax` functions from <math.h>.

 

In C the code seems to work fine and correctly returns the t value for the ray intersection, however in the GLSL depending on the ray position and angle this gives me incorrect results.

 

If anyone has any insight about why the GLSL isn't functioning as expected, I'd be very grateful!

 

(Note: I'm debugging on OSX and have had problems building glslDebugger, so have had to resort to outputting colours to assist in debugging which proves very time consuming - if anyone knows of any Mac based tools which could help, I'd be grateful for those also).

 

 

Share this post


Link to post
Share on other sites
Advertisement

I could be mistaken but I guess that in 3.3 GLSL, division by zero would result in undefined behavior. From page 21 of the official spec (https://www.opengl.org/registry/doc/GLSLangSpec.3.30.6.pdf) :

 

"Similarly, treatment of conditions such as divide by 0 may lead to an unspecified result, but in no case should such a condition lead to the interruption or termination of processing. Generally, there are no signaling NaNs, and operating on NaNs (Not a Number) or infs (positive or negative infinities) gives undefined results."

Share this post


Link to post
Share on other sites

The division by zero is indeed a problem with this approach. E.g. the Ember raytracer just uses FLT_MIN if one of the direction elements is zero (at least in one of their older versions). Another option is to use the Separating Axis Test (SAT) for ray/box instead. This will not clip the ray though, but only provide a boolean answer whether the ray intersects the box. I use this in most of my tree traversals as I found it faster and more robust.

Share this post


Link to post
Share on other sites

Shameless self-promotion: http://www.gamedev.net/page/resources/_/technical/math-and-physics/intersection-math-algorithms-learn-to-derive-r3033 - this is probably more optimal algorithm for the GPU (I've been using (similar version to this) in OpenCL/CUDA for quite a bit of time, I'm sure GLSL is able to handle it too).

 

Don't try to build your code to avoid divisions by 0, especially in world of real time ray tracing, build it in such way so you can handle division by zeroes (use the properties of IEEE floating point math - like INF value or NAN value).

 

Let me try to dig something on my AABB intersection here on local machine.

 

EDIT:

 

To add here - this is my gpu-based Ray-AABB test done when working with KD-tree (and testing against AABB sorrounding whole tree):

float4 v1 = (aabbMin - o) * invRay;
float4 v2 = (aabbMax - o) * invRay;
float4 n = min(v1, v2);
float4 f = max(v1, v2);
float enter = max(n.x, max(n.y, n.z));
float exit = min(f.x, min(f.y, f.z));

if (exit > 0.0f && enter < exit)
{
    // We have a hit
}

Where:

  • aabbMin - is float4 which contains min point in X, Y and Z members and 1.0 in W member
  • aabbMax - is float4 which contains max point in X, Y and Z members and 1.0 in W member 
  • o - is ray origin point with 1.0 in W member
  • invRay - is float4(1.0 / rayDirection.x, 1.0 / rayDirection.y, 1.0 / rayDirection.z, 0.0)

Functions min (and respectively max) are working per component, like:

float4 min(float4 a, float4 b)
{
    return (float4)(min(a.x, b.x), min(a.y, b.y), min(a.z, b.z), 0.0);
}
Edited by Vilem Otte

Share this post


Link to post
Share on other sites

Thanks for the advice everyone!

 

@SapphireG and @Dirk the advice on the handling of infinite and NaN values was very helpful - it caused me to test fundamentally the output of min and max functions with NaN values...this has helped me confirm that my rendering error is not caused by this function.

 

@Vilem - thanks for sharing your code - your Ray AABB test looks very similar to an optimised one I had been testing:

float rayAABBIntersection (Ray ray, AABB aabb)
{
   float t[10];
   t[1] = (aabb.min.x - ray.origin.x) * ray.inverseDirection.x;
   t[2] = (aabb.max.x - ray.origin.x) * ray.inverseDirection.x;
   t[3] = (aabb.min.y - ray.origin.y) * ray.inverseDirection.y;
   t[4] = (aabb.max.y - ray.origin.y) * ray.inverseDirection.y;
   t[5] = (aabb.min.z - ray.origin.z) * ray.inverseDirection.z;
   t[6] = (aabb.max.z - ray.origin.z) * ray.inverseDirection.z;
   t[7] = max(max(min(t[1], t[2]), min(t[3], t[4])), min(t[5], t[6]));
   t[8] = min(min(max(t[1], t[2]), max(t[3], t[4])), max(t[5], t[6]));
   t[9] = (t[8] < 0 || t[7] > t[8]) ? INFINITY : t[7];
   return t[9];
}

but yours looks to be slightly more optimal (using vector operations to remove a few lines by the looks of it)!

I hope you don't mind if I use yours for inspiration! ;)

 

I've since moved to working in OpenGL 4.1 (GLSL #version 410) and by checking the results of `min` and `max` have confirmed they work the same way as my C code when presented with NaN values - (if one of the values is NaN, the other is always returned).

 

I think the error in my ray-tracing must be elsewhere....so I'd like to quickly break down the steps I'm taking to calculate a ray to see if my approach is valid:

  1. I upload a simple cube VBO to the GPU.
  2. In my vert shader I pass the vert-position (before applying a model-view-projection transform) to the fragment shader. This should give me the surface position for a fragment in model space. (Correct?)
  3. I multiply the eye position (in world space) by the inverse model-view matrix to get the eye position in model space and pass this to the fragment shader.
  4. To calculate a ray for each fragment, I use the eye-position as the origin, and the direction is calculated as the fragment surface position, minus the eye position, normalized.

Is this a valid approach for calculating a (model-space) ray per fragment?

 

Many thanks again for all your help, it's greatly appreciated!

Share this post


Link to post
Share on other sites

 

) I'm attempting to find the interval t in which a ray intersects an area-aligned-bounding-box (AABB).

 

I've created a GLSL function (left un-optimised for clarity) based on the slab method described here

 

 

you define z_far thats right the same z_far you define but you store it in shader

 

in vertex shader you pass to fragment vertex position as: aabb_world_matrix * v 

 

in fragment shader

float t = DistanceBetweenPoints(Camera_position, world_vertex_pos) / z_far;

 

where in vertex shader

world_vertex_pos.x = dp43(WM1, Vpos);
world_vertex_pos.y = dp43(WM2, Vpos);
world_vertex_pos.z = dp43(WM3, Vpos);
float dp43(vec4 matrow, vec3 p)
{
return ( (matrow.x*p.x) + (matrow.y*p.y) + (matrow.z*p.z) + matrow.w );
}

where WM1 is vec4 and stores world matrix row 1

where WM2 is vec4 and stores world matrix row 2

where WM3 is vec4 and stores world matrix row 3

where WM4 is vec4 and stores world matrix row 4

 

and Vpos is an attribute that i pass ass vertex buffer anyway its old glsl

 

hypnotized-smiley.gif

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!