diff --git a/SConstruct b/SConstruct index 8d52be54..1e3eaed3 100644 --- a/SConstruct +++ b/SConstruct @@ -517,7 +517,7 @@ plugins += env.SharedLibrary('plugins/obj', ['src/shapes/obj.cpp']) plugins += env.SharedLibrary('plugins/ply', ['src/shapes/ply/ply.cpp', 'src/shapes/ply/ply_parser.cpp'], CPPPATH = env['CPPPATH'] + ['src/shapes/ply']) plugins += env.SharedLibrary('plugins/serialized', ['src/shapes/serialized.cpp']) -#plugins += env.SharedLibrary('plugins/sphere', ['src/shapes/sphere.cpp']) +plugins += env.SharedLibrary('plugins/sphere', ['src/shapes/sphere.cpp']) #plugins += env.SharedLibrary('plugins/cylinder', ['src/shapes/cylinder.cpp']) #plugins += env.SharedLibrary('plugins/hair', ['src/shapes/hair.cpp', 'src/shapes/miterseg.cpp']) #plugins += env.SharedLibrary('plugins/group', ['src/shapes/group.cpp']) diff --git a/doc/integrator.tex b/doc/integrator.tex index e1ace30f..cf77b514 100644 --- a/doc/integrator.tex +++ b/doc/integrator.tex @@ -200,7 +200,7 @@ To avoid having to do this every time \code{Li()} is called, we can override the \code{preprocess} function: \begin{cpp} /// Preprocess function -- called on the initiating machine - void preprocess(const Scene *scene, RenderQueue *queue, + bool preprocess(const Scene *scene, RenderQueue *queue, const RenderJob *job, int sceneResID, int cameraResID, int samplerResID) { SampleIntegrator::preprocess(scene, queue, job, sceneResID, @@ -213,6 +213,8 @@ we can override the \code{preprocess} function: for (int i=0; i<8; ++i) m_maxDist = std::max(m_maxDist, (cameraPosition - sceneAABB.getCorner(i)).length()); + + return true; } \end{cpp} The bottom of this function should be relatively self-explanatory. The diff --git a/include/mitsuba/core/ray.h b/include/mitsuba/core/ray.h index 7e519fc9..3a5e270c 100644 --- a/include/mitsuba/core/ray.h +++ b/include/mitsuba/core/ray.h @@ -185,7 +185,7 @@ struct RayInterval4 { mint = SSEConstants::eps; maxt = SSEConstants::p_inf; } - + inline RayInterval4(const Ray *rays) { for (int i=0; i<4; i++) { mint.f[i] = rays[i].mint; diff --git a/include/mitsuba/render/gkdtree.h b/include/mitsuba/render/gkdtree.h index f52ed8ba..01fe8e94 100644 --- a/include/mitsuba/render/gkdtree.h +++ b/include/mitsuba/render/gkdtree.h @@ -40,8 +40,8 @@ #define MTS_KD_BLOCKSIZE_KD (512*1024/sizeof(KDNode)) #define MTS_KD_BLOCKSIZE_IDX (512*1024/sizeof(uint32_t)) -/// 64 byte temporary storage for intersection computations -#define MTS_KD_INTERSECTION_TEMP 64 +/// 32 byte temporary storage for intersection computations +#define MTS_KD_INTERSECTION_TEMP 32 /// Use a simple hashed 8-entry mailbox per thread //#define MTS_KD_MAILBOX_ENABLED 1 diff --git a/include/mitsuba/render/kdtree.h b/include/mitsuba/render/kdtree.h index 1edd097d..79dab483 100644 --- a/include/mitsuba/render/kdtree.h +++ b/include/mitsuba/render/kdtree.h @@ -90,14 +90,14 @@ public: * use of ray coherence to do this very efficiently. Requires SSE. */ void rayIntersectPacket(const RayPacket4 &packet, - const RayInterval4 &interval, Intersection4 &its) const; + const RayInterval4 &interval, Intersection4 &its, void *temp) const; /** * \brief Fallback for incoherent rays * \sa rayIntesectPacket */ void rayIntersectPacketIncoherent(const RayPacket4 &packet, - const RayInterval4 &interval, Intersection4 &its) const; + const RayInterval4 &interval, Intersection4 &its, void *temp) const; #endif MTS_DECLARE_CLASS() @@ -135,7 +135,7 @@ protected: const TriMesh *mesh = static_cast(shape); return mesh->getTriangles()[idx].getClippedAABB(mesh->getVertexPositions(), aabb); } else { - return shape->getAABB(); + return shape->getClippedAABB(aabb); } } @@ -146,7 +146,7 @@ protected: /// Temporarily holds some intersection information struct IntersectionCache { size_type shapeIndex; - size_type index; + size_type primIndex; Float u, v; }; @@ -173,7 +173,7 @@ protected: return ENo; t = tempT; cache->shapeIndex = shapeIdx; - cache->index = idx; + cache->primIndex = primIdx; cache->u = tempU; cache->v = tempV; return EYes; @@ -181,8 +181,9 @@ protected: } else { const Shape *shape = m_shapes[shapeIndex]; if (shape->rayIntersect(ray, mint, maxt, t, - reinterpret_cast(temp) + 4)) { + reinterpret_cast(temp) + 8)) { cache->shapeIndex = shapeIdx; + cache->primIndex = KNoTriangleFlag; return EYes; } } @@ -193,7 +194,7 @@ protected: if (ta.rayIntersect(ray, mint, maxt, tempU, tempV, tempT)) { t = tempT; cache->shapeIndex = ta.shapeIndex; - cache->index = ta.index; + cache->primIndex = ta.primIndex; cache->u = tempU; cache->v = tempV; return EYes; @@ -202,8 +203,9 @@ protected: uint32_t shapeIndex = ta.shapeIndex; const Shape *shape = m_shapes[shapeIndex]; if (shape->rayIntersect(ray, mint, maxt, t, - reinterpret_cast(temp) + 4)) { + reinterpret_cast(temp) + 8)) { cache->shapeIndex = shapeIndex; + cache->primIndex = KNoTriangleFlag; return EYes; } } @@ -221,17 +223,6 @@ protected: }; #endif - /** - * Fallback for incoherent rays - */ - inline void rayIntersectPacketIncoherent(const Ray *rays, - Intersection *its) const { - for (int i=0; i<4; i++) { - if (!rayIntersect(rays[i], its[i])) - its[i].t = std::numeric_limits::infinity(); - } - } - /// Virtual destructor virtual ~KDTree(); private: diff --git a/include/mitsuba/render/shape.h b/include/mitsuba/render/shape.h index 2e02fa77..882e75ed 100644 --- a/include/mitsuba/render/shape.h +++ b/include/mitsuba/render/shape.h @@ -111,9 +111,6 @@ public: /// Return a string representation std::string toString() const; public: - /// Incident direction in the local frame - Vector wi; - /// Distance traveled along the ray Float t; @@ -129,7 +126,7 @@ public: /// UV surface coordinates Point2 uv; - /// Position partials wrt. to changes in texture-space + /// Position partials wrt. the texture space parameterization Vector dpdu, dpdv; /// Texture coordinate mapping partials wrt. changes in screen-space @@ -138,6 +135,9 @@ public: /// Interpolated vertex color Spectrum color; + /// Incident direction in the local frame + Vector wi; + /// Affected shape const Shape *shape; @@ -197,13 +197,6 @@ public: virtual void fillIntersectionRecord(const Ray &ray, Float t, const void *temp, Intersection &its) const; -#if defined(MTS_SSE) - /// Perform 4 simultaneous intersection tests using SSE -// virtual __m128 rayIntersectPacket(const RayPacket4 &packet, const - // __m128 mint, __m128 maxt, __m128 inactive, - // Intersection4 &its) const; -#endif - /** * \brief Sample a point on the shape * diff --git a/include/mitsuba/render/triaccel.h b/include/mitsuba/render/triaccel.h index d0209cb5..97991997 100644 --- a/include/mitsuba/render/triaccel.h +++ b/include/mitsuba/render/triaccel.h @@ -44,8 +44,8 @@ struct TriAccel { Float c_nu; Float c_nv; - uint32_t index; uint32_t shapeIndex; + uint32_t primIndex; /// Construct from vertex data. Returns '1' if there was a failure inline int load(const Point &A, const Point &B, const Point &C); diff --git a/src/librender/kdtree.cpp b/src/librender/kdtree.cpp index 37a7e38b..46fb6c14 100644 --- a/src/librender/kdtree.cpp +++ b/src/librender/kdtree.cpp @@ -86,7 +86,7 @@ void KDTree::build() { const Point &v2 = positions[tri.idx[2]]; m_triAccel[idx].load(v0, v1, v2); m_triAccel[idx].shapeIndex = i; - m_triAccel[idx].index = j; + m_triAccel[idx].primIndex = j; ++idx; } } else { @@ -128,7 +128,7 @@ bool KDTree::rayIntersect(const Ray &ray, Intersection &its) const { const Shape *shape = m_shapes[cache->shapeIndex]; if (m_triangleFlag[cache->shapeIndex]) { const TriMesh *trimesh = static_cast(shape); - const Triangle &tri = trimesh->getTriangles()[cache->index]; + const Triangle &tri = trimesh->getTriangles()[cache->primIndex]; const Point *vertexPositions = trimesh->getVertexPositions(); const Normal *vertexNormals = trimesh->getVertexNormals(); const Point2 *vertexTexcoords = trimesh->getVertexTexcoords(); @@ -195,7 +195,8 @@ bool KDTree::rayIntersect(const Ray &ray, Intersection &its) const { return true; } else { shape->fillIntersectionRecord(ray, its.t, - reinterpret_cast(temp) + 4, its); + reinterpret_cast(temp) + 8, its); + return true; } } } @@ -231,7 +232,7 @@ static StatsCounter coherentPackets("General", "Coherent ray packets"); static StatsCounter incoherentPackets("General", "Incoherent ray packets"); void KDTree::rayIntersectPacket(const RayPacket4 &packet, - const RayInterval4 &rayInterval, Intersection4 &its) const { + const RayInterval4 &rayInterval, Intersection4 &its, void *temp) const { CoherentKDStackEntry MM_ALIGN16 stack[MTS_KD_MAXDEPTH]; RayInterval4 MM_ALIGN16 interval; @@ -248,9 +249,9 @@ void KDTree::rayIntersectPacket(const RayPacket4 &packet, interval.mint.ps = _mm_max_ps(interval.mint.ps, rayInterval.mint.ps); interval.maxt.ps = _mm_min_ps(interval.maxt.ps, rayInterval.maxt.ps); - __m128 itsFound = _mm_cmpgt_ps(interval.mint.ps, interval.maxt.ps), - masked = itsFound; - if (_mm_movemask_ps(itsFound) == 0xF) + SSEVector itsFound( _mm_cmpgt_ps(interval.mint.ps, interval.maxt.ps)); + __m128 masked = itsFound.ps; + if (_mm_movemask_ps(itsFound.ps) == 0xF) return; while (currNode != NULL) { @@ -295,54 +296,57 @@ void KDTree::rayIntersectPacket(const RayPacket4 &packet, const index_type primEnd = currNode->getPrimEnd(); if (EXPECT_NOT_TAKEN(primStart != primEnd)) { -#ifdef MTS_USE_TRIACCEL4 - const int count = m_packedTriangles[primStart].indirectionCount; - primStart = m_packedTriangles[primStart].indirectionIndex; - primEnd = primStart + count; -#endif - __m128 - searchStart = _mm_max_ps(rayInterval.mint.ps, - _mm_mul_ps(interval.mint.ps, SSEConstants::om_eps.ps)), - searchEnd = _mm_min_ps(rayInterval.maxt.ps, - _mm_mul_ps(interval.maxt.ps, SSEConstants::op_eps.ps)); + SSEVector + searchStart(_mm_max_ps(rayInterval.mint.ps, + _mm_mul_ps(interval.mint.ps, SSEConstants::om_eps.ps))), + searchEnd(_mm_min_ps(rayInterval.maxt.ps, + _mm_mul_ps(interval.maxt.ps, SSEConstants::op_eps.ps))); for (index_type entry=primStart; entry != primEnd; entry++) { const TriAccel &kdTri = m_triAccel[m_indices[entry]]; if (EXPECT_TAKEN(kdTri.k != KNoTriangleFlag)) { - itsFound = _mm_or_ps(itsFound, - kdTri.rayIntersectPacket(packet, searchStart, searchEnd, masked, its)); + itsFound.ps = _mm_or_ps(itsFound.ps, + kdTri.rayIntersectPacket(packet, searchStart.ps, searchEnd.ps, masked, its)); } else { -#if 0 - /* Not a triangle - invoke the shape's intersection routine */ - __m128 hasIts = m_shapes[kdTri.shapeIndex]->rayIntersectPacket(packet, - searchStart, searchEnd, masked, its); - itsFound = _mm_or_ps(itsFound, hasIts); - its.primIndex.pi = mux_epi32(pstoepi32(hasIts), - load1_epi32(kdTri.index), its.primIndex.pi); - its.shapeIndex.pi = mux_epi32(pstoepi32(hasIts), - load1_epi32(kdTri.shapeIndex), its.shapeIndex.pi); -#endif + const Shape *shape = m_shapes[kdTri.shapeIndex]; + + for (int i=0; i<4; ++i) { + Ray ray; + for (int axis=0; axis<3; axis++) { + ray.o[axis] = packet.o[axis].f[i]; + ray.d[axis] = packet.d[axis].f[i]; + } + Float t; + + if (shape->rayIntersect(ray, searchStart.f[i], searchEnd.f[i], t, + reinterpret_cast(temp) + + i * MTS_KD_INTERSECTION_TEMP + 8)) { + its.t.f[i] = t; + its.shapeIndex.i[i] = kdTri.shapeIndex; + its.primIndex.i[i] = KNoTriangleFlag; + itsFound.i[i] = 0xFFFFFFFF; + } + } } - searchEnd = _mm_min_ps(searchEnd, its.t.ps); + searchEnd.ps = _mm_min_ps(searchEnd.ps, its.t.ps); } } /* Abort if the tree has been traversed or if intersections have been found for all four rays */ - if (_mm_movemask_ps(itsFound) == 0xF || --stackIndex < 0) + if (_mm_movemask_ps(itsFound.ps) == 0xF || --stackIndex < 0) break; /* Pop from the stack */ currNode = stack[stackIndex].node; interval = stack[stackIndex].interval; - masked = _mm_or_ps(itsFound, + masked = _mm_or_ps(itsFound.ps, _mm_cmpgt_ps(interval.mint.ps, interval.maxt.ps)); } } void KDTree::rayIntersectPacketIncoherent(const RayPacket4 &packet, - const RayInterval4 &rayInterval, Intersection4 &its4) const { - uint8_t temp[MTS_KD_INTERSECTION_TEMP]; + const RayInterval4 &rayInterval, Intersection4 &its4, void *temp) const { ++incoherentPackets; for (int i=0; i<4; i++) { @@ -354,12 +358,13 @@ void KDTree::rayIntersectPacketIncoherent(const RayPacket4 &packet, } ray.mint = rayInterval.mint.f[i]; ray.maxt = rayInterval.maxt.f[i]; - if (ray.mint < ray.maxt && rayIntersectHavran(ray, ray.mint, ray.maxt, its4.t.f[i], temp)) { + if (ray.mint < ray.maxt && rayIntersectHavran(ray, ray.mint, + ray.maxt, its4.t.f[i], reinterpret_cast(temp) + i * MTS_KD_INTERSECTION_TEMP)) { const IntersectionCache *cache = reinterpret_cast(temp); its4.u.f[i] = cache->u; - its4.v.f[i] = cache->u; + its4.v.f[i] = cache->v; its4.shapeIndex.i[i] = cache->shapeIndex; - its4.primIndex.i[i] = cache->index; + its4.primIndex.i[i] = cache->primIndex; } } } diff --git a/src/librender/preview.cpp b/src/librender/preview.cpp index 318c9d25..11efca76 100644 --- a/src/librender/preview.cpp +++ b/src/librender/preview.cpp @@ -137,6 +137,7 @@ void PreviewWorker::processCoherent(const WorkUnit *workUnit, WorkResult *workRe const SSEVector MM_ALIGN16 yOffset(0.0f, 0.0f, 1.0f, 1.0f); const int pixelOffset[] = {0, 1, width, width+1}; const __m128 clamping = _mm_set1_ps(1/(m_minDist*m_minDist)); + uint8_t temp[MTS_KD_INTERSECTION_TEMP*4]; const __m128 camTL[3] = { _mm_set1_ps(m_cameraTL.x), @@ -232,9 +233,9 @@ void PreviewWorker::processCoherent(const WorkUnit *workUnit, WorkResult *workRe primRay4.signs[0][0] = primSignsX ? 1 : 0; primRay4.signs[1][0] = primSignsY ? 1 : 0; primRay4.signs[2][0] = primSignsZ ? 1 : 0; - m_kdtree->rayIntersectPacket(primRay4, itv4, its4); + m_kdtree->rayIntersectPacket(primRay4, itv4, its4, temp); } else { - m_kdtree->rayIntersectPacketIncoherent(primRay4, itv4, its4); + m_kdtree->rayIntersectPacketIncoherent(primRay4, itv4, its4, temp); } numRays += 4; @@ -288,7 +289,6 @@ void PreviewWorker::processCoherent(const WorkUnit *workUnit, WorkResult *workRe const Shape *shape = (*m_shapes)[its4.shapeIndex.i[idx]]; const BSDF *bsdf = shape->getBSDF(); - if (EXPECT_TAKEN(primIndex != KNoTriangleFlag)) { const TriMesh *mesh = static_cast(shape); const Triangle &t = mesh->getTriangles()[primIndex]; @@ -343,13 +343,12 @@ void PreviewWorker::processCoherent(const WorkUnit *workUnit, WorkResult *workRe its.dpdv = t0.dpdv * alpha + t1.dpdv * beta + t2.dpdv * gamma; } } else { -#if 0 Ray ray( Point(primRay4.o[0].f[idx], primRay4.o[1].f[idx], primRay4.o[2].f[idx]), Vector(primRay4.d[0].f[idx], primRay4.d[1].f[idx], primRay4.d[2].f[idx]) ); - shape->rayIntersect(ray, its); -#endif + shape->fillIntersectionRecord(ray, its4.t.f[idx], temp + + + idx * MTS_KD_INTERSECTION_TEMP + 8, its); } wo.x = nSecD[0].f[idx]; wo.y = nSecD[1].f[idx]; wo.z = nSecD[2].f[idx]; @@ -427,9 +426,9 @@ void PreviewWorker::processCoherent(const WorkUnit *workUnit, WorkResult *workRe secRay4.signs[0][0] = secSignsX ? 1 : 0; secRay4.signs[1][0] = secSignsY ? 1 : 0; secRay4.signs[2][0] = secSignsZ ? 1 : 0; - m_kdtree->rayIntersectPacket(secRay4, secItv4, secIts4); + m_kdtree->rayIntersectPacket(secRay4, secItv4, secIts4, temp); } else { - m_kdtree->rayIntersectPacketIncoherent(secRay4, secItv4, secIts4); + m_kdtree->rayIntersectPacketIncoherent(secRay4, secItv4, secIts4, temp); } for (int idx=0; idx<4; ++idx) { diff --git a/src/librender/scene.cpp b/src/librender/scene.cpp index 183d0365..815327af 100644 --- a/src/librender/scene.cpp +++ b/src/librender/scene.cpp @@ -239,7 +239,7 @@ void Scene::configure() { Float maxExtents = std::max(extents.x, extents.y); Float distance = maxExtents/(2.0f * std::tan(45 * .5f * M_PI/180)); - props.setTransform("toWorld", Transform::translate(Vector(center.x, center.y, aabb.getMinimum().x - distance))); + props.setTransform("toWorld", Transform::translate(Vector(center.x, center.y, aabb.getMinimum().z - distance))); props.setFloat("fov", 45.0f); m_camera = static_cast (PluginManager::getInstance()->createObject(Camera::m_theClass, props)); diff --git a/src/shapes/sphere.cpp b/src/shapes/sphere.cpp index 00e91003..bd5b2ddf 100644 --- a/src/shapes/sphere.cpp +++ b/src/shapes/sphere.cpp @@ -29,136 +29,116 @@ MTS_NAMESPACE_BEGIN */ class Sphere : public Shape { -private: - Point m_center; - Float m_radius; public: Sphere(const Properties &props) : Shape(props) { - m_radius = props.getFloat("radius", 1.0f); // Negative radius -> inside-out sphere - if (m_objectToWorld.hasScale()) - Log(EError, "The scale needs to be specified using the 'radius' parameter!"); + m_objectToWorld = props.getTransform("toWorld", Transform()); + m_center = m_objectToWorld(Point(0,0,0)); + /// non-uniform scales are not supported! + m_radius = m_objectToWorld(Vector(1,0,0)).length(); + m_worldToObject = m_objectToWorld.inverse(); + m_invSurfaceArea = 1/(4*M_PI*m_radius*m_radius); } - + Sphere(Stream *stream, InstanceManager *manager) : Shape(stream, manager) { + m_objectToWorld = Transform(stream); m_radius = stream->readFloat(); - configure(); + m_center = Point(stream); + m_worldToObject = m_objectToWorld.inverse(); + m_invSurfaceArea = 1/(4*M_PI*m_radius*m_radius); } - void configure() { - Shape::configure(); + void serialize(Stream *stream, InstanceManager *manager) const { + Shape::serialize(stream, manager); + m_objectToWorld.serialize(stream); + stream->writeFloat(m_radius); + m_center.serialize(stream); + } - m_surfaceArea = 4*M_PI*m_radius*m_radius; - m_invSurfaceArea = 1.0f / m_surfaceArea; - m_center = m_objectToWorld(Point(0,0,0)); + AABB getAABB() const { + AABB aabb; Float absRadius = std::abs(m_radius); - m_aabb.min = m_center - Vector(absRadius, absRadius, absRadius); - m_aabb.max = m_center + Vector(absRadius, absRadius, absRadius); - m_bsphere.center = m_center; - m_bsphere.radius = absRadius; + aabb.min = m_center - Vector(absRadius); + aabb.max = m_center + Vector(absRadius); + return aabb; } - bool rayIntersect(const Ray &ray, Float start, Float end, Float &t) const { + Float getSurfaceArea() const { + return 4*M_PI*m_radius*m_radius; + } + + bool rayIntersect(const Ray &ray, Float mint, Float maxt, Float &t, void *tmp) const { + Vector ro = ray.o - m_center; + /* Transform into the local coordinate system and normalize */ - double nearT, farT; - const double ox = (double) ray.o.x - (double) m_center.x, - oy = (double) ray.o.y - (double) m_center.y, - oz = (double) ray.o.z - (double) m_center.z; - const double dx = ray.d.x, dy = ray.d.y, dz = ray.d.z; - const double A = dx*dx + dy*dy + dz*dz; - const double B = 2 * (dx*ox + dy*oy + dz*oz); - const double C = ox*ox + oy*oy + oz*oz - m_radius * m_radius; + Float A = ray.d.x*ray.d.x + ray.d.y*ray.d.y + ray.d.z*ray.d.z; + Float B = 2 * (ray.d.x*ro.x + ray.d.y*ro.y + ray.d.z*ro.z); + Float C = ro.x*ro.x + ro.y*ro.y + + ro.z*ro.z - m_radius*m_radius; - if (!solveQuadraticDouble(A, B, C, nearT, farT)) + Float nearT, farT; + if (!solveQuadratic(A, B, C, nearT, farT)) return false; - if (nearT > end || farT < start) + if (nearT > maxt || farT < mint) return false; - if (nearT < start) { - if (farT > end) + if (nearT < mint) { + if (farT > maxt) return false; - t = (Float) farT; + t = farT; } else { - t = (Float) nearT; + t = nearT; } + return true; } - bool rayIntersect(const Ray &ray, Intersection &its) const { - if (!rayIntersect(ray, ray.mint, ray.maxt, its.t)) - return false; - its.p = ray(its.t); - - Point local = m_worldToObject(its.p); - Float absRadius = std::abs(m_radius); - Float theta = std::acos(local.z / absRadius); + void fillIntersectionRecord(const Ray &ray, Float t, + const void *temp, Intersection &its) const { + its.t = t; + its.p = ray(t); + Vector local = normalize(m_worldToObject(its.p - m_center)); + Float theta = std::acos(local.z); Float phi = std::atan2(local.y, local.x); + if (phi < 0) phi += 2*M_PI; - its.uv.x = theta / M_PI; - its.uv.y = phi / (2*M_PI); - its.shape = this; + + its.uv.x = phi * (0.5 * INV_PI); + its.uv.y = theta * INV_PI; + its.dpdu = m_objectToWorld(Vector(-local.y, local.x, 0) * (2*M_PI)); Float zrad = std::sqrt(local.x*local.x + local.y*local.y); - its.geoFrame.n = Normal(normalize(its.p - m_center)); - if (m_radius < 0) - its.geoFrame.n *= -1; - if (zrad > 0) { - Float invZRad = 1.0f / zrad; - Float cosPhi = local.x * invZRad; - Float sinPhi = local.y * invZRad; - its.dpdu = m_objectToWorld(Vector(-local.y, local.x, 0) * (2*M_PI)); + Float invZRad = 1.0f / zrad, + cosPhi = local.x * invZRad, + sinPhi = local.y * invZRad; its.dpdv = m_objectToWorld(Vector(local.z * cosPhi, local.z * sinPhi, - - absRadius * std::sin(theta)) * M_PI); - its.geoFrame.s = normalize(its.dpdu - its.geoFrame.n - * dot(its.geoFrame.n, its.dpdu)); - its.geoFrame.t = cross(its.geoFrame.n, its.geoFrame.s); + -std::sin(theta)) * M_PI); + its.geoFrame.s = normalize(its.dpdu); + its.geoFrame.t = normalize(its.dpdv); + its.geoFrame.n = normalize(cross(its.dpdv, its.dpdu)); } else { // avoid a singularity - Float cosPhi = 0, sinPhi = 1; + const Float cosPhi = 0, sinPhi = 1; its.dpdv = m_objectToWorld(Vector(local.z * cosPhi, local.z * sinPhi, - -absRadius*std::sin(theta))* M_PI); - its.dpdu = cross(its.dpdv, its.geoFrame.n); - coordinateSystem(its.geoFrame.n, its.geoFrame.s, its.geoFrame.t); + -std::sin(theta)) * M_PI); + its.geoFrame = Frame(normalize(its.p - m_center)); + if (m_radius < 0) + its.geoFrame.n *= -1; } its.shFrame = its.geoFrame; its.wi = its.toLocal(-ray.d); + its.shape = this; its.hasUVPartials = false; - - return true; } -#if defined(MTS_SSE) - /* SSE-accelerated packet tracing is not supported for spheres at the moment */ - __m128 rayIntersectPacket(const RayPacket4 &packet, const - __m128 start, __m128 end, __m128 inactive, Intersection4 &its) const { - SSEVector result(_mm_setzero_ps()), mint(start), maxt(end), mask(inactive); - Float t; - - for (int i=0; i<4; i++) { - Ray ray; - for (int axis=0; axis<3; axis++) { - ray.o[axis] = packet.o[axis].f[i]; - ray.d[axis] = packet.d[axis].f[i]; - } - if (mask.i[i] != 0) - continue; - if (rayIntersect(ray, mint.f[i], maxt.f[i], t)) { - result.i[i] = 0xFFFFFFFF; - its.t.f[i] = t; - } - } - return result.ps; - } -#endif - Float sampleArea(ShapeSamplingRecord &sRec, const Point2 &sample) const { Vector v = squareToSphere(sample); - sRec.n = m_objectToWorld(Normal(v)); - sRec.p = m_objectToWorld(Point(v * m_radius)); - return m_invSurfaceArea; + sRec.n = Normal(v); + sRec.p = Point(v * m_radius) + m_center; + return 1.0f / (4*M_PI*m_radius*m_radius); } /** @@ -193,7 +173,7 @@ public: Ray ray(p, d); Float t; - if (!rayIntersect(ray, 0, std::numeric_limits::infinity(), t)) { + if (!rayIntersect(ray, 0, std::numeric_limits::infinity(), t, NULL)) { // This can happen sometimes due to roundoff errors - just fail to // generate a sample in this case. return 0; @@ -223,28 +203,25 @@ public: return squareToConePdf(cosThetaMax); } - void serialize(Stream *stream, InstanceManager *manager) const { - Shape::serialize(stream, manager); - stream->writeFloat(m_radius); - } - std::string toString() const { std::ostringstream oss; oss << "Sphere[" << endl << " radius = " << m_radius << ", " << endl << " center = " << m_center.toString() << ", " << endl - << " objectToWorld = " << indent(m_objectToWorld.toString()) << "," << endl - << " aabb = " << m_aabb.toString() << "," << endl - << " bsphere = " << m_bsphere.toString() << "," << endl << " bsdf = " << indent(m_bsdf.toString()) << "," << endl << " luminaire = " << indent(m_luminaire.toString()) << "," << endl - << " subsurface = " << indent(m_subsurface.toString()) << "," << endl - << " surfaceArea = " << m_surfaceArea << endl + << " subsurface = " << indent(m_subsurface.toString()) << endl << "]"; return oss.str(); } MTS_DECLARE_CLASS() +private: + Transform m_objectToWorld; + Transform m_worldToObject; + Point m_center; + Float m_radius; + Float m_invSurfaceArea; }; MTS_IMPLEMENT_CLASS_S(Sphere, false, Shape)