Wednesday, January 5, 2011

Octtree's for Two

Next time you're out having dinner with that special someone, why not strike up a romantic conversation about Octtrees? For just such an occasion, here is what I've figured out about the subject.

Hors d'Ĺ“uvres

An octtree can be thought of as a 3 dimensional binary tree. Typically in a tree the OO mentality would suggest that we create objects for each node, with pointers to the child nodes. But we're using C, and want to avoid all of that hassle. If only there was a clever way to avoid all of that setup. With an array.

First of all, we want to malloc() some memory for this. Each node in an octtree has 8 child nodes, so for an Octtree with 1 layer, there are 9 nodes (1 for the root, and 8 children). For 2 layers there are 73 nodes, and for 3 layers 585. But we're playing with Octtrees, surely there is a trick for figuring out these numbers. If we look at those numbers in Octal, things become clearer. $9_{10}$ in octal is $11_8$. $73_{10}$ in octal is $111_8$, $585_{10}$ in octal is $1111_8$. So we can safely say a $6_{10}$ layer octtree has $1111111_8$ nodes in it. One octal number takes 3 bits to represent, so we can cheat.

int ot_size(int depth)
return 011111111111 & (0xFFFFFFFF >> (32-(depth*3 + 1)));

Code Explanation

Prefixing a zero is how you write in Octal in C. We're trimming down the 0xFFF's down to the length we want, so we can get enough octal 1's for the size octtree we want. (Depth * 3) is because each octal digit is 3 bits long.

Wait a minute...

Why the plus one? Lets take a tree, With a root layer, 1 node, first layer, 8 nodes and second layer, 64 nodes. For calculating the size of the node, it works best to describe this as a 3 layer octtree because we need 3 octal digits to count the number of nodes. Like the good C programmers we are, our layers are indexed starting from 0. So while the root node technically is a layer, it's a very special case. Later when I'm saying a "3 layer tree", you can pretend I'm saying "3 layers plus a root node".

Unfortunately we have imposed a limit on ourselves with this function, of having octtree's with no more than 10 layers in them. But if you really need an Octtree with more than 1,227,133,513 nodes in it, you'll be able to adjust the function yourself.

If we create a 4 layer octtree, and want the offset in the array of the first element of the 4th layer,
will give it to us. Incidentally, the leaf layer contains 8 to the 4th power nodes, and in each dimension there is 2 to the 4th power of nodes. We can prove this as

$8^4 = 2^4 \times 2^4 \times 2^4$
$4096 = 16 \times 16 \times16$

Main Course

Now we have allocated our octtree, we want to start populating it with data. So we need a uniform method of mapping space to a branch in the tree. This borders into Space filling curves, but we can do something easier than that. As mentioned earlier, an octal digit takes 3 bits to store, one for each spatial dimension. We can use Morton Numbers for this. Here is a pair of functions that will do the job

int spread_bits(int x)
x = (x | (x << 16)) & 0x030000FF;
x = (x | (x << 8)) & 0x0300F00F;
x = (x | (x << 4)) & 0x030C30C3;
x = (x | (x << 2)) & 0x09249249;
return x;

int cell_offset(int x, int y, int z)
return spread_bits(x)
| spread_bits(y) << 1
| spread_bits(z) << 2;

So if our octtree is covering 1000.0f units of space (and is an even cube for simlicity), and is 4 layers deep, then each cell on the leaf node layer is $\dfrac{1000}{2^4} = 62.5$ units along each edge. So to find the node into which to insert our point of data (p), we do this

int x,y,z, offset;
x = (int)(p.x / 62.5f);
y = (int)(p.y / 62.5f);
z = (int)(p.z / 62.5f);

offset = cell_offset(x,y,z);

Now before you go and index into your array, don't forget you need to add the offset of the current layer to that.

int layer = ot_size(3);
int target = layer + offset;

Now you're probably wondering why I'm keeping the layer offset and the cell offset in separate variables. Here's why...


As I've been harping on about, octal numbers take 3 bits. And I've rather conspicuously avoided talking about tree traversal. That is because we are going to traverse the tree from bottom to top. So that node we've just populated with data, we want to find it's parent node, couldn't be easier.

layer = layer >> 3;
offset = offset >> 3;
target = layer+offset;

Want to find the neighboring nodes of the leaf we were just looking at?

offset = layer+cell_offset(x+1, y, z);

After Dinner Mint

These techniques apply to n-dimensional binary trees. Want a Quad-Tree? Morton numbers with 2 variables is the canonical example. A Binary tree? Morton numbers with 1 variable are just the same number. Enjoy!

Fluid Simulation Goodness

Just got the fluid simulation code I've been working on to the point where it works correctly. Until I get something more impressive, enjoy this video of a dot vortex!

If you're interested in playing with this yourself, check out the excellent articles by Dr Michael Gourlay over at the Intel Software site.

I should note, this simulation has 27000 particles in it, and runs at 60fps while using about 60% of one of my cores (on a 3.33Ghz c2d). That's without any optimisation, and full debug symbols compiled in.