clray

annotate src/scene.cc @ 27:8b2f2ad14ae7

semi-fixed the kdtree construction
author John Tsiombikas <nuclear@member.fsf.org>
date Tue, 17 Aug 2010 20:35:00 +0100
parents c740ae431d51
children 97cfd9675310
rev   line source
John@15 1 #include <math.h>
nuclear@25 2 #include <float.h>
nuclear@26 3 #include <assert.h>
nuclear@22 4 #include "scene.h"
nuclear@27 5 #include "ogl.h"
nuclear@6 6
nuclear@26 7
nuclear@27 8 static void draw_kdtree(const KDNode *node, int level = 0);
nuclear@27 9 static bool build_kdtree(KDNode *kd, int level = 0);
nuclear@27 10 static float eval_cost(const std::list<const Face*> &faces, const AABBox &aabb, int axis, float par_sarea = 1.0);
nuclear@26 11 static void free_kdtree(KDNode *node);
nuclear@27 12 static int kdtree_depth(const KDNode *node);
nuclear@27 13 static void print_item_counts(const KDNode *node, int level);
nuclear@26 14
nuclear@26 15
nuclear@26 16 static int accel_param[NUM_ACCEL_PARAMS] = {
nuclear@26 17 75, // max tree depth
nuclear@26 18 0, // max items per node (0 means ignore limit)
nuclear@26 19 5, // estimated traversal cost
nuclear@26 20 15 // estimated interseciton cost
nuclear@26 21 };
nuclear@26 22
nuclear@26 23
nuclear@26 24 void set_accel_param(int p, int v)
nuclear@26 25 {
nuclear@26 26 assert(p >= 0 && p < NUM_ACCEL_PARAMS);
nuclear@26 27 accel_param[p] = v;
nuclear@26 28 }
nuclear@26 29
nuclear@26 30
John@15 31 #define FEQ(a, b) (fabs((a) - (b)) < 1e-8)
John@15 32 bool Face::operator ==(const Face &f) const
John@15 33 {
John@15 34 for(int i=0; i<3; i++) {
John@15 35 for(int j=0; j<3; j++) {
John@15 36 if(!FEQ(v[i].pos[j], f.v[i].pos[j])) {
John@15 37 return false;
John@15 38 }
John@15 39 if(!FEQ(v[i].normal[j], f.v[i].normal[j])) {
John@15 40 return false;
John@15 41 }
John@15 42 }
John@15 43 if(!FEQ(normal[i], f.normal[i])) {
John@15 44 return false;
John@15 45 }
John@15 46 }
John@15 47 return true;
John@15 48 }
John@15 49
nuclear@25 50 float AABBox::calc_surface_area() const
nuclear@25 51 {
nuclear@25 52 float area1 = (max[0] - min[0]) * (max[1] - min[1]);
nuclear@25 53 float area2 = (max[3] - min[3]) * (max[1] - min[1]);
nuclear@25 54 float area3 = (max[0] - min[0]) * (max[3] - min[3]);
nuclear@25 55
nuclear@25 56 return 2.0f * (area1 + area2 + area3);
nuclear@25 57 }
nuclear@25 58
nuclear@26 59 KDNode::KDNode()
nuclear@26 60 {
nuclear@26 61 axis = 0;
nuclear@26 62 pt = 0.0;
nuclear@26 63 left = right = 0;
nuclear@26 64 num_faces = 0;
nuclear@26 65 }
nuclear@26 66
nuclear@25 67
nuclear@24 68 Scene::Scene()
nuclear@24 69 {
nuclear@24 70 facebuf = 0;
nuclear@24 71 num_faces = -1;
nuclear@24 72 kdtree = 0;
nuclear@24 73 }
nuclear@24 74
nuclear@24 75 Scene::~Scene()
nuclear@24 76 {
nuclear@24 77 delete [] facebuf;
nuclear@24 78 }
nuclear@24 79
nuclear@13 80 bool Scene::add_mesh(Mesh *m)
nuclear@13 81 {
nuclear@13 82 // make sure triangles have material ids
nuclear@13 83 for(size_t i=0; i<m->faces.size(); i++) {
nuclear@13 84 m->faces[i].matid = m->matid;
nuclear@13 85 }
nuclear@24 86
nuclear@24 87 try {
nuclear@24 88 meshes.push_back(m);
nuclear@24 89 }
nuclear@24 90 catch(...) {
nuclear@24 91 return false;
nuclear@24 92 }
nuclear@24 93
nuclear@24 94 // invalidate facebuffer and count
nuclear@24 95 delete [] facebuf;
nuclear@24 96 facebuf = 0;
nuclear@24 97 num_faces = -1;
nuclear@24 98
nuclear@13 99 return true;
nuclear@13 100 }
nuclear@13 101
John@14 102 int Scene::get_num_meshes() const
John@14 103 {
John@14 104 return (int)meshes.size();
John@14 105 }
John@14 106
nuclear@13 107 int Scene::get_num_faces() const
nuclear@13 108 {
nuclear@24 109 if(num_faces >= 0) {
nuclear@24 110 return num_faces;
nuclear@24 111 }
nuclear@24 112
nuclear@24 113 num_faces = 0;
nuclear@13 114 for(size_t i=0; i<meshes.size(); i++) {
nuclear@13 115 num_faces += meshes[i]->faces.size();
nuclear@13 116 }
nuclear@13 117 return num_faces;
nuclear@13 118 }
nuclear@13 119
John@14 120 int Scene::get_num_materials() const
John@14 121 {
John@14 122 return (int)matlib.size();
John@14 123 }
John@14 124
John@14 125 Material *Scene::get_materials()
John@14 126 {
John@14 127 if(matlib.empty()) {
John@14 128 return 0;
John@14 129 }
John@14 130 return &matlib[0];
John@14 131 }
John@14 132
John@14 133 const Material *Scene::get_materials() const
John@14 134 {
John@14 135 if(matlib.empty()) {
John@14 136 return 0;
John@14 137 }
John@14 138 return &matlib[0];
John@14 139 }
nuclear@24 140
nuclear@24 141 const Face *Scene::get_face_buffer() const
nuclear@24 142 {
nuclear@24 143 if(facebuf) {
nuclear@24 144 return facebuf;
nuclear@24 145 }
nuclear@24 146
nuclear@24 147 int num_meshes = get_num_meshes();
nuclear@24 148
nuclear@24 149 printf("constructing face buffer with %d faces (out of %d meshes)\n", get_num_faces(), num_meshes);
nuclear@24 150 facebuf = new Face[num_faces];
nuclear@24 151 Face *fptr = facebuf;
nuclear@24 152
nuclear@24 153 for(int i=0; i<num_meshes; i++) {
nuclear@24 154 for(size_t j=0; j<meshes[i]->faces.size(); j++) {
nuclear@24 155 *fptr++ = meshes[i]->faces[j];
nuclear@24 156 }
nuclear@24 157 }
nuclear@24 158 return facebuf;
nuclear@24 159 }
nuclear@24 160
nuclear@24 161
nuclear@27 162 void Scene::draw_kdtree() const
nuclear@27 163 {
nuclear@27 164 glPushAttrib(GL_ENABLE_BIT);
nuclear@27 165 glDisable(GL_LIGHTING);
nuclear@27 166 glDepthMask(0);
nuclear@27 167
nuclear@27 168 glBegin(GL_LINES);
nuclear@27 169 ::draw_kdtree(kdtree, 0);
nuclear@27 170 glEnd();
nuclear@27 171
nuclear@27 172 glDepthMask(1);
nuclear@27 173 glPopAttrib();
nuclear@27 174 }
nuclear@27 175
nuclear@27 176 static float palette[][3] = {
nuclear@27 177 {0, 1, 0},
nuclear@27 178 {0, 1, 0},
nuclear@27 179 {0, 1, 0},
nuclear@27 180 {1, 0, 0},
nuclear@27 181 {1, 0, 0},
nuclear@27 182 {1, 0, 0},
nuclear@27 183 {0, 0, 1},
nuclear@27 184 {0, 0, 1},
nuclear@27 185 {0, 0, 1},
nuclear@27 186 {1, 1, 0},
nuclear@27 187 {1, 1, 0},
nuclear@27 188 {1, 1, 0},
nuclear@27 189 {0, 0, 1},
nuclear@27 190 {0, 0, 1},
nuclear@27 191 {0, 0, 1},
nuclear@27 192 {1, 0, 1},
nuclear@27 193 {1, 0, 1},
nuclear@27 194 {1, 0, 1}
nuclear@27 195 };
nuclear@27 196 static int pal_size = sizeof palette / sizeof *palette;
nuclear@27 197
nuclear@27 198 static void draw_kdtree(const KDNode *node, int level)
nuclear@27 199 {
nuclear@27 200 if(!node) return;
nuclear@27 201
nuclear@27 202 draw_kdtree(node->left, level + 1);
nuclear@27 203 draw_kdtree(node->right, level + 1);
nuclear@27 204
nuclear@27 205 glColor3fv(palette[level % pal_size]);
nuclear@27 206
nuclear@27 207 glVertex3fv(node->aabb.min);
nuclear@27 208 glVertex3f(node->aabb.max[0], node->aabb.min[1], node->aabb.min[2]);
nuclear@27 209 glVertex3f(node->aabb.max[0], node->aabb.min[1], node->aabb.min[2]);
nuclear@27 210 glVertex3f(node->aabb.max[0], node->aabb.max[1], node->aabb.min[2]);
nuclear@27 211 glVertex3f(node->aabb.max[0], node->aabb.max[1], node->aabb.min[2]);
nuclear@27 212 glVertex3f(node->aabb.min[0], node->aabb.max[1], node->aabb.min[2]);
nuclear@27 213 glVertex3f(node->aabb.min[0], node->aabb.max[1], node->aabb.min[2]);
nuclear@27 214 glVertex3fv(node->aabb.min);
nuclear@27 215
nuclear@27 216 glVertex3f(node->aabb.min[0], node->aabb.min[1], node->aabb.max[2]);
nuclear@27 217 glVertex3f(node->aabb.max[0], node->aabb.min[1], node->aabb.max[2]);
nuclear@27 218 glVertex3f(node->aabb.max[0], node->aabb.min[1], node->aabb.max[2]);
nuclear@27 219 glVertex3fv(node->aabb.max);
nuclear@27 220 glVertex3fv(node->aabb.max);
nuclear@27 221 glVertex3f(node->aabb.min[0], node->aabb.max[1], node->aabb.max[2]);
nuclear@27 222 glVertex3f(node->aabb.min[0], node->aabb.max[1], node->aabb.max[2]);
nuclear@27 223 glVertex3f(node->aabb.min[0], node->aabb.min[1], node->aabb.max[2]);
nuclear@27 224
nuclear@27 225 glVertex3fv(node->aabb.min);
nuclear@27 226 glVertex3f(node->aabb.min[0], node->aabb.min[1], node->aabb.max[2]);
nuclear@27 227 glVertex3f(node->aabb.max[0], node->aabb.min[1], node->aabb.min[2]);
nuclear@27 228 glVertex3f(node->aabb.max[0], node->aabb.min[1], node->aabb.max[2]);
nuclear@27 229 glVertex3f(node->aabb.max[0], node->aabb.max[1], node->aabb.min[2]);
nuclear@27 230 glVertex3fv(node->aabb.max);
nuclear@27 231 glVertex3f(node->aabb.min[0], node->aabb.max[1], node->aabb.min[2]);
nuclear@27 232 glVertex3f(node->aabb.min[0], node->aabb.max[1], node->aabb.max[2]);
nuclear@27 233 }
nuclear@27 234
nuclear@27 235 bool Scene::build_kdtree()
nuclear@24 236 {
nuclear@24 237 const Face *faces = get_face_buffer();
nuclear@24 238 int num_faces = get_num_faces();
nuclear@24 239
nuclear@25 240 printf("Constructing kd-tree out of %d faces ...\n", num_faces);
nuclear@25 241
nuclear@27 242 int icost = accel_param[ACCEL_PARAM_COST_INTERSECT];
nuclear@27 243 int tcost = accel_param[ACCEL_PARAM_COST_TRAVERSE];
nuclear@27 244 printf(" max items per leaf: %d\n", accel_param[ACCEL_PARAM_MAX_NODE_ITEMS]);
nuclear@27 245 printf(" SAH parameters - tcost: %d - icost: %d\n", tcost, icost);
nuclear@27 246
nuclear@25 247 free_kdtree(kdtree);
nuclear@25 248 kdtree = new KDNode;
nuclear@25 249
nuclear@25 250 /* Start the construction of the kdtree by adding all faces of the scene
nuclear@25 251 * to the new root node. At the same time calculate the root's AABB.
nuclear@25 252 */
nuclear@25 253 kdtree->aabb.min[0] = kdtree->aabb.min[1] = kdtree->aabb.min[2] = FLT_MAX;
nuclear@25 254 kdtree->aabb.max[0] = kdtree->aabb.max[1] = kdtree->aabb.max[2] = -FLT_MAX;
nuclear@25 255
nuclear@24 256 for(int i=0; i<num_faces; i++) {
nuclear@25 257 const Face *face = faces + i;
nuclear@25 258
nuclear@25 259 // for each vertex of the face ...
nuclear@25 260 for(int j=0; j<3; j++) {
nuclear@25 261 const float *pos = face->v[j].pos;
nuclear@25 262
nuclear@25 263 // for each element (xyz) of the position vector ...
nuclear@25 264 for(int k=0; k<3; k++) {
nuclear@25 265 if(pos[k] < kdtree->aabb.min[k]) {
nuclear@25 266 kdtree->aabb.min[k] = pos[k];
nuclear@25 267 }
nuclear@25 268 if(pos[k] > kdtree->aabb.max[k]) {
nuclear@25 269 kdtree->aabb.max[k] = pos[k];
nuclear@25 270 }
nuclear@25 271 }
nuclear@25 272 }
nuclear@25 273
nuclear@25 274 kdtree->faces.push_back(face); // add the face
nuclear@26 275 kdtree->num_faces++;
nuclear@24 276 }
nuclear@24 277
nuclear@26 278 // calculate the heuristic for the root
nuclear@26 279 kdtree->cost = eval_cost(kdtree->faces, kdtree->aabb, kdtree->axis);
nuclear@26 280
nuclear@25 281 // now proceed splitting the root recursively
nuclear@27 282 if(!::build_kdtree(kdtree)) {
nuclear@27 283 fprintf(stderr, "failed to build kdtree\n");
nuclear@27 284 return false;
nuclear@27 285 }
nuclear@27 286
nuclear@27 287 printf(" tree depth: %d\n", kdtree_depth(kdtree));
nuclear@27 288 print_item_counts(kdtree, 0);
nuclear@27 289 return true;
nuclear@24 290 }
nuclear@24 291
nuclear@27 292 static bool build_kdtree(KDNode *kd, int level)
nuclear@24 293 {
nuclear@26 294 int opt_max_items = accel_param[ACCEL_PARAM_MAX_NODE_ITEMS];
nuclear@27 295 int tcost = accel_param[ACCEL_PARAM_COST_TRAVERSE];
nuclear@27 296
nuclear@27 297 if(kd->num_faces == 0) {
nuclear@27 298 return true;
nuclear@25 299 }
nuclear@25 300
nuclear@27 301 int axis = level % 3;
nuclear@27 302 //float parent_sa = kd->aabb.calc_surface_area();
nuclear@26 303
nuclear@26 304 float best_cost[2], best_sum_cost = FLT_MAX;
nuclear@26 305 float best_split;
nuclear@26 306
nuclear@26 307 std::list<const Face*>::iterator it = kd->faces.begin();
nuclear@26 308 while(it != kd->faces.end()) {
nuclear@26 309 const Face *face = *it++;
nuclear@26 310
nuclear@26 311 for(int i=0; i<3; i++) {
nuclear@26 312 AABBox aabb_left, aabb_right;
nuclear@26 313 const float *split = face->v[i].pos;
nuclear@26 314
nuclear@26 315 aabb_left = aabb_right = kd->aabb;
nuclear@26 316 aabb_left.max[axis] = split[axis];
nuclear@26 317 aabb_right.min[axis] = split[axis];
nuclear@26 318
nuclear@26 319 float left_cost = eval_cost(kd->faces, aabb_left, axis);
nuclear@26 320 float right_cost = eval_cost(kd->faces, aabb_right, axis);
nuclear@27 321 float sum_cost = left_cost + right_cost - tcost; // tcost is added twice
nuclear@26 322
nuclear@26 323 if(sum_cost < best_sum_cost) {
nuclear@26 324 best_cost[0] = left_cost;
nuclear@26 325 best_cost[1] = right_cost;
nuclear@26 326 best_sum_cost = sum_cost;
nuclear@26 327 best_split = split[axis];
nuclear@26 328 }
nuclear@26 329 }
nuclear@26 330 }
nuclear@26 331
nuclear@27 332 printf("current cost: %f, best_cost: %f\n", kd->cost, best_sum_cost);
nuclear@27 333 if(best_sum_cost > kd->cost && (opt_max_items == 0 || kd->num_faces <= opt_max_items)) {
nuclear@27 334 return true; // stop splitting if it doesn't reduce the cost
nuclear@26 335 }
nuclear@26 336 kd->pt = best_split;
nuclear@26 337
nuclear@26 338 // create the two children
nuclear@26 339 KDNode *kdleft, *kdright;
nuclear@26 340 kdleft = new KDNode;
nuclear@26 341 kdright = new KDNode;
nuclear@26 342
nuclear@26 343 kdleft->aabb = kdright->aabb = kd->aabb;
nuclear@26 344
nuclear@26 345 kdleft->aabb.max[axis] = best_split;
nuclear@26 346 kdright->aabb.min[axis] = best_split;
nuclear@26 347
nuclear@26 348 kdleft->cost = best_cost[0];
nuclear@26 349 kdright->cost = best_cost[1];
nuclear@26 350
nuclear@27 351 //kdleft->axis = kdright->axis = (axis + 1) % 3;
nuclear@27 352
nuclear@26 353 it = kd->faces.begin();
nuclear@26 354 while(it != kd->faces.end()) {
nuclear@26 355 const Face *face = *it++;
nuclear@26 356
nuclear@26 357 if(face->v[0].pos[axis] < best_split ||
nuclear@26 358 face->v[1].pos[axis] < best_split ||
nuclear@26 359 face->v[2].pos[axis] < best_split) {
nuclear@26 360 kdleft->faces.push_back(face);
nuclear@26 361 kdleft->num_faces++;
nuclear@26 362 }
nuclear@26 363 if(face->v[0].pos[axis] >= best_split ||
nuclear@26 364 face->v[1].pos[axis] >= best_split ||
nuclear@26 365 face->v[2].pos[axis] >= best_split) {
nuclear@26 366 kdright->faces.push_back(face);
nuclear@26 367 kdright->num_faces++;
nuclear@26 368 }
nuclear@26 369 }
nuclear@27 370 kd->faces.clear(); // only leaves have faces
nuclear@26 371
nuclear@26 372 kd->left = kdleft;
nuclear@26 373 kd->right = kdright;
nuclear@27 374
nuclear@27 375 return build_kdtree(kd->left, level + 1) && build_kdtree(kd->right, level + 1);
nuclear@26 376 }
nuclear@26 377
nuclear@27 378 static float eval_cost(const std::list<const Face*> &faces, const AABBox &aabb, int axis, float par_sarea)
nuclear@26 379 {
nuclear@26 380 int num_inside = 0;
nuclear@26 381 int tcost = accel_param[ACCEL_PARAM_COST_TRAVERSE];
nuclear@26 382 int icost = accel_param[ACCEL_PARAM_COST_INTERSECT];
nuclear@26 383
nuclear@26 384 std::list<const Face*>::const_iterator it = faces.begin();
nuclear@26 385 while(it != faces.end()) {
nuclear@26 386 const Face *face = *it++;
nuclear@26 387
nuclear@26 388 for(int i=0; i<3; i++) {
nuclear@26 389 if(face->v[i].pos[axis] >= aabb.min[axis] && face->v[i].pos[axis] < aabb.max[axis]) {
nuclear@26 390 num_inside++;
nuclear@26 391 break;
nuclear@26 392 }
nuclear@26 393 }
nuclear@26 394 }
nuclear@26 395
nuclear@27 396 float sarea = aabb.calc_surface_area();
nuclear@27 397 if(sarea < 1e-8) {
nuclear@27 398 return FLT_MAX; // heavily penalize 0-area voxels
nuclear@27 399 }
nuclear@27 400
nuclear@27 401 return tcost + (sarea / par_sarea) * num_inside * icost;
nuclear@24 402 }
nuclear@25 403
nuclear@25 404 static void free_kdtree(KDNode *node)
nuclear@25 405 {
nuclear@25 406 if(node) {
nuclear@25 407 free_kdtree(node->left);
nuclear@25 408 free_kdtree(node->right);
nuclear@25 409 delete node;
nuclear@25 410 }
nuclear@25 411 }
nuclear@27 412
nuclear@27 413 static int kdtree_depth(const KDNode *node)
nuclear@27 414 {
nuclear@27 415 if(!node) return 0;
nuclear@27 416
nuclear@27 417 int left = kdtree_depth(node->left);
nuclear@27 418 int right = kdtree_depth(node->right);
nuclear@27 419 return (left > right ? left : right) + 1;
nuclear@27 420 }
nuclear@27 421
nuclear@27 422 static void print_item_counts(const KDNode *node, int level)
nuclear@27 423 {
nuclear@27 424 if(!node) return;
nuclear@27 425
nuclear@27 426 for(int i=0; i<level; i++) {
nuclear@27 427 fputs(" ", stdout);
nuclear@27 428 }
nuclear@27 429 printf("- %d (cost: %f)\n", node->num_faces, node->cost);
nuclear@27 430
nuclear@27 431 print_item_counts(node->left, level + 1);
nuclear@27 432 print_item_counts(node->right, level + 1);
nuclear@27 433 }