00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024 #include "unsafekdtreeaccel.h"
00025 #include "paramset.h"
00026
00027 using namespace lux;
00028
00029
00030 UnsafeKdTreeAccel::
00031 UnsafeKdTreeAccel(const vector<Primitive* > &p,
00032 int icost, int tcost,
00033 float ebonus, int maxp, int maxDepth)
00034 : isectCost(icost), traversalCost(tcost),
00035 maxPrims(maxp), emptyBonus(ebonus) {
00036 vector<Primitive* > prims;
00037 for (u_int i = 0; i < p.size(); ++i)
00038 p[i]->FullyRefine(prims);
00039
00040 curMailboxId = 0;
00041 nMailboxes = prims.size();
00042 mailboxPrims = (MailboxPrim *)AllocAligned(nMailboxes *
00043 sizeof(MailboxPrim));
00044 for (u_int i = 0; i < nMailboxes; ++i)
00045 new (&mailboxPrims[i]) MailboxPrim((const Primitive*&)prims[i]);
00046
00047 nextFreeNode = nAllocedNodes = 0;
00048 if (maxDepth <= 0)
00049 maxDepth =
00050 Round2Int(8 + 1.3f * Log2Int(float(prims.size())));
00051
00052 vector<BBox> primBounds;
00053 primBounds.reserve(prims.size());
00054 for (u_int i = 0; i < prims.size(); ++i) {
00055 BBox b = prims[i]->WorldBound();
00056 bounds = Union(bounds, b);
00057 primBounds.push_back(b);
00058 }
00059
00060 UnsafeBoundEdge *edges[3];
00061 for (int i = 0; i < 3; ++i)
00062 edges[i] = new UnsafeBoundEdge[2*prims.size()];
00063 int *prims0 = new int[prims.size()];
00064 int *prims1 = new int[(maxDepth+1) * prims.size()];
00065
00066 int *primNums = new int[prims.size()];
00067 for (u_int i = 0; i < prims.size(); ++i)
00068 primNums[i] = i;
00069
00070 buildTree(0, bounds, primBounds, primNums,
00071 prims.size(), maxDepth, edges,
00072 prims0, prims1);
00073
00074 delete[] primNums;
00075 for (int i = 0; i < 3; ++i)
00076 delete[] edges[i];
00077 delete[] prims0;
00078 delete[] prims1;
00079 }
00080
00081 UnsafeKdTreeAccel::~UnsafeKdTreeAccel() {
00082 for (u_int i = 0; i < nMailboxes; ++i)
00083 mailboxPrims[i].~MailboxPrim();
00084 FreeAligned(mailboxPrims);
00085 FreeAligned(nodes);
00086 }
00087
00088 void UnsafeKdTreeAccel::buildTree(int nodeNum,
00089 const BBox &nodeBounds,
00090 const vector<BBox> &allPrimBounds, int *primNums,
00091 int nPrims, int depth, UnsafeBoundEdge *edges[3],
00092 int *prims0, int *prims1, int badRefines) {
00093 BOOST_ASSERT(nodeNum == nextFreeNode);
00094
00095 if (nextFreeNode == nAllocedNodes) {
00096 int nAlloc = max(2 * nAllocedNodes, 512);
00097 UnsafeKdAccelNode *n = (UnsafeKdAccelNode *)AllocAligned(nAlloc *
00098 sizeof(UnsafeKdAccelNode));
00099 if (nAllocedNodes > 0) {
00100 memcpy(n, nodes,
00101 nAllocedNodes * sizeof(UnsafeKdAccelNode));
00102 FreeAligned(nodes);
00103 }
00104 nodes = n;
00105 nAllocedNodes = nAlloc;
00106 }
00107 ++nextFreeNode;
00108
00109 if (nPrims <= maxPrims || depth == 0) {
00110 nodes[nodeNum].initLeaf(primNums, nPrims,
00111 mailboxPrims, arena);
00112 return;
00113 }
00114
00115
00116 int bestAxis = -1, bestOffset = -1;
00117 float bestCost = INFINITY;
00118 float oldCost = isectCost * float(nPrims);
00119 Vector d = nodeBounds.pMax - nodeBounds.pMin;
00120 float totalSA = (2.f * (d.x*d.y + d.x*d.z + d.y*d.z));
00121 float invTotalSA = 1.f / totalSA;
00122
00123 int axis;
00124 if (d.x > d.y && d.x > d.z) axis = 0;
00125 else axis = (d.y > d.z) ? 1 : 2;
00126 int retries = 0;
00127 retrySplit:
00128
00129 for (int i = 0; i < nPrims; ++i) {
00130 int pn = primNums[i];
00131 const BBox &bbox = allPrimBounds[pn];
00132 edges[axis][2*i] =
00133 UnsafeBoundEdge(bbox.pMin[axis], pn, true);
00134 edges[axis][2*i+1] =
00135 UnsafeBoundEdge(bbox.pMax[axis], pn, false);
00136 }
00137 sort(&edges[axis][0], &edges[axis][2*nPrims]);
00138
00139 int nBelow = 0, nAbove = nPrims;
00140 for (int i = 0; i < 2*nPrims; ++i) {
00141 if (edges[axis][i].type == UnsafeBoundEdge::END) --nAbove;
00142 float edget = edges[axis][i].t;
00143 if (edget > nodeBounds.pMin[axis] &&
00144 edget < nodeBounds.pMax[axis]) {
00145
00146 int otherAxis[3][2] = { {1, 2}, {0, 2}, {0, 1} };
00147 int otherAxis0 = otherAxis[axis][0];
00148 int otherAxis1 = otherAxis[axis][1];
00149 float belowSA = 2 * (d[otherAxis0] * d[otherAxis1] +
00150 (edget - nodeBounds.pMin[axis]) *
00151 (d[otherAxis0] + d[otherAxis1]));
00152 float aboveSA = 2 * (d[otherAxis0] * d[otherAxis1] +
00153 (nodeBounds.pMax[axis] - edget) *
00154 (d[otherAxis0] + d[otherAxis1]));
00155 float pBelow = belowSA * invTotalSA;
00156 float pAbove = aboveSA * invTotalSA;
00157 float eb = (nAbove == 0 || nBelow == 0) ? emptyBonus : 0.f;
00158 float cost = traversalCost + isectCost * (1.f - eb) *
00159 (pBelow * nBelow + pAbove * nAbove);
00160
00161 if (cost < bestCost) {
00162 bestCost = cost;
00163 bestAxis = axis;
00164 bestOffset = i;
00165 }
00166 }
00167 if (edges[axis][i].type == UnsafeBoundEdge::START) ++nBelow;
00168 }
00169 BOOST_ASSERT(nBelow == nPrims && nAbove == 0);
00170
00171 if (bestAxis == -1 && retries < 2) {
00172 ++retries;
00173 axis = (axis+1) % 3;
00174 goto retrySplit;
00175 }
00176 if (bestCost > oldCost) ++badRefines;
00177 if ((bestCost > 4.f * oldCost && nPrims < 16) ||
00178 bestAxis == -1 || badRefines == 3) {
00179 nodes[nodeNum].initLeaf(primNums, nPrims,
00180 mailboxPrims, arena);
00181 return;
00182 }
00183
00184 int n0 = 0, n1 = 0;
00185 for (int i = 0; i < bestOffset; ++i)
00186 if (edges[bestAxis][i].type == UnsafeBoundEdge::START)
00187 prims0[n0++] = edges[bestAxis][i].primNum;
00188 for (int i = bestOffset+1; i < 2*nPrims; ++i)
00189 if (edges[bestAxis][i].type == UnsafeBoundEdge::END)
00190 prims1[n1++] = edges[bestAxis][i].primNum;
00191
00192 float tsplit = edges[bestAxis][bestOffset].t;
00193 nodes[nodeNum].initInterior(bestAxis, tsplit);
00194 BBox bounds0 = nodeBounds, bounds1 = nodeBounds;
00195 bounds0.pMax[bestAxis] = bounds1.pMin[bestAxis] = tsplit;
00196 buildTree(nodeNum+1, bounds0,
00197 allPrimBounds, prims0, n0, depth-1, edges,
00198 prims0, prims1 + nPrims, badRefines);
00199 nodes[nodeNum].aboveChild = nextFreeNode;
00200 buildTree(nodes[nodeNum].aboveChild, bounds1, allPrimBounds,
00201 prims1, n1, depth-1, edges,
00202 prims0, prims1 + nPrims, badRefines);
00203 }
00204
00205 bool UnsafeKdTreeAccel::Intersect(const Ray &ray,
00206 Intersection *isect) const {
00207
00208 float tmin, tmax;
00209 if (!bounds.IntersectP(ray, &tmin, &tmax))
00210 return false;
00211
00212 int rayId = curMailboxId++;
00213 Vector invDir(1.f/ray.d.x, 1.f/ray.d.y, 1.f/ray.d.z);
00214 #define MAX_TODO 64
00215 KdToDo todo[MAX_TODO];
00216 int todoPos = 0;
00217
00218 bool hit = false;
00219 const UnsafeKdAccelNode *node = &nodes[0];
00220 while (node != NULL) {
00221
00222 if (ray.maxt < tmin) break;
00223
00224
00225
00226 if (!node->IsLeaf()) {
00227
00228
00229 int axis = node->SplitAxis();
00230 float tplane = (node->SplitPos() - ray.o[axis]) *
00231 invDir[axis];
00232
00233 const UnsafeKdAccelNode *firstChild, *secondChild;
00234
00235 int belowFirst = (ray.o[axis] < node->SplitPos()) ||
00236 (ray.o[axis] == node->SplitPos() && ray.d[axis] < 0);
00237 if (belowFirst) {
00238 firstChild = node + 1;
00239 secondChild = &nodes[node->aboveChild];
00240 }
00241 else {
00242 firstChild = &nodes[node->aboveChild];
00243 secondChild = node + 1;
00244 }
00245
00246
00247
00248 if (tplane > tmax || tplane <= 0)
00249 node = firstChild;
00250 else if (tplane < tmin)
00251 node = secondChild;
00252 else {
00253
00254 todo[todoPos].node = secondChild;
00255 todo[todoPos].tmin = tplane;
00256 todo[todoPos].tmax = tmax;
00257 ++todoPos;
00258 node = firstChild;
00259 tmax = tplane;
00260 }
00261 }
00262 else {
00263
00264 u_int nPrimitives = node->nPrimitives();
00265 if (nPrimitives == 1) {
00266 MailboxPrim *mp = node->onePrimitive;
00267
00268 if (mp->lastMailboxId != rayId) {
00269 mp->lastMailboxId = rayId;
00270 if (mp->primitive->Intersect(ray, isect))
00271 hit = true;
00272 }
00273 }
00274 else {
00275 MailboxPrim **prims = node->primitives;
00276 for (u_int i = 0; i < nPrimitives; ++i) {
00277 MailboxPrim *mp = prims[i];
00278
00279 if (mp->lastMailboxId != rayId) {
00280 mp->lastMailboxId = rayId;
00281 if (mp->primitive->Intersect(ray, isect))
00282 hit = true;
00283 }
00284 }
00285 }
00286
00287 if (todoPos > 0) {
00288 --todoPos;
00289 node = todo[todoPos].node;
00290 tmin = todo[todoPos].tmin;
00291 tmax = todo[todoPos].tmax;
00292 }
00293 else
00294 break;
00295 }
00296 }
00297 return hit;
00298 }
00299
00300 bool UnsafeKdTreeAccel::IntersectP(const Ray &ray) const {
00301
00302 float tmin, tmax;
00303 if (!bounds.IntersectP(ray, &tmin, &tmax))
00304 return false;
00305
00306 int rayId = curMailboxId++;
00307 Vector invDir(1.f/ray.d.x, 1.f/ray.d.y, 1.f/ray.d.z);
00308 #define MAX_TODO 64
00309 KdToDo todo[MAX_TODO];
00310 int todoPos = 0;
00311 const UnsafeKdAccelNode *node = &nodes[0];
00312 while (node != NULL) {
00313
00314
00315
00316
00317 if (node->IsLeaf()) {
00318
00319 u_int nPrimitives = node->nPrimitives();
00320 if (nPrimitives == 1) {
00321 MailboxPrim *mp = node->onePrimitive;
00322 if (mp->lastMailboxId != rayId) {
00323 mp->lastMailboxId = rayId;
00324 if (mp->primitive->IntersectP(ray))
00325 return true;
00326 }
00327 }
00328 else {
00329 MailboxPrim **prims = node->primitives;
00330 for (u_int i = 0; i < nPrimitives; ++i) {
00331 MailboxPrim *mp = prims[i];
00332 if (mp->lastMailboxId != rayId) {
00333 mp->lastMailboxId = rayId;
00334 if (mp->primitive->IntersectP(ray))
00335 return true;
00336 }
00337 }
00338 }
00339
00340 if (todoPos > 0) {
00341 --todoPos;
00342 node = todo[todoPos].node;
00343 tmin = todo[todoPos].tmin;
00344 tmax = todo[todoPos].tmax;
00345 }
00346 else
00347 break;
00348 }
00349 else {
00350
00351
00352 int axis = node->SplitAxis();
00353 float tplane = (node->SplitPos() - ray.o[axis]) *
00354 invDir[axis];
00355
00356 const UnsafeKdAccelNode *firstChild, *secondChild;
00357
00358 int belowFirst = (ray.o[axis] < node->SplitPos()) ||
00359 (ray.o[axis] == node->SplitPos() && ray.d[axis] < 0);
00360 if (belowFirst) {
00361 firstChild = node + 1;
00362 secondChild = &nodes[node->aboveChild];
00363 }
00364 else {
00365 firstChild = &nodes[node->aboveChild];
00366 secondChild = node + 1;
00367 }
00368
00369
00370
00371 if (tplane > tmax || tplane <= 0)
00372 node = firstChild;
00373 else if (tplane < tmin)
00374 node = secondChild;
00375 else {
00376
00377 todo[todoPos].node = secondChild;
00378 todo[todoPos].tmin = tplane;
00379 todo[todoPos].tmax = tmax;
00380 ++todoPos;
00381 node = firstChild;
00382 tmax = tplane;
00383 }
00384 }
00385 }
00386 return false;
00387 }
00388
00389 Primitive* UnsafeKdTreeAccel::CreateAccelerator(const vector<Primitive* > &prims,
00390 const ParamSet &ps) {
00391 int isectCost = ps.FindOneInt("intersectcost", 80);
00392 int travCost = ps.FindOneInt("traversalcost", 1);
00393 float emptyBonus = ps.FindOneFloat("emptybonus", 0.5f);
00394 int maxPrims = ps.FindOneInt("maxprims", 1);
00395 int maxDepth = ps.FindOneInt("maxdepth", -1);
00396 return new UnsafeKdTreeAccel(prims, isectCost, travCost,
00397 emptyBonus, maxPrims, maxDepth);
00398 }