Discontinuous ray casting through 3D texture

Hi,

I am doing volume rendering of a 3D texture using a basic ray casting algorithm.

Here is how it works so far : the fragment shader iterates along a ray directed from the eye to the scene, taking samples from the 3D texture at each step. Frag color is generated based upon the data collected along the ray (for example, we may want to display the maximal intensity along the ray).

As long as we just want to scan the whole 3D texture, that’s easy : we send a cube of GL_TRIANGLES through the pipeline, and the vertex shader computes ray entrance location (or ray exit but we use face culling to discard back faces). Then ray exit can be computed easily in the fragment shader based on ray entrance.

But I would like to display only a subset of the whole 3D texture, defined by some arbitrary bounding geometry. The idea is to be able to discard some parts of the 3D texture without having to change the texture itself. So, instead of sending a cube of vertices to the pipeline, one could send, let’s say, a pyramid, and the iteration that takes place in the fragment shader would stay into the bounds of that pyramid.

As long as the 3D texture subset is a convex polyhedron, it’s fine : front and back faces locations can be first stored (as colors) in two separate FBOs, and then passed to the fragment shader as 2D textures.

Where it gets tricky is when the geometry is not a convex polyhedron, because then a single ray can enter and exit the geometry several times :

illustration : [ATTACH=CONFIG]1477[/ATTACH]

So my question is : is there a convenient way in OpenGL to do such a “piecewise ray casting” ?

There are too major issues I could not overcome :
[ul]
[li]How to compute the - possibly multiple - [ray entrance, ray exit] pairs using the GPU ?
[/li][li]How to transmit those pairs to the fragment shader so that, ideally, a single fragment shader execution can process a complete discontinuous ray ?
[/li][li]
[/li][/ul]
And if none or those are possible, is there another approach I could try ?

[QUOTE=MrGruk;1287455]There are too major issues I could not overcome :
[ul]
[li]How to compute the - possibly multiple - [ray entrance, ray exit] pairs using the GPU ?[/li][li]How to transmit those pairs to the fragment shader so that, ideally, a single fragment shader execution can process a complete discontinuous ray ?[/li][/ul][/QUOTE]

you have to capture the “arbitrary geometry” first, this can be done with “transform feedback objects”. after you’ve captured all triangles (into a buffer object on te GPU), check each of them for intersection with the ray in another shader execution, that way you can determine the points where the ray enters / exits the geometry. i’ve no clue how you can determine wheater a certain part of the ray is within / outside the geometry … perhaps by sorting the intersection points after their depth and assume that “alternating” parts are …

outside [intersection1] within [intersection2] outside [intersection3] within … etc

It took me quite a while, but thanks to your help I finally managed to make something that works. So thank you very much !

It is far from perfect. There are still some issues I would like to fix (see below). But it gets the job done, so here are the main steps :
[ul]
[li]Put the bounding geometry in some VBO. It will be used as vertex shader input for both passes
[/li][li]1st pass : storing GL coords of bounding geometry using a feedback buffer
[/li][LIST]
[li]Disable face culling an depth testing so every triangle can be processed
[/li][li]Enable GL_RASTERIZER_DISCARD in order to stop the pipeline before fragment shader execution
[/li][li]create an OpenGL buffer to store transform feedback
[/li][li]bind it using glBindBufferBase(GL_TRANSFORM_FEEDBACK_BUFFER, 0, bufferId);
[/li][li]in vertex shader, just transform vertices to GL coordinates
[/li][li]in geometry shader, compute the orientation of each triangle given the
[/li]current camera direction
[li](there is no fragment shader)
[/li][li]execute shader program between glBeginTransformFeedback(GL_TRIANGLES) and
[/li]glEndTransformFeedback() calls in order to write on feedback buffer instead of doing a regular rendering. Use bounding geometry as vertex shader input.
In the feedback buffer object, store, for each triangle :
[LIST]
[li]GL coordinates and corresponding 3D texture coordinates
[/li][li]triangle orientation (in GL coordinates, and considering current camera)
[/li][/ul]
[/LIST]
[li]2nd pass : actual rendering
[/li][ul]
[li]Enable face culling an depth testing so that, for each ray, only the first entrance point is processed
[/li][li]Disable GL_RASTERIZER_DISCARD
[/li][li]attach the data of the feedback buffer to a texture buffer object (TBO) using glTexBuffer() so it can be used as array input by the fragment shader
[/li][li]in vertex shader, just transform vertices to GL coordinates
[/li][li]in fragment shader :
[/li][LIST]
[li]use “uniform samplerBuffer yourTboName;” to get access to the TBO
[/li][li]parse the triangle list contained in the TBO using texelFetch(), and make the list of intersection points between the ray and the triangles. I don’t use perspective, so this is purely 2D calculations for me
[/li][li]sort that list based on z coordinates of intersection points (from near to far)
[/li][li]parse the list an perform the ray casting algorithm for each [ray entrance, ray exit] pair. To know if an intersection is ray entrance of ray exit, just read it from the TBO
[/li][li]check for consecutive ray entrance or consecutive ray exit points and make sure to handle those properly. It may happen quite often due to float imprecision
[/li][/ul]
[li]execute shader program with bounding geometry as input
[/li][/LIST]
[/LIST]

The clumsy part is the fragment shader of the 2nd pass. There are many issues :
[ul]
[li]Each fragment individually needs to parse all triangles, so performance will drop as the bounding geometry becomes more complex.
[/li][li]It feels so wrong to check for ray - triangles intersections in the fragment shader because - if my understanding of OpenGL pipeline is correct - OpenGL already computes those intersection points internally before running the fragment shader (or discarding fragments). And I assume it does it way faster than I do. Would there be a way to avoid recalculating those intersection points in the fragment shader ?
[/li][li]GLSL is not well-suited for making lists and sorting them. The only way I found to store intersection points is by using a static array in my fragment shader. That means that if there are more intersection points than there is room in the static array, it doesn’t work. But if the static array size is too big (let’s say 256), then the performance drops dramatically due to the only presence of the array (even if most of the array is actually never used - and of course, not involved in the sorting process). Is there a better way I could handle this in GLSL ?
[/li][/ul]

regarding the first pass:
– you dont need to disable depth testing because it takes place after primitive assembly & rasterization

regarding hte second pass:
– you can use a “shader storage buffer” instead of a “texture buffer” (if your system supports GL 4.3)

regarding sorting intersection points by depth:
– you can do that with a “compute shader” for example, in a “1.5th pass”, your input can be a shader storage buffer containing the intersection points (vec4), and an output (shader storage) buffer in which you write the sorted points. the problem with the “static array” size doesnt exist for shader storage buffers because you can query their size on runtime (besides reading / writing to these buffers)

now is depends on what you want to do with these points (didnt understand fully what you’re trying to render). calculations are usually better done on a “per-vertex”-level because there are very likely less vertices than fragments, so there are many fragment shader invocations.

EDIT:
ok, your fragment color is the sum of several points on the ray within the “arbitrary geometry”. maybe you’ve heared about “order-independent transparency” (OIT), there is an example program for that in the book “OpenGL Programming Guide 8th Edition” (Chapter 11 “Memory”). it may help you with our problem …

how it works:
for each pixel on screen you create a “linked list” of fragments belonging to this pixel. in a 2nd pass, you collect all the “hits” for each pixel, calculate the resultion fragment color and write it to the framebuffer.

therefore you need:
– a “big buffer” containing ALL fragment hits
– a integer texture of the screen size
----> contains the location (array index) of the first hit at this pixel in “big buffer”
– a atomic counter buffer (for indexing)

a “fragment hit” looks like this:
–> vec4 fragment color
–> vec4 fragment position / depth
–> uint nexthit (array index)

APPROACH #2:

step 1: render “arbitrary geometry”, capture triangles, discard fragments

step 2: calculate ray x triangle intersection points, capture those into another buffer

step 3: generate line segments that are within the “arbitrary geomatry”

step 4: render line segments, subdivide these in the geometry shader to severyl point primitives, capture generate fragments into “big buffer” for hits

step 5: for each screen pixel, blend the captured fragments together to the “final color”

maybe the need for the “ray casting” dissapears with an OIT approach (?)

regarding the first pass:
– you dont need to disable depth testing because it takes place after primitive assembly & rasterization

Indeed. Thanks for pointing that out.

ok, your fragment color is the sum of several points on the ray within the “arbitrary geometry”.

The way ray points are processed depends on the volume rendering technique used, so it is not always as simple as a sum. For example, I want to be able to :

[ul]
[li]find the maximal intensity point along the ray (my 3D texture is monochromatic), and use it as fragment color. This makes a basic volume rendering, like here (I already do this kind of rendering, but now I want to be able to skip some parts of the 3D texture)
[/li][li]use additional rendering tricks, such as adding some intensity attenuation as the ray goes deeper into the 3D texture
[/li][li]find the nearest point that has its intensity over a given threshold, in order to make surface rendering
[/li][li]etc
[/li][/ul]

maybe you’ve heared about “order-independent transparency” (OIT), there is an example program for that in the book “OpenGL Programming Guide 8th Edition” (Chapter 11 “Memory”). it may help you with our problem …

how it works:
for each pixel on screen you create a “linked list” of fragments belonging to this pixel. in a 2nd pass, you collect all the “hits” for each pixel, calculate the resultion fragment color and write it to the framebuffer.

therefore you need:
– a “big buffer” containing ALL fragment hits
– a integer texture of the screen size
----> contains the location (array index) of the first hit at this pixel in “big buffer”
– a atomic counter buffer (for indexing)

a “fragment hit” looks like this:
–> vec4 fragment color
–> vec4 fragment position / depth
–> uint nexthit (array index)

This looks extremely promising. I think this is exactly what I was looking for. Thanks a lot ! I will definitely give it a try and provide some feedback.

APPROACH #2:

step 1: render “arbitrary geometry”, capture triangles, discard fragments

step 2: calculate ray x triangle intersection points, capture those into another buffer

step 3: generate line segments that are within the “arbitrary geomatry”

step 4: render line segments, subdivide these in the geometry shader to severyl point primitives, capture generate fragments into “big buffer” for hits

step 5: for each screen pixel, blend the captured fragments together to the “final color”

maybe the need for the “ray casting” dissapears with an OIT approach (?)

This sounds great too, but I need to keep fine control of what I do with the points sampled along the ray, and I’m afraid blending will not provide enough flexibility. But the geometry shader trick is very interesting. I will probably give it a try too, to see how performance compares with the “ray casting loop inside fragment shader” version.