# Ray / Octree traversal: Parametric algorithm implementation

This is my implementation of a Ray/Octree traversal algorithm using an adapted version of the Revelles 2000 paper An efficient parametric algorithm for octree traversal. I’ve made it compliant with IEEE double precision standard. This implementation is minimal and very bare-bones, but should be a good starting point if you want to play with Voxel raycasting.

If you plan on using this code for a project, please include a reference to this website, and my name. The details of the Creative Commons license I use can be found here.

```// Ray-Octree intersection based on
//
// An Efficient Parametric Algorithm for Octree Traversal (2000)
// by J. Revelles , C. Ureña , M. Lastra
//
// Implementation by Jeroen Baert - www.forceflow.be - [email protected]
//
// which can be found at: http://creativecommons.org/licenses/by-nc-sa/3.0/
//
// TL:DR? You must attribute me, cannot use this for commercial purposes,
//
// Thanks,
//
// Jeroen
//

unsigned char a; // because an unsigned char is 8 bits

int first_node(double tx0, double ty0, double tz0, double txm, double tym, double tzm){
unsigned char answer = 0;	// initialize to 00000000
// select the entry plane and set bits
if(tx0 > ty0){
if(tx0 > tz0){ // PLANE YZ
if(tym < tx0) answer|=2;	// set bit at position 1
if(tzm < tx0) answer|=1;	// set bit at position 0 			return (int) answer; 		} 	} 	else { 		if(ty0 > tz0){ // PLANE XZ
if(txm < ty0) answer|=4;	// set bit at position 2
if(tzm < ty0) answer|=1;	// set bit at position 0
}
}
// PLANE XY
if(txm < tz0) answer|=4;	// set bit at position 2
if(tym < tz0) answer|=2;	// set bit at position 1
}

int new_node(double txm, int x, double tym, int y, double tzm, int z){
if(txm < tym){
if(txm < tzm){return x;}  // YZ plane
}
else{
if(tym < tzm){return y;} // XZ plane
}
return z; // XY plane;
}

void proc_subtree (double tx0, double ty0, double tz0, double tx1, double ty1, double tz1, Node* node){
float txm, tym, tzm;
int currNode;

if(tx1 < 0 || ty1 < 0 || tz1 < 0) return; 	if(node->terminal){
cout << "Reached leaf node " << node->debug_ID << endl;
return;
}
else{ cout << "Reached node " << node->debug_ID << endl;} 	txm = 0.5*(tx0 + tx1); 	tym = 0.5*(ty0 + ty1); 	tzm = 0.5*(tz0 + tz1); 	currNode = first_node(tx0,ty0,tz0,txm,tym,tzm); 	do{ 		switch (currNode) 		{ 		case 0: {  			proc_subtree(tx0,ty0,tz0,txm,tym,tzm,node->children[a]);
currNode = new_node(txm,4,tym,2,tzm,1);
break;}
case 1: {
proc_subtree(tx0,ty0,tzm,txm,tym,tz1,node->children[1^a]);
currNode = new_node(txm,5,tym,3,tz1,8);
break;}
case 2: {
proc_subtree(tx0,tym,tz0,txm,ty1,tzm,node->children[2^a]);
currNode = new_node(txm,6,ty1,8,tzm,3);
break;}
case 3: {
proc_subtree(tx0,tym,tzm,txm,ty1,tz1,node->children[3^a]);
currNode = new_node(txm,7,ty1,8,tz1,8);
break;}
case 4: {
proc_subtree(txm,ty0,tz0,tx1,tym,tzm,node->children[4^a]);
currNode = new_node(tx1,8,tym,6,tzm,5);
break;}
case 5: {
proc_subtree(txm,ty0,tzm,tx1,tym,tz1,node->children[5^a]);
currNode = new_node(tx1,8,tym,7,tz1,8);
break;}
case 6: {
proc_subtree(txm,tym,tz0,tx1,ty1,tzm,node->children[6^a]);
currNode = new_node(tx1,8,ty1,8,tzm,7);
break;}
case 7: {
proc_subtree(txm,tym,tzm,tx1,ty1,tz1,node->children[7^a]);
currNode = 8;
break;}
}
} while (currNodesize[0] - ray.origin[0];
ray.direction[0] = - ray.direction[0];
a |= 4 ; //bitwise OR (latest bits are XYZ)
}
if(ray.direction[1] < 0){ 		ray.origin[1] = octree->size[1] - ray.origin[1];
ray.direction[1] = - ray.direction[1];
a |= 2 ;
}
if(ray.direction[2] < 0){ 		ray.origin[2] = octree->size[2] - ray.origin[2];
ray.direction[2] = - ray.direction[2];
a |= 1 ;
}

double divx = 1 / ray.direction[0]; // IEEE stability fix
double divy = 1 / ray.direction[1];
double divz = 1 / ray.direction[2];

double tx0 = (octree->min[0] - ray.origin[0]) * divx;
double tx1 = (octree->max[0] - ray.origin[0]) * divx;
double ty0 = (octree->min[1] - ray.origin[1]) * divy;
double ty1 = (octree->max[1] - ray.origin[1]) * divy;
double tz0 = (octree->min[2] - ray.origin[2]) * divz;
double tz1 = (octree->max[2] - ray.origin[2]) * divz;

if( max(max(tx0,ty0),tz0) < min(min(tx1,ty1),tz1) ){ 		proc_subtree(tx0,ty0,tz0,tx1,ty1,tz1,octree->root);
}
}```

This site uses Akismet to reduce spam. Learn how your comment data is processed.