clray: ca445da26588 src/scene.cc

clray

view src/scene.cc @ 37:ca445da26588

out of resources?!

author	John Tsiombikas <nuclear@member.fsf.org>
date	Fri, 27 Aug 2010 01:00:04 +0100
parents	4dec8853bf75
children	4bcf78e572d6

line source

1 #include <stdlib.h>

2 #include <math.h>

3 #include <float.h>

4 #include <assert.h>

5 #include <map>

6 #include "scene.h"

7 #include "ogl.h"

10 #define MIN(a, b) ((a) < (b) ? (a) : (b))

11 #define MAX(a, b) ((a) > (b) ? (a) : (b))

14 static void flatten_kdtree(const KDNode *node, KDNodeGPU *kdbuf, int *count, std::map<const KDNode*, int> *flatmap);

15 static void fix_child_links(const KDNode *node, KDNodeGPU *kdbuf, const std::map<const KDNode*, int> &flatmap);

16 static void draw_kdtree(const KDNode *node, int level = 0);

17 static bool build_kdtree(KDNode *kd, const Face *faces, int level = 0);

18 static float eval_cost(const Face *faces, const int *face_idx, int num_faces, const AABBox &aabb, int axis);

19 static void free_kdtree(KDNode *node);

20 static void print_item_counts(const KDNode *node, int level);

23 static int accel_param[NUM_ACCEL_PARAMS] = {

24 40, // max tree depth

25 0, // max items per node (0 means ignore limit)

26 5, // estimated traversal cost

27 15 // estimated interseciton cost

28 };

31 void set_accel_param(int p, int v)

32 {

33 assert(p >= 0 && p < NUM_ACCEL_PARAMS);

34 accel_param[p] = v;

35 }

38 #define FEQ(a, b) (fabs((a) - (b)) < 1e-8)

39 bool Face::operator ==(const Face &f) const

40 {

41 for(int i=0; i<3; i++) {

42 for(int j=0; j<3; j++) {

43 if(!FEQ(v[i].pos[j], f.v[i].pos[j])) {

44 return false;

45 }

46 if(!FEQ(v[i].normal[j], f.v[i].normal[j])) {

47 return false;

48 }

49 }

50 if(!FEQ(normal[i], f.normal[i])) {

51 return false;

52 }

53 }

54 return true;

55 }

57 float AABBox::calc_surface_area() const

58 {

59 float area1 = (max[0] - min[0]) * (max[1] - min[1]);

60 float area2 = (max[3] - min[3]) * (max[1] - min[1]);

61 float area3 = (max[0] - min[0]) * (max[3] - min[3]);

63 return 2.0f * (area1 + area2 + area3);

64 }

66 KDNode::KDNode()

67 {

68 left = right = 0;

69 cost = 0.0;

70 }

73 Scene::Scene()

74 {

75 facebuf = 0;

76 num_faces = -1;

77 kdtree = 0;

78 kdbuf = 0;

79 }

81 Scene::~Scene()

82 {

83 delete [] facebuf;

84 delete [] kdbuf;

85 free_kdtree(kdtree);

86 }

88 bool Scene::add_mesh(Mesh *m)

89 {

90 // make sure triangles have material ids

91 for(size_t i=0; i<m->faces.size(); i++) {

92 m->faces[i].matid = m->matid;

93 }

95 try {

96 meshes.push_back(m);

97 }

98 catch(...) {

99 return false;

100 }

101

102 // invalidate facebuffer and count

103 delete [] facebuf;

104 facebuf = 0;

105 num_faces = -1;

106

107 return true;

108 }

109

110 int Scene::get_num_meshes() const

111 {

112 return (int)meshes.size();

113 }

114

115 int Scene::get_num_faces() const

116 {

117 if(num_faces >= 0) {

118 return num_faces;

119 }

120

121 num_faces = 0;

122 for(size_t i=0; i<meshes.size(); i++) {

123 num_faces += meshes[i]->faces.size();

124 }

125 return num_faces;

126 }

127

128 int Scene::get_num_materials() const

129 {

130 return (int)matlib.size();

131 }

132

133 int Scene::get_num_kdnodes() const

134 {

135 return kdtree_nodes(kdtree);

136 }

137

138 Material *Scene::get_materials()

139 {

140 if(matlib.empty()) {

141 return 0;

142 }

143 return &matlib[0];

144 }

145

146 const Material *Scene::get_materials() const

147 {

148 if(matlib.empty()) {

149 return 0;

150 }

151 return &matlib[0];

152 }

153

154 const Face *Scene::get_face_buffer() const

155 {

156 if(facebuf) {

157 return facebuf;

158 }

159

160 int num_meshes = get_num_meshes();

161

162 printf("constructing face buffer with %d faces (out of %d meshes)\n", get_num_faces(), num_meshes);

163 facebuf = new Face[num_faces];

164 Face *fptr = facebuf;

165

166 for(int i=0; i<num_meshes; i++) {

167 for(size_t j=0; j<meshes[i]->faces.size(); j++) {

168 *fptr++ = meshes[i]->faces[j];

169 }

170 }

171 return facebuf;

172 }

173

174 static int find_node_index(const KDNode *node, const std::map<const KDNode*, int> &flatmap);

175

176 static bool validate_flat_tree(const KDNode *node, const KDNodeGPU *kdbuf, int count, const std::map<const KDNode*, int> &flatmap)

177 {

178 if(!node) return true;

179

180 int idx = find_node_index(node, flatmap);

181 int left_idx = find_node_index(node->left, flatmap);

182 int right_idx = find_node_index(node->right, flatmap);

183

184 if(kdbuf[idx].left != left_idx) {

185 fprintf(stderr, "%d's left should be %d. it's: %d\n", idx, left_idx, kdbuf[idx].left);

186 return false;

187 }

188 if(kdbuf[idx].right != right_idx) {

189 fprintf(stderr, "%d's right should be %d. it's: %d\n", idx, right_idx, kdbuf[idx].right);

190 return false;

191 }

192 return validate_flat_tree(node->left, kdbuf, count, flatmap) && validate_flat_tree(node->right, kdbuf, count, flatmap);

193 }

194

195 const KDNodeGPU *Scene::get_kdtree_buffer() const

196 {

197 if(kdbuf) {

198 return kdbuf;

199 }

200

201 if(!kdtree) {

202 ((Scene*)this)->build_kdtree();

203 }

204

205 /* we'll associate kdnodes with flattened array indices through this map as

206 * we add them so that the second child-index pass, will be able to find

207 * the children indices of each node.

208 */

209 std::map<const KDNode*, int> flatmap;

210 flatmap[0] = -1;

211

212 int num_nodes = get_num_kdnodes();

213 kdbuf = new KDNodeGPU[num_nodes];

214

215 int count = 0;

216

217 // first arrange the kdnodes into an array (flatten)

218 flatten_kdtree(kdtree, kdbuf, &count, &flatmap);

219

220 // then fix all the left/right links to point to the correct nodes

221 fix_child_links(kdtree, kdbuf, flatmap);

222

223 assert(validate_flat_tree(kdtree, kdbuf, count, flatmap));

224

225 return kdbuf;

226 }

227

228 static void flatten_kdtree(const KDNode *node, KDNodeGPU *kdbuf, int *count, std::map<const KDNode*, int> *flatmap)

229 {

230 int idx = (*count)++;

231

232 // copy the node

233 kdbuf[idx].aabb = node->aabb;

234 for(size_t i=0; i<node->face_idx.size(); i++) {

235 kdbuf[idx].face_idx[i] = node->face_idx[i];

236 }

237 kdbuf[idx].num_faces = node->face_idx.size();

238

239 // update the node* -> array-position mapping

240 (*flatmap)[node] = idx;

241

242 // recurse to the left/right (if we're not in a leaf node)

243 if(node->left) {

244 assert(node->right);

245

246 flatten_kdtree(node->left, kdbuf, count, flatmap);

247 flatten_kdtree(node->right, kdbuf, count, flatmap);

248 }

249 }

250

251 static int find_node_index(const KDNode *node, const std::map<const KDNode*, int> &flatmap)

252 {

253 std::map<const KDNode*, int>::const_iterator it = flatmap.find(node);

254 assert(it != flatmap.end());

255 return it->second;

256 }

257

258 static void fix_child_links(const KDNode *node, KDNodeGPU *kdbuf, const std::map<const KDNode*, int> &flatmap)

259 {

260 if(!node) return;

261

262 int idx = find_node_index(node, flatmap);

263

264 kdbuf[idx].left = find_node_index(node->left, flatmap);

265 if(!node->left && kdbuf[idx].left != -1) abort();

266 kdbuf[idx].right = find_node_index(node->right, flatmap);

267 if(!node->right && kdbuf[idx].right != -1) abort();

268

269 fix_child_links(node->left, kdbuf, flatmap);

270 fix_child_links(node->right, kdbuf, flatmap);

271 }

272

273 void Scene::draw_kdtree() const

274 {

275 glPushAttrib(GL_ENABLE_BIT);

276 glDisable(GL_LIGHTING);

277 glDepthMask(0);

278

279 glBegin(GL_LINES);

280 ::draw_kdtree(kdtree, 0);

281 glEnd();

282

283 glDepthMask(1);

284 glPopAttrib();

285 }

286

287 static float palette[][3] = {

288 {0, 1, 0},

289 {1, 0, 0},

290 {0, 0, 1},

291 {1, 1, 0},

292 {0, 0, 1},

293 {1, 0, 1}

294 };

295 static int pal_size = sizeof palette / sizeof *palette;

296

297 static void draw_kdtree(const KDNode *node, int level)

298 {

299 if(!node) return;

300

301 draw_kdtree(node->left, level + 1);

302 draw_kdtree(node->right, level + 1);

303

304 glColor3fv(palette[level % pal_size]);

305

306 glVertex3fv(node->aabb.min);

307 glVertex3f(node->aabb.max[0], node->aabb.min[1], node->aabb.min[2]);

308 glVertex3f(node->aabb.max[0], node->aabb.min[1], node->aabb.min[2]);

309 glVertex3f(node->aabb.max[0], node->aabb.max[1], node->aabb.min[2]);

310 glVertex3f(node->aabb.max[0], node->aabb.max[1], node->aabb.min[2]);

311 glVertex3f(node->aabb.min[0], node->aabb.max[1], node->aabb.min[2]);

312 glVertex3f(node->aabb.min[0], node->aabb.max[1], node->aabb.min[2]);

313 glVertex3fv(node->aabb.min);

314

315 glVertex3f(node->aabb.min[0], node->aabb.min[1], node->aabb.max[2]);

316 glVertex3f(node->aabb.max[0], node->aabb.min[1], node->aabb.max[2]);

317 glVertex3f(node->aabb.max[0], node->aabb.min[1], node->aabb.max[2]);

318 glVertex3fv(node->aabb.max);

319 glVertex3fv(node->aabb.max);

320 glVertex3f(node->aabb.min[0], node->aabb.max[1], node->aabb.max[2]);

321 glVertex3f(node->aabb.min[0], node->aabb.max[1], node->aabb.max[2]);

322 glVertex3f(node->aabb.min[0], node->aabb.min[1], node->aabb.max[2]);

323

324 glVertex3fv(node->aabb.min);

325 glVertex3f(node->aabb.min[0], node->aabb.min[1], node->aabb.max[2]);

326 glVertex3f(node->aabb.max[0], node->aabb.min[1], node->aabb.min[2]);

327 glVertex3f(node->aabb.max[0], node->aabb.min[1], node->aabb.max[2]);

328 glVertex3f(node->aabb.max[0], node->aabb.max[1], node->aabb.min[2]);

329 glVertex3fv(node->aabb.max);

330 glVertex3f(node->aabb.min[0], node->aabb.max[1], node->aabb.min[2]);

331 glVertex3f(node->aabb.min[0], node->aabb.max[1], node->aabb.max[2]);

332 /*if(!node->left) return;

333

334 AABBox *bleft = &node->left->aabb;

335

336 int axis = level % 3;

337 switch(axis) {

338 case 0:

339 glVertex3f(bleft->max[0], bleft->min[1], bleft->min[2]);

340 glVertex3f(bleft->max[0], bleft->max[1], bleft->min[2]);

341 glVertex3f(bleft->max[0], bleft->max[1], bleft->min[2]);

342 glVertex3f(bleft->max[0], bleft->max[1], bleft->max[2]);

343 glVertex3f(bleft->max[0], bleft->max[1], bleft->max[2]);

344 glVertex3f(bleft->max[0], bleft->min[1], bleft->max[2]);

345 glVertex3f(bleft->max[0], bleft->min[1], bleft->max[2]);

346 glVertex3f(bleft->max[0], bleft->min[1], bleft->min[2]);

347 break;

348

349 case 1:

350 glVertex3f(bleft->min[0], bleft->min[1], bleft->max[2]);

351 glVertex3f(bleft->min[0], bleft->max[1], bleft->max[2]);

352 glVertex3f(bleft->min[0], bleft->max[1], bleft->max[2]);

353 glVertex3f(bleft->max[0], bleft->max[1], bleft->max[2]);

354 glVertex3f(bleft->max[0], bleft->max[1], bleft->max[2]);

355 glVertex3f(bleft->max[0], bleft->min[1], bleft->max[2]);

356 glVertex3f(bleft->max[0], bleft->min[1], bleft->max[2]);

357 glVertex3f(bleft->min[0], bleft->min[1], bleft->max[2]);

358 break;

359

360 case 2:

361 glVertex3f(bleft->min[0], bleft->max[1], bleft->min[2]);

362 glVertex3f(bleft->max[0], bleft->max[1], bleft->min[2]);

363 glVertex3f(bleft->max[0], bleft->max[1], bleft->min[2]);

364 glVertex3f(bleft->max[0], bleft->max[1], bleft->max[2]);

365 glVertex3f(bleft->max[0], bleft->max[1], bleft->max[2]);

366 glVertex3f(bleft->min[0], bleft->max[1], bleft->max[2]);

367 glVertex3f(bleft->min[0], bleft->max[1], bleft->max[2]);

368 glVertex3f(bleft->min[0], bleft->max[1], bleft->min[2]);

369 break;

370

371 default:

372 break;

373 }*/

374 }

375

376 bool Scene::build_kdtree()

377 {

378 assert(kdtree == 0);

379

380 const Face *faces = get_face_buffer();

381 int num_faces = get_num_faces();

382

383 printf("Constructing kd-tree out of %d faces ...\n", num_faces);

384

385 int icost = accel_param[ACCEL_PARAM_COST_INTERSECT];

386 int tcost = accel_param[ACCEL_PARAM_COST_TRAVERSE];

387 printf(" max items per leaf: %d\n", accel_param[ACCEL_PARAM_MAX_NODE_ITEMS]);

388 printf(" SAH parameters - tcost: %d - icost: %d\n", tcost, icost);

389

390 free_kdtree(kdtree);

391 kdtree = new KDNode;

392

393 /* Start the construction of the kdtree by adding all faces of the scene

394 * to the new root node. At the same time calculate the root's AABB.

395 */

396 kdtree->aabb.min[0] = kdtree->aabb.min[1] = kdtree->aabb.min[2] = FLT_MAX;

397 kdtree->aabb.max[0] = kdtree->aabb.max[1] = kdtree->aabb.max[2] = -FLT_MAX;

398

399 for(int i=0; i<num_faces; i++) {

400 const Face *face = faces + i;

401

402 // for each vertex of the face ...

403 for(int j=0; j<3; j++) {

404 const float *pos = face->v[j].pos;

405

406 // for each element (xyz) of the position vector ...

407 for(int k=0; k<3; k++) {

408 if(pos[k] < kdtree->aabb.min[k]) {

409 kdtree->aabb.min[k] = pos[k];

410 }

411 if(pos[k] > kdtree->aabb.max[k]) {

412 kdtree->aabb.max[k] = pos[k];

413 }

414 }

415 }

416

417 kdtree->face_idx.push_back(i); // add the face

418 }

419

420 // calculate the heuristic for the root

421 kdtree->cost = eval_cost(faces, &kdtree->face_idx[0], kdtree->face_idx.size(), kdtree->aabb, 0);

422

423 // now proceed splitting the root recursively

424 if(!::build_kdtree(kdtree, faces)) {

425 fprintf(stderr, "failed to build kdtree\n");

426 return false;

427 }

428

429 printf(" tree depth: %d\n", kdtree_depth(kdtree));

430 print_item_counts(kdtree, 0);

431 return true;

432 }

433

434 struct Split {

435 int axis;

436 float pos;

437 float sum_cost;

438 float cost_left, cost_right;

439 };

440

441 static void find_best_split(const KDNode *node, int axis, const Face *faces, Split *split)

442 {

443 Split best_split;

444 best_split.sum_cost = FLT_MAX;

445

446 for(size_t i=0; i<node->face_idx.size(); i++) {

447 const Face *face = faces + node->face_idx[i];

448

449 float splitpt[2];

450 splitpt[0] = MIN(face->v[0].pos[axis], MIN(face->v[1].pos[axis], face->v[2].pos[axis]));

451 splitpt[1] = MAX(face->v[0].pos[axis], MAX(face->v[1].pos[axis], face->v[2].pos[axis]));

452

453 for(int j=0; j<2; j++) {

454 AABBox aabb_left, aabb_right;

455 aabb_left = aabb_right = node->aabb;

456 aabb_left.max[axis] = splitpt[j];

457 aabb_right.min[axis] = splitpt[j];

458

459 float left_cost = eval_cost(faces, &node->face_idx[0], node->face_idx.size(), aabb_left, axis);

460 float right_cost = eval_cost(faces, &node->face_idx[0], node->face_idx.size(), aabb_right, axis);

461 float sum_cost = left_cost + right_cost - accel_param[ACCEL_PARAM_COST_TRAVERSE]; // tcost is added twice

462

463 if(sum_cost < best_split.sum_cost) {

464 best_split.cost_left = left_cost;

465 best_split.cost_right = right_cost;

466 best_split.sum_cost = sum_cost;

467 best_split.pos = splitpt[j];

468 }

469 }

470 }

471

472 assert(split);

473 *split = best_split;

474 split->axis = axis;

475 }

476

477 static bool build_kdtree(KDNode *kd, const Face *faces, int level)

478 {

479 int opt_max_depth = accel_param[ACCEL_PARAM_MAX_TREE_DEPTH];

480 int opt_max_items = accel_param[ACCEL_PARAM_MAX_NODE_ITEMS];

481

482 if(kd->face_idx.empty() || level >= opt_max_depth) {

483 return true;

484 }

485

486 Split best_split;

487 best_split.sum_cost = FLT_MAX;

488

489 for(int i=0; i<1; i++) {

490 Split split;

491 find_best_split(kd, i, faces, &split);

492

493 if(split.sum_cost < best_split.sum_cost) {

494 best_split = split;

495 }

496 }

497

498 kd->axis = best_split.axis;

499

500 //printf("current cost: %f, best_cost: %f\n", kd->cost, best_sum_cost);

501 if(best_split.sum_cost > kd->cost && (opt_max_items == 0 || (int)kd->face_idx.size() <= opt_max_items)) {

502 return true; // stop splitting if it doesn't reduce the cost

503 }

504

505 // create the two children

506 KDNode *kdleft, *kdright;

507 kdleft = new KDNode;

508 kdright = new KDNode;

509

510 kdleft->aabb = kdright->aabb = kd->aabb;

511

512 kdleft->aabb.max[kd->axis] = best_split.pos;

513 kdright->aabb.min[kd->axis] = best_split.pos;

514

515 kdleft->cost = best_split.cost_left;

516 kdright->cost = best_split.cost_right;

517

518 // TODO would it be much better if we actually split faces that straddle the splitting plane?

519 for(size_t i=0; i<kd->face_idx.size(); i++) {

520 int fidx = kd->face_idx[i];

521 const Face *face = faces + fidx;

522

523 if(face->v[0].pos[kd->axis] < best_split.pos ||

524 face->v[1].pos[kd->axis] < best_split.pos ||

525 face->v[2].pos[kd->axis] < best_split.pos) {

526 kdleft->face_idx.push_back(fidx);

527 }

528 if(face->v[0].pos[kd->axis] >= best_split.pos ||

529 face->v[1].pos[kd->axis] >= best_split.pos ||

530 face->v[2].pos[kd->axis] >= best_split.pos) {

531 kdright->face_idx.push_back(fidx);

532 }

533 }

534 kd->face_idx.clear(); // only leaves have faces

535

536 kd->left = kdleft;

537 kd->right = kdright;

538

539 return build_kdtree(kd->left, faces, level + 1) && build_kdtree(kd->right, faces, level + 1);

540 }

541

542 static float eval_cost(const Face *faces, const int *face_idx, int num_faces, const AABBox &aabb, int axis)

543 {

544 int num_inside = 0;

545 int tcost = accel_param[ACCEL_PARAM_COST_TRAVERSE];

546 int icost = accel_param[ACCEL_PARAM_COST_INTERSECT];

547

548 for(int i=0; i<num_faces; i++) {

549 const Face *face = faces + face_idx[i];

550

551 for(int j=0; j<3; j++) {

552 if(face->v[j].pos[axis] >= aabb.min[axis] && face->v[j].pos[axis] < aabb.max[axis]) {

553 num_inside++;

554 break;

555 }

556 }

557 }

558

559 float sarea = aabb.calc_surface_area();

560 if(sarea < 1e-6) {

561 return FLT_MAX; // heavily penalize 0-area voxels

562 }

563

564 return tcost + sarea * num_inside * icost;

565 }

566

567 static void free_kdtree(KDNode *node)

568 {

569 if(node) {

570 free_kdtree(node->left);

571 free_kdtree(node->right);

572 delete node;

573 }

574 }

575

576 int kdtree_depth(const KDNode *node)

577 {

578 if(!node) return 0;

579

580 int left = kdtree_depth(node->left);

581 int right = kdtree_depth(node->right);

582 return (left > right ? left : right) + 1;

583 }

584

585 int kdtree_nodes(const KDNode *node)

586 {

587 if(!node) return 0;

588 return kdtree_nodes(node->left) + kdtree_nodes(node->right) + 1;

589 }

590

591 static void print_item_counts(const KDNode *node, int level)

592 {

593 if(!node) return;

594

595 for(int i=0; i<level; i++) {

596 fputs(" ", stdout);

597 }

598 printf("- %d (cost: %f)\n", (int)node->face_idx.size(), node->cost);

599

600 print_item_counts(node->left, level + 1);

601 print_item_counts(node->right, level + 1);

602 }