33 }; |
33 }; |
34 |
34 |
35 class SplitPos |
35 class SplitPos |
36 { |
36 { |
37 public: |
37 public: |
38 float pos; |
38 Float pos; |
39 int lnum, rnum; |
39 int lnum, rnum; |
40 SplitPos(): pos(0.0), lnum(0), rnum(0) {}; |
40 SplitPos(): pos(0.0), lnum(0), rnum(0) {}; |
41 SplitPos(float &aPos): pos(aPos), lnum(0), rnum(0) {}; |
41 SplitPos(Float &aPos): pos(aPos), lnum(0), rnum(0) {}; |
42 friend bool operator<(const SplitPos& a, const SplitPos& b) |
42 friend bool operator<(const SplitPos& a, const SplitPos& b) |
43 { return a.pos < b.pos; }; |
43 { return a.pos < b.pos; }; |
44 friend bool operator==(const SplitPos& a, const SplitPos& b) |
44 friend bool operator==(const SplitPos& a, const SplitPos& b) |
45 { return a.pos == b.pos; }; |
45 { return a.pos == b.pos; }; |
46 }; |
46 }; |
52 // stack element for kd-tree traversal |
52 // stack element for kd-tree traversal |
53 class StackElem |
53 class StackElem |
54 { |
54 { |
55 public: |
55 public: |
56 KdNode* node; /* pointer to far child */ |
56 KdNode* node; /* pointer to far child */ |
57 float t; /* the entry/exit signed distance */ |
57 Float t; /* the entry/exit signed distance */ |
58 Vector3 pb; /* the coordinates of entry/exit point */ |
58 Vector3 pb; /* the coordinates of entry/exit point */ |
59 StackElem(KdNode *anode, const float &at, const Vector3 &apb): |
59 StackElem(KdNode *anode, const Float &at, const Vector3 &apb): |
60 node(anode), t(at), pb(apb) {}; |
60 node(anode), t(at), pb(apb) {}; |
61 }; |
61 }; |
62 |
62 |
63 // ---------------------------------------- |
63 // ---------------------------------------- |
64 |
64 |
79 if (shapebb.H.z > bbox.H.z) bbox.H.z = shapebb.H.z; |
79 if (shapebb.H.z > bbox.H.z) bbox.H.z = shapebb.H.z; |
80 } |
80 } |
81 }; |
81 }; |
82 |
82 |
83 Shape *Container::nearest_intersection(const Shape *origin_shape, const Ray &ray, |
83 Shape *Container::nearest_intersection(const Shape *origin_shape, const Ray &ray, |
84 float &nearest_distance) |
84 Float &nearest_distance) |
85 { |
85 { |
86 Shape *nearest_shape = NULL; |
86 Shape *nearest_shape = NULL; |
87 ShapeList::iterator shape; |
87 ShapeList::iterator shape; |
88 for (shape = shapes.begin(); shape != shapes.end(); shape++) |
88 for (shape = shapes.begin(); shape != shapes.end(); shape++) |
89 if (*shape != origin_shape && (*shape)->intersect(ray, nearest_distance)) |
89 if (*shape != origin_shape && (*shape)->intersect(ray, nearest_distance)) |
181 } |
181 } |
182 } |
182 } |
183 } |
183 } |
184 |
184 |
185 // choose best split pos |
185 // choose best split pos |
186 const float K = 1.4; // constant, K = cost of traversal / cost of ray-triangle intersection |
186 const Float K = 1.4; // constant, K = cost of traversal / cost of ray-triangle intersection |
187 float SAV = 2*(bbox.w()*bbox.h() + bbox.w()*bbox.d() + bbox.h()*bbox.d()); // surface area of node |
187 Float SAV = 2*(bbox.w()*bbox.h() + bbox.w()*bbox.d() + bbox.h()*bbox.d()); // surface area of node |
188 float cost = SAV * (K + shapes->size()); // initial cost = non-split cost |
188 Float cost = SAV * (K + shapes->size()); // initial cost = non-split cost |
189 SplitPos *splitpos = NULL; |
189 SplitPos *splitpos = NULL; |
190 bool leaf = true; |
190 bool leaf = true; |
191 BBox lbb = bbox; |
191 BBox lbb = bbox; |
192 BBox rbb = bbox; |
192 BBox rbb = bbox; |
193 for (spl = splitlist.begin(); spl != splitlist.end(); spl++) |
193 for (spl = splitlist.begin(); spl != splitlist.end(); spl++) |
194 { |
194 { |
195 // calculate SAH cost of this split |
195 // calculate SAH cost of this split |
196 lbb.H.cell[axis] = spl->pos; |
196 lbb.H.cell[axis] = spl->pos; |
197 rbb.L.cell[axis] = spl->pos; |
197 rbb.L.cell[axis] = spl->pos; |
198 float SAL = 2*(lbb.w()*lbb.h() + lbb.w()*lbb.d() + lbb.h()*lbb.d()); |
198 Float SAL = 2*(lbb.w()*lbb.h() + lbb.w()*lbb.d() + lbb.h()*lbb.d()); |
199 float SAR = 2*(rbb.w()*rbb.h() + rbb.w()*rbb.d() + rbb.h()*rbb.d()); |
199 Float SAR = 2*(rbb.w()*rbb.h() + rbb.w()*rbb.d() + rbb.h()*rbb.d()); |
200 float splitcost = K + SAL/SAV*(K+spl->lnum) + SAR/SAV*(K+spl->rnum); |
200 Float splitcost = K + SAL/SAV*(K+spl->lnum) + SAR/SAV*(K+spl->rnum); |
201 |
201 |
202 if (splitcost < cost) |
202 if (splitcost < cost) |
203 { |
203 { |
204 leaf = false; |
204 leaf = false; |
205 cost = splitcost; |
205 cost = splitcost; |
235 F << "f " << f+1 << " " << f+2 << " " << f+3 << " " << f+4 << endl; |
235 F << "f " << f+1 << " " << f+2 << " " << f+3 << " " << f+4 << endl; |
236 f += 4; |
236 f += 4; |
237 #endif |
237 #endif |
238 |
238 |
239 split = splitpos->pos; |
239 split = splitpos->pos; |
240 float lnum = splitpos->lnum; |
240 Float lnum = splitpos->lnum; |
241 float rnum = splitpos->rnum; |
241 Float rnum = splitpos->rnum; |
242 |
242 |
243 // split this node |
243 // split this node |
244 delete shapes; |
244 delete shapes; |
245 children = new KdNode[2]; |
245 children = new KdNode[2]; |
246 int state = 0; |
246 int state = 0; |
313 built = true; |
313 built = true; |
314 } |
314 } |
315 |
315 |
316 /* algorithm by Vlastimil Havran, Heuristic Ray Shooting Algorithms, appendix C */ |
316 /* algorithm by Vlastimil Havran, Heuristic Ray Shooting Algorithms, appendix C */ |
317 Shape *KdTree::nearest_intersection(const Shape *origin_shape, const Ray &ray, |
317 Shape *KdTree::nearest_intersection(const Shape *origin_shape, const Ray &ray, |
318 float &nearest_distance) |
318 Float &nearest_distance) |
319 { |
319 { |
320 float a, b; /* entry/exit signed distance */ |
320 Float a, b; /* entry/exit signed distance */ |
321 float t; /* signed distance to the splitting plane */ |
321 Float t; /* signed distance to the splitting plane */ |
322 |
322 |
323 /* if we have no tree, fall back to naive test */ |
323 /* if we have no tree, fall back to naive test */ |
324 if (!built) |
324 if (!built) |
325 return Container::nearest_intersection(origin_shape, ray, nearest_distance); |
325 return Container::nearest_intersection(origin_shape, ray, nearest_distance); |
326 |
326 |
350 exPt = st.back(); |
350 exPt = st.back(); |
351 /* loop until a leaf is found */ |
351 /* loop until a leaf is found */ |
352 while (!node->isLeaf()) |
352 while (!node->isLeaf()) |
353 { |
353 { |
354 /* retrieve position of splitting plane */ |
354 /* retrieve position of splitting plane */ |
355 float splitVal = node->getSplit(); |
355 Float splitVal = node->getSplit(); |
356 short axis = node->getAxis(); |
356 short axis = node->getAxis(); |
357 |
357 |
358 if (enPt->pb[axis] <= splitVal) |
358 if (enPt->pb[axis] <= splitVal) |
359 { |
359 { |
360 if (exPt->pb[axis] <= splitVal) |
360 if (exPt->pb[axis] <= splitVal) |
399 |
399 |
400 /* current node is the leaf . . . empty or full */ |
400 /* current node is the leaf . . . empty or full */ |
401 /* "intersect ray with each object in the object list, discarding " |
401 /* "intersect ray with each object in the object list, discarding " |
402 "those lying before stack[enPt].t or farther than stack[exPt].t" */ |
402 "those lying before stack[enPt].t or farther than stack[exPt].t" */ |
403 Shape *nearest_shape = NULL; |
403 Shape *nearest_shape = NULL; |
404 float dist = exPt->t; |
404 Float dist = exPt->t; |
405 ShapeList::iterator shape; |
405 ShapeList::iterator shape; |
406 for (shape = node->shapes->begin(); shape != node->shapes->end(); shape++) |
406 for (shape = node->shapes->begin(); shape != node->shapes->end(); shape++) |
407 if (*shape != origin_shape && (*shape)->intersect(ray, dist) |
407 if (*shape != origin_shape && (*shape)->intersect(ray, dist) |
408 && dist >= enPt->t) |
408 && dist >= enPt->t) |
409 nearest_shape = *shape; |
409 nearest_shape = *shape; |