Hi,

I was looking for a library that can extract surface from volumetric voxel data. I just found polyvox. I am not sure still it will work in my case.
My Problem:
I have an array of voxel data which I represent as a 3d grid points. Those points represents the density data from CT scan data. This volumetric data represents a static object , though could be rotated (camera). Currently I render the points (900000) using just glVertex3f.
I have another object which is made of a point cloud (5000) and I try to find the collision between objects and I can exactly find which objects are collided. Currently,in the case of collision I inactivate the voxels (not to render the points).
I would like to improve the rendering part by extracting surface from voxel data of the static object and should be able to modify the voxel and the surface.

My major concern are:
- I should be able to render it in real time (30 fps)
- I should also be able to modify the voxel data hence the surface (30fps)

Could you please advice me I could use PolyVox for my work?

My Platform:
Windows XP(SP3)
OpenGl
C++

I am looking forward to some reply from someone .

Thanking you.

Sincerely
Gajananan

Hi,

Well, there is a good chance PolyVox can be useful for you. The library will indeed take volume data and convert it to triangles using the marching cubes algorithm. However, the format of the volumes is different to that produced by CT scanners. In CT data the voxel value corresponds to the material density, where as in PolyVox volumes the voxel value is simply an identifier for a particular material. 0 is always empty spce, and the other values you can define yourself.

Given a CT scan of the human body, you will basically divide the voxels into different materials using one or more thresholds. Anything below -800HU is air in a CT scan, so in the corresponding PolyVox volume you would set this to zero. You could then set all other voxels to '1', or you might want to use more thesholds to seperate bone, fat, etc. You can have up to 255 of these different materials. However, I would start by just dividing into 'air' and 'body'.

PolyVox will then split the volume into a number of cubes (regions), and when you run the surface extraction you will get a mesh for each region. The mesh is basically just a index/vertex list which you can the render using OpenGL or whatever else you like.

So what are the catches? Well, I would imagine your volumes are bigger then the 'castle' which you probably saw in the demo. That castle volume was 256^3, whereas I expect your volumes are at least 512^3. This means more polygons, and I don't know whether this will be a problem. I think it will still work, as most graphics cards are actually limited by the number of 'batches' they can render, rather than the number of polygons.

As for modifying the volume, again it is fine on the smaller volume but I can't be certain about your larger one. There are plenty of optimisations to be done here though...

I should point out that PolyVox does not do any kind of collision detection - it just generates the meshes. In the castle demo the physics is implemented using the 'Bullet' physics engine, but I wouldn't recommend that for you. If you just want collision (rather than full physics) then it's easier and faster to just directly check each of the 5000 points in your cloud against the voxels. Actually, the Bullet physics is one of the slowest parts of my demo, especially when modifying the volume.

Also, there is an extreamly basic OpenGL demo in the PolyVox source code folder. It just generates a volume and then a mesh (no modification or nice materials) but maybe it's a start. You wil need GLUT to run it.

Hi,

Thanking you for your prompt reply!



In my CT scan data, currently I have 85 * 109 * 98 resolution and values of 0, 50,100 ,255 only.



I have a CT scan data of a human tooth. According to your explanation, first I can set
- 0 for Voxel value 0
- 1 for Voxel value 50,100, 255
I am not sure about setting more thresholds.



As I said, currently I have 85 * 109 * 98, which I cut from a high resolution data. I will switch to at least 256 ^ 3. In that case, PolyVox would not have any issue.



Sorry I was not clear in my previous post. Actually I have done the collision detection part using the technique called 3d Occupancy Grid. This is an integer based operation and much faster.



Yes , I have already built PolyVox and the example and run it. Currently I am trying to understand how it is implemented. I am already using glut , no problem then.

Later , I tried to load 64 ^ 3 (croped from my original data) in OpenGl example, but I am confused with how we should set the "volIter.setVoxel"
Does this mean the Iso Value used in MC ?

/Gajananan

Great, well PolyVox requires the volumes to be cubic and to have side lengths which are a power of two. The OpenGL demo uses a 128^3 volume so I guess that will be fine for your purposes.



Yep, that is correct. Latar, you can try setting 1=50, 2=100, and 3=255 if you want to. But start simple



Yep, should be fine.



Ok, well it was originally based on this article by Paul Bourke. I have since created more complex (but hopefully faster) implementation of the algorithm, but the original one is still in there as a reference. Look in sources/PolyVoxImpl/ReferenceSurfaceExtractor.cpp and the main function is called 'extractReferenceSurfaceImpl()'.

This function takes a single region of the volume and generates the corresponding mesh (a mesh is actually called an IndexedSurfacePatch as I recall). So you call that function for each region of your volume, and you can call it again later if the data changes.

The main() function in main.cpp contains a nested loop which calls extractReferenceSurface() for each region in the volume. The meshes are stored, and then rendered later in the display() function.



Yeah... I haven't quite finalised the way the iterators should work. It currently uses a function called 'moveForwardInRegionXYZ()' which is supposed to be fast but is also a little bit confusing. The easiest way is just to use the functions volIter.setPosition() to point the iterator at a particulr voxel, and volIter.setVoxel() to set it's value. Put these in a nested loop which goes from 0-63 in x, y, and z. Does that make sense?

And let me know if you have any problems - it's very likely there are bugs. But it's great that someone else is using the library

I tested it by writing my own OpenGl example then I compared the rendered surface with the image from another surface rendering tool . Results seems fine as I can see the shape of the tooth. Since my example too basic with out light model I am yet to get the desired output.
Meantime I also tested with the following code.I think it works fine.

if(vector.dv > 50)
{
volIter.setVoxel(1);

}
else if(vector.dv == 100)
{
volIter.setVoxel(2);
}
else if(vector.dv == 255)
{
volIter.setVoxel(3);
}
else
{
volIter.setVoxel(0);
}

I still not understood with the idea how Iso value is implemented for MC. Usually we should be able to set different Iso values and see different Iso surfaces. I am not sure Polyvox does this.

/Gajananan

Great, if you look at the SurfaceVertex class you will see that the vertices already contain normals, so it should be easy to add lighting.



The code above is what I had in mind but, thinking about it more carefully, PolyVox won't give quite the right behaviour for your situation. I have to go to work now, but I'll come back tonight and give you more details about how it decides when to generate a surface.

Ok, sorry for the delay...

Right, PolyVox was originally designed for a game environment - the idea was that voxels with a value of 0 would represent empty space where other values would represent various solid materials such as wood, stone, etc. Crucially, it was anticipated that the camera would only ever be in areas of empty space. Consider the diagram below for a voxel world:



As you can see, PolyVox generates a mesh for the white-yellow boundary, and for the white-green boundary, but not for the yellow-green boundary. Even if a mesh was generated for the yellow-green boundary, the user would never be able to see it because it was hidden behing the other one.

In your case you want something slightly different. Because your data came from a density field, you probably never have the case where material 0 (in white) is directly adjacent to material 2 (in green). Your data will probably look more like the following:



It seems that, in your case, you don't want to find the boundaries for 0-1, 0-2, etc. Instead you wan 0-1, 1-2, etc. Do you see what I mean?

Unfortunatly this is slightly beyond what I designed PolyVox to do. However, it shouldn't be too difficult to modify the source code to do what you need. ReferenceSurfaceExtractor.cpp contains the following code:

This basically checks each vertex to determine whether they are empty space (value == 0) or not.If I were you, I would start by playing with the logic in this section.

Hope that helps!

Hi,

Thanking you for your detailed reply. I could not work for last two days due to an illness. I just read your reply.


I understand your point. You picture about my data is correct. But I am not sure If it is enough to have 0-1 boundary only, I need to think about 1-2 when I cut the model and try to go inside.In my model , 0,50,100,255 is the order of density from outside to inside.


I would try this!

/Gajananan

Hi

I tried to load the tooth model [128 * 3]in PolyVox example! this is the rotated view of the second image.




You can compare this with an image that I got after volume rendering of it.


You can also see what I get when I render the sphere:


My question is : why it is showing like flat, I mean no 3D effect in tooth and sphere cases ?

/gajananan

Because the OpenGL demo is very basic, and doesn't apply any lighting or materials to the mesh. If lighting was applied, then parts of the surface which are directly facing the light would be brighter than those which were at an angle, and the result would be a more 3D image.

I forget exactly how you set up lighting in OpenGL but I'm sure there will be articles about it online (you only need a basic approach). As well as setting up the light you will also need to provide OpenGL with the surface normals of the mesh so it can perform the lighting calculations. PolyVox does generate these normals for you, but I'll have to let you know later how to access them as I don't have the code in foront of me right now...

Hope that helps!