diff --git a/include/mitsuba/render/gkdtree.h b/include/mitsuba/render/gkdtree.h index 4142c9ed..9abfedbc 100644 --- a/include/mitsuba/render/gkdtree.h +++ b/include/mitsuba/render/gkdtree.h @@ -512,7 +512,7 @@ public: m_stopPrims = 4; m_maxBadRefines = 3; m_exactPrimThreshold = 65536; - m_maxDepth = 1; + m_maxDepth = 0; m_retract = true; m_parallelBuild = true; } @@ -898,6 +898,13 @@ public: expLeavesVisited /= rootSA; expPrimitivesIntersected /= rootSA; sahCost /= rootSA; + + /* Slightly enlarge the bounding box + (necessary e.g. when the scene is planar) */ + m_aabb.min -= (m_aabb.max-m_aabb.min) * Epsilon + + Vector(Epsilon, Epsilon, Epsilon); + m_aabb.max += (m_aabb.max-m_aabb.min) * Epsilon + + Vector(Epsilon, Epsilon, Epsilon); Log(EDebug, "Detailed kd-tree statistics:"); Log(EDebug, " Inner nodes : %i", ctx.innerNodeCount); @@ -923,42 +930,12 @@ public: /** * \brief Intersect a ray against all primitives stored in the kd-tree */ - bool rayIntersect(const Ray &ray, Intersection &its) { - uint32_t temp[MTS_KD_INTERSECTION_TEMP]; - its.t = std::numeric_limits::infinity(); - Float mint, maxt; - - if (m_aabb.rayIntersect(ray, mint, maxt)) { - if (ray.mint > mint) mint = ray.mint; - if (ray.maxt < maxt) maxt = ray.maxt; - - if (EXPECT_TAKEN(maxt > mint)) { - if (rayIntersect(ray, mint, maxt, its.t, temp)) { - cast()->fillIntersectionDetails(ray, its.t, temp, its); - return true; - } - } - } - return false; - } + bool rayIntersect(const Ray &ray, Intersection &its); /** * \brief Test a ray for intersection against all primitives stored in the kd-tree */ - bool rayIntersect(const Ray &ray) { - uint32_t temp[MTS_KD_INTERSECTION_TEMP]; - Float mint, maxt, t; - - if (m_aabb.rayIntersect(ray, mint, maxt)) { - if (ray.mint > mint) mint = ray.mint; - if (ray.maxt < maxt) maxt = ray.maxt; - - if (EXPECT_TAKEN(maxt > mint)) - if (rayIntersect(ray, mint, maxt, t, temp)) - return true; - } - return false; - } + bool rayIntersect(const Ray &ray); protected: /// Primitive classification during tree-construction @@ -1302,7 +1279,7 @@ protected: * \brief Hashed mailbox implementation */ struct HashedMailbox { - HashedMailbox() { + inline HashedMailbox() { memset(entries, 0xFF, sizeof(index_type)*MTS_KD_MAILBOX_SIZE); } @@ -2455,7 +2432,7 @@ protected: }; /// Ray traversal stack entry for incoherent ray tracing - struct KDStackEntry { + struct KDStackEntryHavran { /* Pointer to the far child */ const KDNode * __restrict node; /* Distance traveled along the ray (entry or exit) */ @@ -2467,16 +2444,17 @@ protected: }; /** - * \brief Internal kd-tree traversal implementation + * \brief Internal kd-tree traversal implementation (Havran variant) */ - template FINLINE bool rayIntersect(const Ray &ray, + template FINLINE bool rayIntersectHavran(const Ray &ray, Float mint, Float maxt, Float &t, void *temp) { - static const int prevAxisTable[] = { 2, 0, 1 }; - static const int nextAxisTable[] = { 1, 2, 0 }; - - KDStackEntry stack[MTS_KD_MAXDEPTH*2]; + KDStackEntryHavran stack[MTS_KD_MAXDEPTH]; const KDNode * __restrict farChild, * __restrict currNode = m_nodes; + #if 0 + static const int prevAxisTable[] = { 2, 0, 1 }; + static const int nextAxisTable[] = { 1, 2, 0 }; + #endif #if defined(MTS_KD_MAILBOX_ENABLED) HashedMailbox mailbox; @@ -2540,14 +2518,18 @@ protected: stack[exPt].prev = tmp; stack[exPt].t = t; stack[exPt].node = farChild; -#if 1 + + #if 1 + /* Faster than the original code with the + prevAxis & nextAxis table */ + stack[exPt].p = ray(t); + #else const int nextAxis = nextAxisTable[axis]; const int prevAxis = prevAxisTable[axis]; stack[exPt].p[nextAxis] = ray.o[nextAxis] + t*ray.d[nextAxis]; stack[exPt].p[prevAxis] = ray.o[prevAxis] + t*ray.d[prevAxis]; -#else - stack[exPt].p = ray(t); -#endif + #endif + stack[exPt].p[axis] = splitVal; } @@ -2560,6 +2542,8 @@ protected: /* Floating-point arithmetic.. - use both absolute and relative epsilons when looking for intersections in the subinterval */ +// const Float searchStart = std::max(mint, stack[enPt].t * m_eps - eps); +// Float searchEnd = std::min(maxt, stack[exPt].t * p_eps + eps); const Float searchStart = std::max(mint, stack[enPt].t * m_eps - eps); Float searchEnd = std::min(maxt, stack[exPt].t * p_eps + eps); @@ -2570,7 +2554,7 @@ protected: const index_type primIdx = m_indices[entry]; #if defined(MTS_KD_MAILBOX_ENABLED) - if (mailbox.get(primIdx)) + if (mailbox.contains(primIdx)) continue; #endif @@ -2585,15 +2569,13 @@ protected: } #if defined(MTS_KD_MAILBOX_ENABLED) - if (result == ENever) + if (result == ENever) mailbox.put(primIdx); #endif } - if (foundIntersection) { - t = searchEnd; + if (foundIntersection) return true; - } /* Pop from the stack and advance to the next node on the interval */ enPt = exPt; @@ -2604,6 +2586,33 @@ protected: return false; } + struct KDStackEntry { + const KDNode * __restrict node; + Float mint, maxt; + }; + + /** + * \brief Internal kd-tree traversal implementation (Havran variant) + */ + template FINLINE bool rayIntersect(const Ray &ray, + Float mint, Float maxt, Float &t, void *temp) { + KDStackEntry stack[MTS_KD_MAXDEPTH]; + int stackPos = 0; + KDNode *node = m_nodes; + + while (node != NULL) { + if (EXPECT_TAKEN(!node->isLeaf())) { + float split = (Float) node->getSplit(); + int axis = node->getAxis(); + float t = (split - ray.o[axis]) * ray.dRcp[axis]; + + const KDNode * __restrict first, * __restrict second; + + } + } + return false; + } + private: KDNode *m_nodes; index_type *m_indices; @@ -2622,6 +2631,44 @@ private: BuildInterface m_interface; }; + +template bool GenericKDTree::rayIntersect( + const Ray &ray, Intersection &its) { + + uint32_t temp[MTS_KD_INTERSECTION_TEMP]; + its.t = std::numeric_limits::infinity(); + Float mint, maxt; + + if (m_aabb.rayIntersect(ray, mint, maxt)) { + if (ray.mint > mint) mint = ray.mint; + if (ray.maxt < maxt) maxt = ray.maxt; + + if (EXPECT_TAKEN(maxt > mint)) { + if (rayIntersectHavran(ray, mint, maxt, its.t, temp)) { + cast()->fillIntersectionDetails(ray, its.t, temp, its); + return true; + } + } + } + return false; +} + +template bool GenericKDTree::rayIntersect( + const Ray &ray) { + uint32_t temp[MTS_KD_INTERSECTION_TEMP]; + Float mint, maxt, t; + + if (m_aabb.rayIntersect(ray, mint, maxt)) { + if (ray.mint > mint) mint = ray.mint; + if (ray.maxt < maxt) maxt = ray.maxt; + + if (EXPECT_TAKEN(maxt > mint)) + if (rayIntersectHavran(ray, mint, maxt, t, temp)) + return true; + } + return false; +} + MTS_NAMESPACE_END #endif /* __KDTREE_GENERIC_H */ diff --git a/include/mitsuba/render/kdtree.h b/include/mitsuba/render/kdtree.h index 3d0aada5..66e57fa6 100644 --- a/include/mitsuba/render/kdtree.h +++ b/include/mitsuba/render/kdtree.h @@ -30,10 +30,10 @@ * whole following block: */ #ifdef MTS_SSE -#define MTS_USE_TRIACCEL4 1 -#define MTS_USE_TRIACCEL 1 +//#define MTS_USE_TRIACCEL4 1 +//#define MTS_USE_TRIACCEL 1 #else -#define MTS_USE_TRIACCEL 1 +//#define MTS_USE_TRIACCEL 1 #endif #if defined(MTS_HAS_COHERENT_RT) diff --git a/src/libcore/triangle.cpp b/src/libcore/triangle.cpp index 3d1ebc5a..6738635b 100644 --- a/src/libcore/triangle.cpp +++ b/src/libcore/triangle.cpp @@ -35,7 +35,7 @@ bool Triangle::rayIntersect(const Vertex *buffer, const Ray &ray, /* if determinant is near zero, ray lies in plane of triangle */ Float det = dot(edge1, pvec); - if (det > -Epsilon && det < Epsilon) + if (det > -1e-8f && det < 1e-8f) return false; Float inv_det = 1.0f / det; diff --git a/src/librender/kdtree_traversal.cpp b/src/librender/kdtree_traversal.cpp index 4d7d35f9..c870d963 100644 --- a/src/librender/kdtree_traversal.cpp +++ b/src/librender/kdtree_traversal.cpp @@ -186,7 +186,7 @@ bool KDTree::rayIntersect(const Ray &ray, Intersection &its, Float mint, Float m Moeller and Trumbore (without having done any pre-computation) */ const KDTriangle &kdTri = m_triangles[m_indices[entry]]; if (EXPECT_TAKEN(kdTri.index != KNoTriangleFlag)) { - const TriMesh *mesh = m_meshes[kdTri.shapeIndex]; + const TriMesh *mesh = (const TriMesh *) m_shapes[kdTri.shapeIndex]; const Triangle &tri = mesh->getTriangles()[kdTri.index]; if (tri.rayIntersect(mesh->getVertexBuffer(), ray, u, v, t)) { if (t >= searchStart && t < searchEnd) { @@ -202,7 +202,7 @@ bool KDTree::rayIntersect(const Ray &ray, Intersection &its, Float mint, Float m } } } else { - if (m_shapes[kdTri.shapeIndex].rayIntersect(ray, searchStart, searchEnd, t)) { + if (m_shapes[kdTri.shapeIndex]->rayIntersect(ray, searchStart, searchEnd, t)) { if (shadowRay) return true; foundIntersection = true; diff --git a/src/tests/test_kd.cpp b/src/tests/test_kd.cpp index 0ec6fa2c..386cf537 100644 --- a/src/tests/test_kd.cpp +++ b/src/tests/test_kd.cpp @@ -107,7 +107,6 @@ public: Float mint, Float maxt, Float &t, void *tmp) { Float tempT, tempU, tempV; if (m_triangles[idx].rayIntersect(m_vertexBuffer, ray, tempU, tempV, tempT)) { - cout << "Got one!" << endl; if (tempT >= mint && tempT <= maxt) { index_type *indexPtr = reinterpret_cast(tmp); Float *floatPtr = reinterpret_cast(indexPtr + 1); @@ -129,17 +128,17 @@ public: Float t, const void *tmp, Intersection &its) const { its.p = ray(t); - const index_type *indexPtr = reinterpret_cast(tmp); - const Float *floatPtr = reinterpret_cast(indexPtr + 1); + //const index_type *indexPtr = reinterpret_cast(tmp); + //const Float *floatPtr = reinterpret_cast(indexPtr + 1); - const Triangle &tri = m_triangles[*indexPtr]; - const Vertex &v0 = m_vertexBuffer[tri.idx[0]]; - const Vertex &v1 = m_vertexBuffer[tri.idx[1]]; - const Vertex &v2 = m_vertexBuffer[tri.idx[2]]; - - const Float u = *floatPtr++, v = *floatPtr++; - const Vector b(1 - u - v, u, v); + //const Triangle &tri = m_triangles[*indexPtr]; + //const Vertex &v0 = m_vertexBuffer[tri.idx[0]]; + //const Vertex &v1 = m_vertexBuffer[tri.idx[1]]; + //const Vertex &v2 = m_vertexBuffer[tri.idx[2]]; + //const Float u = *floatPtr++, v = *floatPtr++; + //const Vector b(1 - u - v, u, v); +/* its.uv = v0.uv * b.x + v1.uv * b.y + v2.uv * b.z; its.dpdu = v0.dpdu * b.x + v1.dpdu * b.y + v2.dpdu * b.z; its.dpdv = v0.dpdv * b.x + v1.dpdv * b.y + v2.dpdv * b.z; @@ -150,6 +149,7 @@ public: * dot(its.shFrame.n, its.dpdu)); its.shFrame.t = cross(its.shFrame.n, its.shFrame.s); its.wi = its.toLocal(-ray.d); +*/ its.hasUVPartials = false; } @@ -172,18 +172,19 @@ public: tree.build(); BSphere bsphere(mesh->getBSphere()); + /* ref oldTree = new KDTree(); oldTree->addShape(mesh); oldTree->build(); + */ for (int j=0; j<3; ++j) { ref random = new Random(); ref timer = new Timer(); - size_t nRays = 100, nIntersections = 0, nIntersectionsBF = 0, nIntersectionsOld = 0; + size_t nRays = 10000000, nIntersections = 0; Log(EInfo, "Bounding sphere: %s", bsphere.toString().c_str()); Log(EInfo, "Shooting " SIZE_T_FMT " rays ..", nRays); - uint32_t tmp[8]; for (size_t i=0; inextFloat(), random->nextFloat()), @@ -193,25 +194,37 @@ public: Ray r(p1, normalize(p2-p1)); Intersection its; - for (uint32_t j=0; jrayIntersect(r, its)) - nIntersectionsOld++; } - Log(EInfo, "KD: Found " SIZE_T_FMT " intersections in %i ms", + Log(EInfo, "New: Found " SIZE_T_FMT " intersections in %i ms", nIntersections, timer->getMilliseconds()); - Log(EInfo, "BF: Found " SIZE_T_FMT " intersections in %i ms", - nIntersectionsBF, timer->getMilliseconds()); - Log(EInfo, "Old: Found " SIZE_T_FMT " intersections in %i ms", - nIntersectionsOld, timer->getMilliseconds()); - Log(EInfo, "%.3f MRays/s", + Log(EInfo, "New: %.3f MRays/s", nRays / (timer->getMilliseconds() * (Float) 1000)); +/* + random = new Random(); + timer->reset(); + nIntersections=0; + + for (size_t i=0; inextFloat(), random->nextFloat()), + sample2(random->nextFloat(), random->nextFloat()); + Point p1 = bsphere.center + squareToSphere(sample1) * bsphere.radius; + Point p2 = bsphere.center + squareToSphere(sample2) * bsphere.radius; + Ray r(p1, normalize(p2-p1)); + Intersection its; + + if (oldTree->rayIntersect(r, its)) + nIntersections++; + } + + Log(EInfo, "Old: Found " SIZE_T_FMT " intersections in %i ms", + nIntersections, timer->getMilliseconds()); + Log(EInfo, "Old: %.3f MRays/s", + nRays / (timer->getMilliseconds() * (Float) 1000)); +*/ } } };