clray: df239a52a091 src/scene.cc

clray

view src/scene.cc @ 55:df239a52a091

extensive render stats for the CPU raytracer

author	John Tsiombikas <nuclear@member.fsf.org>
date	Sat, 11 Sep 2010 03:00:21 +0100
parents	54a96b738afe
children	3d13924b22e6

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"

9 #define CHECK_AABB(aabb) \

10 assert(aabb.max[0] >= aabb.min[0] && aabb.max[1] >= aabb.min[1] && aabb.max[2] >= aabb.min[2])

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

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

17 static int flatten_kdtree(const KDNode *node, KDNodeGPU *kdbuf, int *count);

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

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

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

21 static void free_kdtree(KDNode *node);

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

25 static int accel_param[NUM_ACCEL_PARAMS] = {

26 64, // max tree depth

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

28 5, // estimated traversal cost

29 15 // estimated interseciton cost

30 };

33 void set_accel_param(int p, int v)

34 {

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

36 accel_param[p] = v;

37 }

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

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

42 {

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

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

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

46 return false;

47 }

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

49 return false;

50 }

51 }

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

53 return false;

54 }

55 }

56 return true;

57 }

59 float AABBox::calc_surface_area() const

60 {

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

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

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

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

66 }

68 KDNode::KDNode()

69 {

70 left = right = 0;

71 cost = 0.0;

72 }

75 Scene::Scene()

76 {

77 facebuf = 0;

78 num_faces = -1;

79 kdtree = 0;

80 kdbuf = 0;

81 }

83 Scene::~Scene()

84 {

85 delete [] facebuf;

86 delete [] kdbuf;

87 free_kdtree(kdtree);

88 }

90 bool Scene::add_mesh(Mesh *m)

91 {

92 // make sure triangles have material ids

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

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

95 }

97 try {

98 meshes.push_back(m);

99 }

100 catch(...) {

101 return false;

102 }

103

104 // invalidate facebuffer and count

105 delete [] facebuf;

106 facebuf = 0;

107 num_faces = -1;

108

109 return true;

110 }

111

112 bool Scene::add_light(const Light &lt)

113 {

114 try {

115 lights.push_back(lt);

116 }

117 catch(...) {

118 return false;

119 }

120 return true;

121 }

122

123 int Scene::get_num_meshes() const

124 {

125 return (int)meshes.size();

126 }

127

128 int Scene::get_num_lights() const

129 {

130 return (int)lights.size();

131 }

132

133 int Scene::get_num_faces() const

134 {

135 if(num_faces >= 0) {

136 return num_faces;

137 }

138

139 num_faces = 0;

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

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

142 }

143 return num_faces;

144 }

145

146 int Scene::get_num_materials() const

147 {

148 return (int)matlib.size();

149 }

150

151 int Scene::get_num_kdnodes() const

152 {

153 return kdtree_nodes(kdtree);

154 }

155

156 Mesh **Scene::get_meshes()

157 {

158 if(meshes.empty()) {

159 return 0;

160 }

161 return &meshes[0];

162 }

163

164 const Mesh * const *Scene::get_meshes() const

165 {

166 if(meshes.empty()) {

167 return 0;

168 }

169 return &meshes[0];

170 }

171

172 Light *Scene::get_lights()

173 {

174 if(lights.empty()) {

175 return 0;

176 }

177 return &lights[0];

178 }

179

180 const Light *Scene::get_lights() const

181 {

182 if(lights.empty()) {

183 return 0;

184 }

185 return &lights[0];

186 }

187

188 Material *Scene::get_materials()

189 {

190 if(matlib.empty()) {

191 return 0;

192 }

193 return &matlib[0];

194 }

195

196 const Material *Scene::get_materials() const

197 {

198 if(matlib.empty()) {

199 return 0;

200 }

201 return &matlib[0];

202 }

203

204 const Face *Scene::get_face_buffer() const

205 {

206 if(facebuf) {

207 return facebuf;

208 }

209

210 int num_meshes = get_num_meshes();

211

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

213 facebuf = new Face[num_faces];

214 Face *fptr = facebuf;

215

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

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

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

219 }

220 }

221 return facebuf;

222 }

223

224 const KDNodeGPU *Scene::get_kdtree_buffer() const

225 {

226 if(kdbuf) {

227 return kdbuf;

228 }

229

230 if(!kdtree) {

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

232 }

233

234 int num_nodes = get_num_kdnodes();

235 kdbuf = new KDNodeGPU[num_nodes];

236

237 int count = 0;

238

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

240 flatten_kdtree(kdtree, kdbuf, &count);

241

242 return kdbuf;

243 }

244

245 static int flatten_kdtree(const KDNode *node, KDNodeGPU *kdbuf, int *count)

246 {

247 const size_t max_node_items = sizeof kdbuf[0].face_idx / sizeof kdbuf[0].face_idx[0];

248 int idx = (*count)++;

249

250 // copy the node

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

252 kdbuf[idx].num_faces = 0;

253

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

255 if(i >= max_node_items) {

256 fprintf(stderr, "WARNING too many faces per leaf node!\n");

257 break;

258 }

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

260 kdbuf[idx].num_faces++;

261 }

262

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

264 if(node->left) {

265 assert(node->right);

266

267 kdbuf[idx].left = flatten_kdtree(node->left, kdbuf, count);

268 kdbuf[idx].right = flatten_kdtree(node->right, kdbuf, count);

269 } else {

270 kdbuf[idx].left = kdbuf[idx].right = -1;

271 }

272

273 return idx;

274 }

275

276 void Scene::draw_kdtree() const

277 {

278 glPushAttrib(GL_ENABLE_BIT);

279 glDisable(GL_LIGHTING);

280 glDepthMask(0);

281

282 glBegin(GL_LINES);

283 ::draw_kdtree(kdtree, 0);

284 glEnd();

285

286 glDepthMask(1);

287 glPopAttrib();

288 }

289

290 static float palette[][3] = {

291 {0, 1, 0},

292 {1, 0, 0},

293 {0, 0, 1},

294 {1, 1, 0},

295 {0, 0, 1},

296 {1, 0, 1}

297 };

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

299

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

301 {

302 if(!node) return;

303

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

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

306

307 glColor3fv(palette[level % pal_size]);

308

309 glVertex3fv(node->aabb.min);

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

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

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

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

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

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

316 glVertex3fv(node->aabb.min);

317

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

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

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

321 glVertex3fv(node->aabb.max);

322 glVertex3fv(node->aabb.max);

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

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

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

326

327 glVertex3fv(node->aabb.min);

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

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

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

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

332 glVertex3fv(node->aabb.max);

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

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

335 }

336

337 bool Scene::build_kdtree()

338 {

339 assert(kdtree == 0);

340

341 const Face *faces = get_face_buffer();

342 int num_faces = get_num_faces();

343

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

345

346 int icost = accel_param[ACCEL_PARAM_COST_INTERSECT];

347 int tcost = accel_param[ACCEL_PARAM_COST_TRAVERSE];

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

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

350

351 free_kdtree(kdtree);

352 kdtree = new KDNode;

353

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

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

356 */

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

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

359

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

361 const Face *face = faces + i;

362

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

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

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

366

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

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

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

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

371 }

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

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

374 }

375 }

376 }

377

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

379 }

380

381 CHECK_AABB(kdtree->aabb);

382

383 // calculate the heuristic for the root

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

385

386 // now proceed splitting the root recursively

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

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

389 return false;

390 }

391

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

393 print_item_counts(kdtree, 0);

394 return true;

395 }

396

397 struct Split {

398 int axis;

399 float pos;

400 float sum_cost;

401 float cost_left, cost_right;

402 };

403

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

405 {

406 Split best_split;

407 best_split.sum_cost = FLT_MAX;

408

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

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

411

412 float splitpt[2];

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

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

415

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

417 if(splitpt[j] <= node->aabb.min[axis] || splitpt[j] >= node->aabb.max[axis]) {

418 continue;

419 }

420

421 AABBox aabb_left, aabb_right;

422 aabb_left = aabb_right = node->aabb;

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

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

425

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

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

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

429

430 if(sum_cost < best_split.sum_cost) {

431 best_split.cost_left = left_cost;

432 best_split.cost_right = right_cost;

433 best_split.sum_cost = sum_cost;

434 best_split.pos = splitpt[j];

435 }

436 }

437 }

438

439 assert(split);

440 *split = best_split;

441 split->axis = axis;

442 }

443

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

445 {

446 int opt_max_depth = accel_param[ACCEL_PARAM_MAX_TREE_DEPTH];

447 int opt_max_items = accel_param[ACCEL_PARAM_MAX_NODE_ITEMS];

448

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

450 return true;

451 }

452

453 Split best_split;

454 best_split.axis = -1;

455 best_split.sum_cost = FLT_MAX;

456

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

458 Split split;

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

460

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

462 best_split = split;

463 }

464 }

465

466 if(best_split.axis == -1) {

467 return true; // can't split any more, only 0-area splits available

468 }

469

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

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

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

473 }

474

475 kd->axis = best_split.axis;

476

477 // create the two children

478 KDNode *kdleft, *kdright;

479 kdleft = new KDNode;

480 kdright = new KDNode;

481

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

483

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

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

486

487 kdleft->cost = best_split.cost_left;

488 kdright->cost = best_split.cost_right;

489

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

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

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

493 const Face *face = faces + fidx;

494

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

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

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

498 kdleft->face_idx.push_back(fidx);

499 }

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

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

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

503 kdright->face_idx.push_back(fidx);

504 }

505 }

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

507

508 kd->left = kdleft;

509 kd->right = kdright;

510

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

512 }

513

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

515 {

516 int num_inside = 0;

517 int tcost = accel_param[ACCEL_PARAM_COST_TRAVERSE];

518 int icost = accel_param[ACCEL_PARAM_COST_INTERSECT];

519

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

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

522

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

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

525 num_inside++;

526 break;

527 }

528 }

529 }

530

531 float dx = aabb.max[0] - aabb.min[0];

532 float dy = aabb.max[1] - aabb.min[1];

533 float dz = aabb.max[2] - aabb.min[2];

534

535 if(dx < 0.0) {

536 fprintf(stderr, "FOO DX = %f\n", dx);

537 abort();

538 }

539 if(dy < 0.0) {

540 fprintf(stderr, "FOO DX = %f\n", dy);

541 abort();

542 }

543 if(dz < 0.0) {

544 fprintf(stderr, "FOO DX = %f\n", dz);

545 abort();

546 }

547

548 if(dx < 1e-6 || dy < 1e-6 || dz < 1e-6) {

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

550 }

551

552 float sarea = 2.0 * (dx + dy + dz);//aabb.calc_surface_area();

553 return tcost + sarea * num_inside * icost;

554 }

555

556 static void free_kdtree(KDNode *node)

557 {

558 if(node) {

559 free_kdtree(node->left);

560 free_kdtree(node->right);

561 delete node;

562 }

563 }

564

565 int kdtree_depth(const KDNode *node)

566 {

567 if(!node) return 0;

568

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

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

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

572 }

573

574 int kdtree_nodes(const KDNode *node)

575 {

576 if(!node) return 0;

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

578 }

579

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

581 {

582 if(!node) return;

583

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

585 fputs(" ", stdout);

586 }

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

588

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

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

591 }