What I am actually doing is for my kdtree program. I have an object
which I want to divide spatially using kdtrees. There is bounding box
for an entire object. I also built a bounding sphere for every
triangle and then based on the centers of these bounding spheres, I
will calculate a median. At this point the box must be split into left
and right child boxes. Now, the splitting axis which is used to
subdivide a bounding box will depend on depth.
axis = depth % 3
suppose axis ==0 then the box will be split along x axis and the
median center is taken as the split point. If axis == 1 then split
along y axis, if axis == 2 then split along z axis. To cacluate
median, I need to sort the list of bounding spheres within a
particular kd tree node/box. I will sort this list (along x, y ,
z)depending upon value of axis.
My code -
In the below a triangle is
typedef struct triangle_struct
{
int v0, v1, v2; /* 3 indices into vertex list that identify a
triangle */
}triangle;
typedef vector vertex; /* vertex is a vector of doubles */
Object is
typedef struct object_struct
{
int ntri, nvert; /* number of triangles and no of vertices
respectively */
triangle *tri; /* list of triangles */
vertex *vert; /* vertex list */
} object;
kdtree.h
----------------------------------------
#ifndef kd_tree
#define kd_tree
#include "common.h"
#include "reader.h"
typedef struct bounding_sphere_struct
{
vector cen; /*center of bounding sphere */
double r; /* radius */
int tid; /* id of triangle inside this sphere */
triangle tr; /* triangle in this sphere */
}bsphere;
typedef struct kdnode_struct
{
vector maxB, minB; /* maximum bound of the kdnode/box ,
minimum bound*/
vector split; /* point along which this box will be split */
struct kdnode_struct *left, *right; /* left child right child
*/
bsphere* bhead; /* list of bounding spheres inside this box */
int maxspheres; /* Maximum number of spheres in this box */
}kdnode;
int compute_bounding_sphere(bsphere *bs, vector a, vector b, vector
c);
kdnode* init_kdtree(object *, kdnode *);
void build_kdtree(kdnode *kd, int depth);
#endif
kdtree.c
--------------------------------------------
#include "kdtree.h"
int axis = 0;
/* this function runs correctly, as i have checked it */
int compute_bounding_sphere(bsphere *bs, vector a, vector b, vector c)
{
double dotABAB, dotABAC, dotACAC, d,s ,t ;
vector tp, temp1, temp2;
vector_sub(&b, &a, &temp1);
dotABAB = vector_dot( &temp1, &temp1);
vector_sub(&c, &a, &temp2);
dotABAC = vector_dot(&temp1, &temp2);
dotACAC = vector_dot(&temp2, &temp2);
d = 2.0 *(dotABAB*dotACAC - dotABAC*dotABAC);
if(fabs(d) <= 0)
return FAILURE;
s = (dotABAB*dotACAC - dotACAC*dotABAC) / d;
t = (dotACAC*dotABAB - dotABAB*dotABAC) / d;
tp = c;
/* s controls height over AC, t over AB, (1-s-t) over BC */
if (s <= 0.0)
{
bs->cen.x = 0.5*(a.x + c.x);
bs->cen.y = 0.5* (a.y + c.y);
bs->cen.z = 0.5* (a.z + c.z);
}
else
if (t <= 0.0)
{
bs->cen.x = 0.5 *(a.x + b.x);
bs->cen.y = 0.5 * (a.y + b.y);
bs->cen.z = 0.5* (a.z + b.z);
}
else
if (s + t >= 1.0)
{
bs->cen.x = 0.5*(b.x + c.x);
bs->cen.y = 0.5 * (b.y + c.y);
bs->cen.z = 0.5*(b.z + c.z);
}
else
{
bs->cen.x = a.x + s*(b.x - a.x) + t*(c.x - a.x);
bs->cen.y = a.y + s*(b.y - a.y)+ t*(c.y - a.y);
bs->cen.z = a.z + s* (b.z - a.z)+t*(c.z - a.z);
}
vector_sub(&(bs->cen), &tp, &temp1);
bs->r = sqrt( vector_dot(&temp1, &temp1));
return SUCCESS;
}
/* this is incomplete */
void sort_x(kdnode *kd)
{
}
/* this is incomplete */
void sort_y(kdnode *kd)
{
}
/* this is incomplete */
void sort_z(kdnode *kd)
{
}
/* this is incomplete as well but i will make it a recursive function
probably
void build_kdtree(kdnode *kd, int depth)
{
int median;
if(kd->maxspheres <= MINSPHERES || depth < DEPTHLIMIT)
{
kd->left = kd->right = NULL;
return;
}
else
{
if(axis == 0) /* sort kdnode->bhead along x axis, find the
median and use it to split */
{
}
if(axis == 1)/* sort kdnode->bhead along x axis, find the
median and use it to split */
{
}
if(axis == 2)/* sort kdnode->bhead along x axis, find the
median and use it to split */
{
}
build_kdtree(kdnode->left, depth+1);
build_kdtree(kdnode->right, depth+1);
}
/* Initializes the root node of the kdtree */
kdnode* init_kdtree(object *o, kdnode *kdroot)
{
int i,flag;
flag = FAILURE;
kdroot = malloc(sizeof(kdnode));
kdroot->left = kdroot->right = NULL;
kdroot->maxB.x = kdroot->maxB.y = kdroot->maxB.z = -DBL_MAX;
kdroot->minB.x = kdroot->minB.y = kdroot->minB.z = DBL_MAX;
for(i = 0; i<o->nvert; i++)
{
if(o->vert[i].x <kdroot-minB.x)
kdroot->minB.x = o->vert[i].x;
if(o->vert[i].x kdroot->maxB.x)
kdroot->maxB.x = o->vert[i].x;
if(o->vert[i].y < kdroot->minB.y)
kdroot->minB.y = o->vert[i].y;
if(o->vert[i].y kdroot->maxB.y)
kdroot->maxB.y = o->vert[i].y;
if(o->vert[i].z < kdroot->minB.z)
kdroot->minB.z = o->vert[i].z;
if(o->vert[i].z kdroot->maxB.z)
kdroot->maxB.z = o->vert[i].z;
}
/*printf("MAXB: %f %f %f MINB: %f %f %f \n", kdroot-
>maxB.x, kdroot->maxB.y, kdroot->maxB.z, kdroot->minB.x, kdroot-
minB.y, kdroot->minB.z); */
kdroot->bhead = calloc(sizeof(bsphere), o->ntri);
kdroot->maxspheres = o->ntri;
for(i=0; i<o->ntri; i++)
{
flag = compute_bounding_sphere(&(kdroot->bhead[i]), o-
>vert[o->tri[i].v0], o->vert[o->tri[i].v1], o->vert[o->tri[i].v2]);
if(flag == FAILURE)
{
printf("Error in calculation of bounding
spheres\n");
exit(FAILURE);
}
kdroot->bhead[i].tr = o->tri[i];
kdroot->bhead[i].tid = i;
/*printf("Flag:%d Center:%f %f %f radius: %f\n",flag,
kdroot->bhead[i].cen.x, kdroot->bhead[i].cen.y, kdroot-
>bhead[i].cen.z,kdroot->bhead[i].r);*/
}
build_kdtree(kdroot, 0);
return kdroot;
}
