4

In regards of Marching Cubes, i have some doubts regarding its algorithm and implementation. I have gone through the excellent Paul Bourke article of Marching Cubes and also the available source code on the site, yet, i still encountered some problems in term of understanding as well as how to implement the algo in my own way. The questions are as below:

  • Gridcell size - I have read that gridcell size affects the quality of the produced 3D model. For example if i have a stack of xray image set with size of (200*200*200), so, a slab of gridcells will be constructed from 2 adjacent slices of image. Thus, the total of the gridcells in a slab would be (200-1)*(200-1) with each of the gridcell corner correspond to the pixel value/density of the image. Is this correct?? Besides, how do we implement different size for gridcell??

  • Voxel size - I have read the few references of Marching Cubes and i cant seem to find how voxel size is taken care of in the algorithm. Please correct me if i am wrong, in my case, the gap between the adjacent layer of images is 1 mil in size; thus, how do i take care of those in Marching Cubes algo or is it a dead end?? Is it taken care of as the size of Gridcell?? (Assumption: the size of one pixel in its xy coordinate is 19 micron while the gap/z is 25.4 micron/1 mil in length)

  • Coordinates of the Gridcell corner(Vertices Coordinates of the cube) - I am trying to assign the coordinates of the corners of gridcells with index i j k by nested looping of the image set dimension (200*200*200). Is this correct?? Is there any faster way to do it??

Note: i have seen implementation of MC in VTK but it is quite hard for me to digest it as it depends on some other VTK classes.

vincent911001
  • 523
  • 1
  • 6
  • 20

1 Answers1

5

Lots of questions. I am going to try to give some pointers. First of 200^3 is a pretty small dataset for ct! what about 1024^3? :)

Marching cubes is built for regular grids. So wether data is defined at cube vertices or centers really does not matter: just shift by half the cube size! If you have irregular data use something else or resample to a regular grid first.

You also seem to be missing the "marching" part: The idea is to find one cube with the surface and flood fill out from there. Cubes that are all outside or all inside stop the search. That way most cubes in your huge regular grid do not even need to be looked at.

Scaling to real units should be a last step. Think of the input volume as normalized to 1x1x1. Then scale the output vertices to physical units. The data you have is the data you have. Any resampling should be done before during reconstruction or filtering. It has no place in the geometry stage.

I am not sure I understand the last question, but one thing that is really important for further processing is to create a connected, indexed mesh. One important trick there is to just keep a kind of hash table of the previous slice/line/neighbor. So you can quickly look up already created vertices and reuse their index. The result should be a connected mesh with unique vertices. This you can then use in any kind of geometry processing.

starmole
  • 4,974
  • 1
  • 28
  • 48
  • Hi, ya, the image size given is pretty small, but in reality, i would have an image set with size for about 2000*2000*201 something. Regarding the gridcell, is the eight corners of the gridcell comprised of the pixel values from adjacent slices(2 pixels from n slice + 2 pixels from n+1 slice)?? – vincent911001 Jun 02 '15 at 05:39
  • really thanks for your advice, it really does clear some of my doubts. Thanks a lot – vincent911001 Jun 02 '15 at 05:44
  • @vincent911001: never mind slice size. your cube is always x,y,z - x+1,y+1+z+1. So 4 values from slice z at(x,y)-(x+1,y+1) and 4 values from slice z+1 at (x,y)-(x+1,y+1) is your cube. Any scaling should be post or pre processing. – starmole Jun 02 '15 at 05:51
  • Hi, starmole, thanks for your clear description, it does really help. BTW, is the vertices of the gridcell created using the looping counter index(For example: i, j, k that is used to loop through whole volume)?? – vincent911001 Jun 02 '15 at 05:54
  • One simple way is whenever you create a vertex in cube (x,y,z) plus (1/0,1/0,1/0): look it up in a hash table indexed by uniquekey(x,y,z,0/1,0/1,0/1). If it is not there create it, otherwise use the index. Now the first trick is to realize, that for every x,y,z you only need 3 possible vertices: down x, down y, down z. The second trick is to know that when you go through tour grid top down, left right, front bottom, that table is very small. You just check if ( havevertex(x,y,z,down/left/back) ) use it else addit(). And the hasvertex is just a hash table. – starmole Jun 02 '15 at 06:18
  • Thanks a lot starmole, if i understand your advice correctly, it is to put or store the interpolated vertices generated from MC into the hashtable in order to be able to reuse it or export it into some kind of geometry file format, is it correct?? Before i get the output of vertices or facet, i would need 8 corners of vertices to be processed in MC in order to determine which edges is intersected and after that apply interpolation to get the vertices, right?? what i would like to know is how to produce the initial vertices of the gridcells?? – vincent911001 Jun 02 '15 at 06:37
  • what i mean is if you create vertices for the first cube, the neighboring cubes should reuse them (indexed mesh). Do one cube at a time. But before making a new vertex at an edge (x,y,z) - (x+dx, y+dy, z+dz) look it up in the hash first - if it is there output that index. If not add it. That way you end up with a connected mesh that is easy to process. – starmole Jun 02 '15 at 06:50
  • Hi, starmole, i guess i know what are you trying to tell me, thanks a lot for your suggestion btw. – vincent911001 Jun 02 '15 at 06:54
  • NP. Marching cubes is one algorithm I really love. :) So I love talking about it! Good luck with it! And if you have more questions feel free to ask. – starmole Jun 02 '15 at 07:02
  • Sure, i will keep you posted if i bump into any problem while implementing it. Thanks a lot btw. – vincent911001 Jun 02 '15 at 07:35
  • Let us [continue this discussion in chat](http://chat.stackoverflow.com/rooms/79469/discussion-between-vincent911001-and-starmole). – vincent911001 Jun 03 '15 at 00:27