From 98b366d5ddd43eaa9ee476a12ae364bedda2f94e Mon Sep 17 00:00:00 2001 From: Benedikt Bitterli Date: Thu, 29 Nov 2018 15:13:46 -0500 Subject: [PATCH] Implement non-exponential transport --- src/core/grids/Grid.hpp | 5 +- src/core/grids/VdbGrid.cpp | 48 ++++--- src/core/grids/VdbGrid.hpp | 5 +- src/core/integrators/TraceBase.cpp | 22 ++-- src/core/integrators/TraceBase.hpp | 3 + .../bidirectional_path_tracer/LightPath.cpp | 18 ++- .../bidirectional_path_tracer/PathVertex.cpp | 10 +- .../bidirectional_path_tracer/PathVertex.hpp | 1 + .../integrators/light_tracer/LightTracer.cpp | 2 +- .../integrators/path_tracer/PathTracer.cpp | 1 + .../integrators/photon_map/PhotonTracer.cpp | 10 +- src/core/io/Scene.cpp | 5 + src/core/io/Scene.hpp | 3 + src/core/media/AtmosphericMedium.cpp | 58 ++++----- src/core/media/AtmosphericMedium.hpp | 9 +- src/core/media/ExponentialMedium.cpp | 73 ++++------- src/core/media/ExponentialMedium.hpp | 9 +- src/core/media/HomogeneousMedium.cpp | 98 +++------------ src/core/media/HomogeneousMedium.hpp | 8 +- src/core/media/Medium.cpp | 19 ++- src/core/media/Medium.hpp | 10 +- src/core/media/VoxelMedium.cpp | 64 ++++------ src/core/media/VoxelMedium.hpp | 7 +- .../transmittances/DavisTransmittance.cpp | 62 +++++++++ .../transmittances/DavisTransmittance.hpp | 31 +++++ .../DavisWeinsteinTransmittance.cpp | 119 ++++++++++++++++++ .../DavisWeinsteinTransmittance.hpp | 34 +++++ .../DoubleExponentialTransmittance.cpp | 68 ++++++++++ .../DoubleExponentialTransmittance.hpp | 35 ++++++ .../transmittances/ErlangTransmittance.cpp | 69 ++++++++++ .../transmittances/ErlangTransmittance.hpp | 34 +++++ .../ExponentialTransmittance.cpp | 57 +++++++++ .../ExponentialTransmittance.hpp | 32 +++++ .../InterpolatedTransmittance.cpp | 85 +++++++++++++ .../InterpolatedTransmittance.hpp | 39 ++++++ .../transmittances/LinearTransmittance.cpp | 79 ++++++++++++ .../transmittances/LinearTransmittance.hpp | 33 +++++ .../transmittances/PulseTransmittance.cpp | 108 ++++++++++++++++ .../transmittances/PulseTransmittance.hpp | 37 ++++++ .../transmittances/QuadraticTransmittance.cpp | 68 ++++++++++ .../transmittances/QuadraticTransmittance.hpp | 31 +++++ src/core/transmittances/Transmittance.hpp | 62 +++++++++ .../transmittances/TransmittanceFactory.cpp | 26 ++++ .../transmittances/TransmittanceFactory.hpp | 17 +++ 44 files changed, 1327 insertions(+), 287 deletions(-) create mode 100644 src/core/transmittances/DavisTransmittance.cpp create mode 100644 src/core/transmittances/DavisTransmittance.hpp create mode 100644 src/core/transmittances/DavisWeinsteinTransmittance.cpp create mode 100644 src/core/transmittances/DavisWeinsteinTransmittance.hpp create mode 100644 src/core/transmittances/DoubleExponentialTransmittance.cpp create mode 100644 src/core/transmittances/DoubleExponentialTransmittance.hpp create mode 100644 src/core/transmittances/ErlangTransmittance.cpp create mode 100644 src/core/transmittances/ErlangTransmittance.hpp create mode 100644 src/core/transmittances/ExponentialTransmittance.cpp create mode 100644 src/core/transmittances/ExponentialTransmittance.hpp create mode 100644 src/core/transmittances/InterpolatedTransmittance.cpp create mode 100644 src/core/transmittances/InterpolatedTransmittance.hpp create mode 100644 src/core/transmittances/LinearTransmittance.cpp create mode 100644 src/core/transmittances/LinearTransmittance.hpp create mode 100644 src/core/transmittances/PulseTransmittance.cpp create mode 100644 src/core/transmittances/PulseTransmittance.hpp create mode 100644 src/core/transmittances/QuadraticTransmittance.cpp create mode 100644 src/core/transmittances/QuadraticTransmittance.hpp create mode 100644 src/core/transmittances/Transmittance.hpp create mode 100644 src/core/transmittances/TransmittanceFactory.cpp create mode 100644 src/core/transmittances/TransmittanceFactory.hpp diff --git a/src/core/grids/Grid.hpp b/src/core/grids/Grid.hpp index 33ae5e17..66220bb0 100644 --- a/src/core/grids/Grid.hpp +++ b/src/core/grids/Grid.hpp @@ -20,9 +20,8 @@ class Grid : public JsonSerializable virtual Box3f bounds() const; virtual float density(Vec3f p) const = 0; - virtual Vec3f transmittance(PathSampleGenerator &sampler, Vec3f p, Vec3f w, float t0, float t1, Vec3f sigmaT) const = 0; - virtual Vec2f inverseOpticalDepth(PathSampleGenerator &sampler, Vec3f p, Vec3f w, float t0, float t1, - float sigmaT, float xi) const = 0; + virtual float opticalDepth(PathSampleGenerator &sampler, Vec3f p, Vec3f w, float t0, float t1) const = 0; + virtual Vec2f inverseOpticalDepth(PathSampleGenerator &sampler, Vec3f p, Vec3f w, float t0, float t1, float xi) const = 0; }; } diff --git a/src/core/grids/VdbGrid.cpp b/src/core/grids/VdbGrid.cpp index 7f2996a0..388c0a70 100644 --- a/src/core/grids/VdbGrid.cpp +++ b/src/core/grids/VdbGrid.cpp @@ -237,7 +237,7 @@ float VdbGrid::density(Vec3f p) const return gridAt(_grid->tree(), p); } -Vec3f VdbGrid::transmittance(PathSampleGenerator &sampler, Vec3f p, Vec3f w, float t0, float t1, Vec3f sigmaT) const +float VdbGrid::opticalDepth(PathSampleGenerator &sampler, Vec3f p, Vec3f w, float t0, float t1) const { auto accessor = _grid->getConstAccessor(); @@ -249,7 +249,7 @@ Vec3f VdbGrid::transmittance(PathSampleGenerator &sampler, Vec3f p, Vec3f w, flo integral += accessor.getValue(voxel)*(tb - ta); return false; }); - return std::exp(-integral*sigmaT); + return integral; } else if (_integrationMethod == IntegrationMethod::ExactLinear) { VdbRaymarcher dda; @@ -261,27 +261,24 @@ Vec3f VdbGrid::transmittance(PathSampleGenerator &sampler, Vec3f p, Vec3f w, flo fa = fb; return false; }); - return std::exp(-integral*sigmaT); + return integral; } else if (_integrationMethod == IntegrationMethod::ResidualRatio) { VdbRaymarcher dda; float scale = _supergridSubsample; float invScale = 1.0f/scale; - sigmaT *= scale; - - float sigmaTc = sigmaT.max(); auto superAccessor = _superGrid->getConstAccessor(); UniformSampler &generator = sampler.uniformGenerator(); float controlIntegral = 0.0f; - Vec3f Tr(1.0f); + float Tr = 1.0f; dda.march(DdaRay(p*invScale + 0.5f, w), t0*invScale, t1*invScale, superAccessor, [&](openvdb::Coord voxel, float ta, float tb) { openvdb::Vec2s v = superAccessor.getValue(voxel); float muC = v.x(); float muR = v.y(); - muR *= sigmaTc; + muR *= scale; controlIntegral += muC*(tb - ta); @@ -289,12 +286,12 @@ Vec3f VdbGrid::transmittance(PathSampleGenerator &sampler, Vec3f p, Vec3f w, flo ta -= BitManip::normalizedLog(generator.nextI())/muR; if (ta >= tb) break; - Tr *= 1.0f - sigmaT*((gridAt(accessor, p + w*ta*scale) - muC)/muR); + Tr *= 1.0f - scale*((gridAt(accessor, p + w*ta*scale) - muC)/muR); } return false; }); - return std::exp(-controlIntegral*sigmaT)*Tr; + return controlIntegral - std::log(Tr); } else { float ta = t0; float fa = gridAt(accessor, p + w*t0); @@ -308,12 +305,11 @@ Vec3f VdbGrid::transmittance(PathSampleGenerator &sampler, Vec3f p, Vec3f w, flo fa = fb; dT = _stepSize; } while (ta < t1); - return std::exp(-integral*sigmaT); + return integral; } } -Vec2f VdbGrid::inverseOpticalDepth(PathSampleGenerator &sampler, Vec3f p, Vec3f w, float t0, float t1, - float sigmaT, float xi) const +Vec2f VdbGrid::inverseOpticalDepth(PathSampleGenerator &sampler, Vec3f p, Vec3f w, float t0, float t1, float tau) const { auto accessor = _grid->getConstAccessor(); @@ -324,9 +320,9 @@ Vec2f VdbGrid::inverseOpticalDepth(PathSampleGenerator &sampler, Vec3f p, Vec3f Vec2f result(t1, 0.0f); bool exited = !dda.march(DdaRay(p + 0.5f, w), t0, t1, accessor, [&](openvdb::Coord voxel, float ta, float tb) { float v = accessor.getValue(voxel); - float delta = v*sigmaT*(tb - ta); - if (integral + delta >= xi) { - result = Vec2f(ta + (tb - ta)*(xi - integral)/delta, v); + float delta = v*(tb - ta); + if (integral + delta >= tau) { + result = Vec2f(ta + (tb - ta)*(tau - integral)/delta, v); return true; } integral += delta; @@ -341,11 +337,11 @@ Vec2f VdbGrid::inverseOpticalDepth(PathSampleGenerator &sampler, Vec3f p, Vec3f Vec2f result(t1, 0.0f); bool exited = !dda.march(DdaRay(p + 0.5f, w), t0, t1, accessor, [&](openvdb::Coord /*voxel*/, float ta, float tb) { float fb = gridAt(accessor, p + tb*w); - float delta = (fb + fa)*0.5f*sigmaT*(tb - ta); - if (integral + delta >= xi) { - float a = (fb - fa)*sigmaT; - float b = fa*sigmaT; - float c = (integral - xi)/(tb - ta); + float delta = (fb + fa)*0.5f*(tb - ta); + if (integral + delta >= tau) { + float a = (fb - fa); + float b = fa; + float c = (integral - tau)/(tb - ta); float mantissa = max(b*b - 2.0f*a*c, 0.0f); float x1 = (-b + std::sqrt(mantissa))/a; result = Vec2f(ta + (tb - ta)*x1, fa + (fb - fa)*x1); @@ -364,11 +360,11 @@ Vec2f VdbGrid::inverseOpticalDepth(PathSampleGenerator &sampler, Vec3f p, Vec3f do { float tb = min(ta + dT, t1); float fb = gridAt(accessor, p + w*tb); - float delta = (fa + fb)*sigmaT*0.5f*(tb - ta); - if (integral + delta >= xi) { - float a = (fb - fa)*sigmaT; - float b = fa*sigmaT; - float c = (integral - xi)/(tb - ta); + float delta = (fa + fb)*0.5f*(tb - ta); + if (integral + delta >= tau) { + float a = (fb - fa); + float b = fa; + float c = (integral - tau)/(tb - ta); float mantissa = max(b*b - 2.0f*a*c, 0.0f); float x1 = (-b + std::sqrt(mantissa))/a; return Vec2f(ta + (tb - ta)*x1, fa + (fb - fa)*x1); diff --git a/src/core/grids/VdbGrid.hpp b/src/core/grids/VdbGrid.hpp index f35dfee5..0cdbbb0c 100644 --- a/src/core/grids/VdbGrid.hpp +++ b/src/core/grids/VdbGrid.hpp @@ -68,9 +68,8 @@ class VdbGrid : public Grid virtual Box3f bounds() const override; float density(Vec3f p) const override; - Vec3f transmittance(PathSampleGenerator &sampler, Vec3f p, Vec3f w, float t0, float t1, Vec3f sigmaT) const override; - Vec2f inverseOpticalDepth(PathSampleGenerator &sampler, Vec3f p, Vec3f w, float t0, float t1, - float sigmaT, float xi) const override; + float opticalDepth(PathSampleGenerator &sampler, Vec3f p, Vec3f w, float t0, float t1) const override; + Vec2f inverseOpticalDepth(PathSampleGenerator &sampler, Vec3f p, Vec3f w, float t0, float t1, float tau) const override; }; } diff --git a/src/core/integrators/TraceBase.cpp b/src/core/integrators/TraceBase.cpp index dfba54e5..0cfeb82f 100644 --- a/src/core/integrators/TraceBase.cpp +++ b/src/core/integrators/TraceBase.cpp @@ -108,7 +108,7 @@ inline Vec3f TraceBase::generalizedShadowRayImpl(PathSampleGenerator &sampler, pdfForward *= forward; pdfBackward *= backward; } else { - throughput *= medium->transmittance(sampler, ray); + throughput *= medium->transmittance(sampler, ray, startsOnSurface, endsOnSurface); } } if (info.primitive == nullptr || info.primitive == endCap) @@ -125,10 +125,11 @@ inline Vec3f TraceBase::generalizedShadowRayImpl(PathSampleGenerator &sampler, } Vec3f TraceBase::generalizedShadowRay(PathSampleGenerator &sampler, Ray &ray, const Medium *medium, - const Primitive *endCap, int bounce) const + const Primitive *endCap, bool startsOnSurface, bool endsOnSurface, int bounce) const { float dummyA, dummyB; - return generalizedShadowRayImpl(sampler, ray, medium, endCap, bounce, false, false, dummyA, dummyB); + return generalizedShadowRayImpl(sampler, ray, medium, endCap, bounce, + startsOnSurface, endsOnSurface, dummyA, dummyB); } Vec3f TraceBase::generalizedShadowRayAndPdfs(PathSampleGenerator &sampler, Ray &ray, const Medium *medium, @@ -147,6 +148,7 @@ Vec3f TraceBase::attenuatedEmission(PathSampleGenerator &sampler, IntersectionTemporary &data, IntersectionInfo &info, int bounce, + bool startsOnSurface, Ray &ray, Vec3f *transmittance) { @@ -162,7 +164,7 @@ Vec3f TraceBase::attenuatedEmission(PathSampleGenerator &sampler, info.w = ray.dir(); light.intersectionInfo(data, info); - Vec3f shadow = generalizedShadowRay(sampler, ray, medium, &light, bounce); + Vec3f shadow = generalizedShadowRay(sampler, ray, medium, &light, startsOnSurface, true, bounce); if (transmittance) *transmittance = shadow; if (shadow == 0.0f) @@ -192,7 +194,7 @@ bool TraceBase::volumeLensSample(const Camera &camera, ray.setPrimaryRay(false); ray.setFarT(lensSample.dist); - Vec3f transmittance = generalizedShadowRay(sampler, ray, medium, nullptr, bounce); + Vec3f transmittance = generalizedShadowRay(sampler, ray, medium, nullptr, false, true, bounce); if (transmittance == 0.0f) return false; @@ -231,7 +233,7 @@ bool TraceBase::surfaceLensSample(const Camera &camera, ray.setPrimaryRay(false); ray.setFarT(sample.dist); - Vec3f transmittance = generalizedShadowRay(*event.sampler, ray, medium, nullptr, bounce); + Vec3f transmittance = generalizedShadowRay(*event.sampler, ray, medium, nullptr, true, true, bounce); if (transmittance == 0.0f) return false; @@ -270,7 +272,7 @@ Vec3f TraceBase::lightSample(const Primitive &light, IntersectionTemporary data; IntersectionInfo info; - Vec3f e = attenuatedEmission(*event.sampler, light, medium, sample.dist, data, info, bounce, ray, transmittance); + Vec3f e = attenuatedEmission(*event.sampler, light, medium, sample.dist, data, info, bounce, true, ray, transmittance); if (e == 0.0f) return Vec3f(0.0f); @@ -306,7 +308,7 @@ Vec3f TraceBase::bsdfSample(const Primitive &light, IntersectionTemporary data; IntersectionInfo info; - Vec3f e = attenuatedEmission(*event.sampler, light, medium, -1.0f, data, info, bounce, ray, nullptr); + Vec3f e = attenuatedEmission(*event.sampler, light, medium, -1.0f, data, info, bounce, true, ray, nullptr); if (e == Vec3f(0.0f)) return Vec3f(0.0f); @@ -338,7 +340,7 @@ Vec3f TraceBase::volumeLightSample(PathSampleGenerator &sampler, IntersectionTemporary data; IntersectionInfo info; - Vec3f e = attenuatedEmission(sampler, light, medium, lightSample.dist, data, info, bounce, ray, nullptr); + Vec3f e = attenuatedEmission(sampler, light, medium, lightSample.dist, data, info, bounce, false, ray, nullptr); if (e == 0.0f) return Vec3f(0.0f); @@ -366,7 +368,7 @@ Vec3f TraceBase::volumePhaseSample(const Primitive &light, IntersectionTemporary data; IntersectionInfo info; - Vec3f e = attenuatedEmission(sampler, light, medium, -1.0f, data, info, bounce, ray, nullptr); + Vec3f e = attenuatedEmission(sampler, light, medium, -1.0f, data, info, bounce, false, ray, nullptr); if (e == Vec3f(0.0f)) return Vec3f(0.0f); diff --git a/src/core/integrators/TraceBase.hpp b/src/core/integrators/TraceBase.hpp index 15ad8901..84e10c3e 100644 --- a/src/core/integrators/TraceBase.hpp +++ b/src/core/integrators/TraceBase.hpp @@ -64,6 +64,7 @@ class TraceBase IntersectionTemporary &data, IntersectionInfo &info, int bounce, + bool startsOnSurface, Ray &ray, Vec3f *transmittance); @@ -148,6 +149,8 @@ class TraceBase Ray &ray, const Medium *medium, const Primitive *endCap, + bool startsOnSurface, + bool endsOnSurface, int bounce) const; Vec3f generalizedShadowRayAndPdfs(PathSampleGenerator &sampler, Ray &ray, diff --git a/src/core/integrators/bidirectional_path_tracer/LightPath.cpp b/src/core/integrators/bidirectional_path_tracer/LightPath.cpp index 4a5a96a9..933d36d8 100644 --- a/src/core/integrators/bidirectional_path_tracer/LightPath.cpp +++ b/src/core/integrators/bidirectional_path_tracer/LightPath.cpp @@ -96,19 +96,25 @@ void LightPath::toAreaMeasure() float LightPath::misWeight(const LightPath &camera, const LightPath &emitter, const PathEdge &edge, int s, int t, float *ratios) { - float *pdfForward = reinterpret_cast(alloca((s + t)*sizeof(float))); - float *pdfBackward = reinterpret_cast(alloca((s + t)*sizeof(float))); - bool *connectable = reinterpret_cast(alloca((s + t)*sizeof(bool))); + float *pdfForward = reinterpret_cast(alloca((s + t)*sizeof(float))); + float *pdfBackward = reinterpret_cast(alloca((s + t)*sizeof(float))); + bool *connectable = reinterpret_cast(alloca((s + t)*sizeof(bool))); + const PathVertex **vertices = reinterpret_cast(alloca((s + t)*sizeof(const PathVertex *))); + + if (!camera[t - 1].segmentConnectable(emitter[s - 1])) + return 0.0f; for (int i = 0; i < s; ++i) { pdfForward [i] = emitter[i].pdfForward(); pdfBackward[i] = emitter[i].pdfBackward(); connectable[i] = !emitter[i].isDirac(); + vertices [i] = &emitter[i]; } for (int i = 0; i < t; ++i) { pdfForward [s + t - (i + 1)] = camera[i].pdfBackward(); pdfBackward[s + t - (i + 1)] = camera[i].pdfForward(); connectable[s + t - (i + 1)] = !camera[i].isDirac(); + vertices [s + t - (i + 1)] = &camera[i]; } connectable[s - 1] = connectable[s] = true; @@ -137,7 +143,7 @@ float LightPath::misWeight(const LightPath &camera, const LightPath &emitter, ratios[s] = 1.0f; for (int i = s + 1; i < s + t; ++i) { pi *= pdfForward[i - 1]/pdfBackward[i - 1]; - if (connectable[i - 1] && connectable[i]) { + if (connectable[i - 1] && connectable[i] && vertices[i - 1]->segmentConnectable(*vertices[i])) { weight += pi; if (ratios) ratios[i] = pi; @@ -149,7 +155,7 @@ float LightPath::misWeight(const LightPath &camera, const LightPath &emitter, pi = 1.0f; for (int i = s - 1; i >= 1; --i) { pi *= pdfBackward[i]/pdfForward[i]; - if (connectable[i - 1] && connectable[i]) { + if (connectable[i - 1] && connectable[i] && vertices[i - 1]->segmentConnectable(*vertices[i])) { weight += pi; if (ratios) ratios[i] = pi; @@ -166,7 +172,7 @@ float LightPath::misWeight(const LightPath &camera, const LightPath &emitter, } else { if (ratios) ratios[0] = 0.0f; - } + } return 1.0f/weight; } diff --git a/src/core/integrators/bidirectional_path_tracer/PathVertex.cpp b/src/core/integrators/bidirectional_path_tracer/PathVertex.cpp index f1d56298..3bbe7889 100644 --- a/src/core/integrators/bidirectional_path_tracer/PathVertex.cpp +++ b/src/core/integrators/bidirectional_path_tracer/PathVertex.cpp @@ -160,7 +160,7 @@ bool PathVertex::sampleNextVertex(const TraceableScene &scene, TraceBase &tracer hitSurface = mediumRecord.mediumSample.exited; edgePdfForward = mediumRecord.mediumSample.pdf; Ray reverseRay(mediumRecord.mediumSample.p, -state.ray.dir(), 0.0f, mediumRecord.mediumSample.t); - edgePdfBackward = state.medium->pdf(state.sampler, reverseRay, onSurface()); + edgePdfBackward = state.medium->pdf(state.sampler, reverseRay, hitSurface, onSurface()); weight *= mediumRecord.mediumSample.weight; if (hitSurface && !didHit) return false; @@ -327,6 +327,14 @@ void PathVertex::evalPdfs(const PathVertex *prev, const PathEdge *prevEdge, cons }} } +bool PathVertex::segmentConnectable(const PathVertex &next) const +{ + if (onSurface() || next.onSurface()) + return true; + + return !_medium->isDirac(); +} + void PathVertex::pointerFixup() { // Yuck. It's best not to ask. diff --git a/src/core/integrators/bidirectional_path_tracer/PathVertex.hpp b/src/core/integrators/bidirectional_path_tracer/PathVertex.hpp index b8d78d43..062bdd91 100644 --- a/src/core/integrators/bidirectional_path_tracer/PathVertex.hpp +++ b/src/core/integrators/bidirectional_path_tracer/PathVertex.hpp @@ -132,6 +132,7 @@ class PathVertex Vec3f eval(const Vec3f &d, bool adjoint) const; void evalPdfs(const PathVertex *prev, const PathEdge *prevEdge, const PathVertex &next, const PathEdge &nextEdge, float *forward, float *backward) const; + bool segmentConnectable(const PathVertex &next) const; void pointerFixup(); diff --git a/src/core/integrators/light_tracer/LightTracer.cpp b/src/core/integrators/light_tracer/LightTracer.cpp index c609696a..88910721 100644 --- a/src/core/integrators/light_tracer/LightTracer.cpp +++ b/src/core/integrators/light_tracer/LightTracer.cpp @@ -29,7 +29,7 @@ void LightTracer::traceSample(PathSampleGenerator &sampler) Ray shadowRay(point.p, splat.d); shadowRay.setFarT(splat.dist); - Vec3f transmission = generalizedShadowRay(sampler, shadowRay, medium, nullptr, 0); + Vec3f transmission = generalizedShadowRay(sampler, shadowRay, medium, nullptr, true, true, 0); if (transmission != 0.0f) { Vec3f value = throughput*transmission*splat.weight *light->evalDirectionalEmission(point, DirectionSample(splat.d)); diff --git a/src/core/integrators/path_tracer/PathTracer.cpp b/src/core/integrators/path_tracer/PathTracer.cpp index 17d84472..10ba1d66 100644 --- a/src/core/integrators/path_tracer/PathTracer.cpp +++ b/src/core/integrators/path_tracer/PathTracer.cpp @@ -50,6 +50,7 @@ Vec3f PathTracer::traceSample(Vec2u pixel, PathSampleGenerator &sampler) while ((didHit || medium) && bounce < _settings.maxBounces) { bool hitSurface = true; if (medium) { + mediumSample.continuedWeight = throughput; if (!medium->sampleDistance(sampler, ray, state, mediumSample)) return emission; throughput *= mediumSample.weight; diff --git a/src/core/integrators/photon_map/PhotonTracer.cpp b/src/core/integrators/photon_map/PhotonTracer.cpp index c3749938..4d1d4117 100644 --- a/src/core/integrators/photon_map/PhotonTracer.cpp +++ b/src/core/integrators/photon_map/PhotonTracer.cpp @@ -128,7 +128,7 @@ static bool evalBeam1D(const PhotonBeam &beam, PathSampleGenerator &sampler, con mediumQuery.setFarT(t); beamEstimate += medium->sigmaT(hitPoint)*invSinTheta/(2.0f*radius) *medium->phaseFunction(hitPoint)->eval(beam.dir, -ray.dir()) - *medium->transmittance(sampler, mediumQuery)*beam.power; + *medium->transmittance(sampler, mediumQuery, true, false)*beam.power; return true; } @@ -149,7 +149,7 @@ static bool evalPlane0D(const PhotonPlane0D &p, PathSampleGenerator &sampler, co mediumQuery.setFarT(t); beamEstimate += sqr(medium->sigmaT(hitPoint))*std::abs(invDet) *medium->phaseFunction(hitPoint)->eval(p.d1, -ray.dir()) - *medium->transmittance(sampler, mediumQuery)*p.power; + *medium->transmittance(sampler, mediumQuery, true, false)*p.power; return true; } @@ -187,7 +187,7 @@ bool evalPlane1D(const PhotonPlane1D &p, PathSampleGenerator &sampler, const Ray Ray mediumQuery(ray); mediumQuery.setFarT(t); - controlVariate -= medium->transmittance(sampler, mediumQuery)*(maxT - minT); + controlVariate -= medium->transmittance(sampler, mediumQuery, true, false)*(maxT - minT); } beamEstimate += sqr(medium->sigmaT(v2))*medium->phaseFunction(v2)->eval(p.d1, -ray.dir())*p.power*controlVariate; @@ -289,7 +289,7 @@ Vec3f PhotonTracer::traceSensorPath(Vec2u pixel, const KdTree &surfaceTr mediumQuery.setFarT(t); estimate += (3.0f*INV_PI*sqr(1.0f - distSq/p.radiusSq))/p.radiusSq *medium->phaseFunction(p.pos)->eval(p.dir, -ray.dir()) - *medium->transmittance(sampler, mediumQuery)*p.power; + *medium->transmittance(sampler, mediumQuery, true, false)*p.power; }; auto beamContribution = [&](uint32 photonIndex, const Vec3pf *bounds, float tMin, float tMax) { const PhotonBeam &beam = beams[photonIndex]; @@ -343,7 +343,7 @@ Vec3f PhotonTracer::traceSensorPath(Vec2u pixel, const KdTree &surfaceTr result += throughput*estimate; } - throughput *= medium->transmittance(sampler, ray); + throughput *= medium->transmittance(sampler, ray, true, true); } if (!didHit || !includeSurfaces) break; diff --git a/src/core/io/Scene.cpp b/src/core/io/Scene.cpp index 241263b5..65a5cd39 100644 --- a/src/core/io/Scene.cpp +++ b/src/core/io/Scene.cpp @@ -91,6 +91,11 @@ std::shared_ptr fetchObject(const std::vector> &list, cons } } +std::shared_ptr Scene::fetchTransmittance(JsonPtr value) const +{ + return instantiate(value, *this); +} + std::shared_ptr Scene::fetchPhase(JsonPtr value) const { return instantiate(value, *this); diff --git a/src/core/io/Scene.hpp b/src/core/io/Scene.hpp index ffdccbfa..2a8e8eb7 100644 --- a/src/core/io/Scene.hpp +++ b/src/core/io/Scene.hpp @@ -13,6 +13,8 @@ #include "ImageIO.hpp" #include "Path.hpp" +#include "transmittances/Transmittance.hpp" + #include "phasefunctions/PhaseFunction.hpp" #include "integrators/Integrator.hpp" @@ -66,6 +68,7 @@ class Scene : public JsonSerializable virtual void loadResources() override; virtual void saveResources() override; + std::shared_ptr fetchTransmittance(JsonPtr value) const; std::shared_ptr fetchPhase(JsonPtr value) const; std::shared_ptr fetchMedium(JsonPtr value) const; std::shared_ptr fetchGrid(JsonPtr value) const; diff --git a/src/core/media/AtmosphericMedium.cpp b/src/core/media/AtmosphericMedium.cpp index 29e37aa5..16bd588b 100644 --- a/src/core/media/AtmosphericMedium.cpp +++ b/src/core/media/AtmosphericMedium.cpp @@ -110,10 +110,10 @@ inline float AtmosphericMedium::densityIntegral(float h, float t0, float t1) con return (SQRT_PI*0.5f/s)*std::exp((-h*h + _radius*_radius)*s*s)*Erf::erfDifference(s*t0, s*t1); } -inline float AtmosphericMedium::inverseOpticalDepth(double h, double t0, double sigmaT, double xi) const +inline float AtmosphericMedium::inverseOpticalDepth(double h, double t0, double tau) const { double s = _effectiveFalloffScale; - double inner = std::erf(s*t0) - 2.0*double(INV_SQRT_PI)*std::exp(s*s*(h - _radius)*(h + _radius))*s*std::log(xi)/sigmaT; + double inner = std::erf(s*t0) + 2.0*double(INV_SQRT_PI)*std::exp(s*s*(h - _radius)*(h + _radius))*s*tau; if (inner >= 1.0) return Ray::infinity(); @@ -134,24 +134,26 @@ bool AtmosphericMedium::sampleDistance(PathSampleGenerator &sampler, const Ray & float maxT = ray.farT() + t0; if (_absorptionOnly) { sample.t = ray.farT(); - sample.weight = std::exp(-_sigmaT*densityIntegral(h, t0, maxT)); + Vec3f tau = densityIntegral(h, t0, maxT)*_sigmaT; + sample.weight = _transmittance->eval(tau, state.firstScatter, true); sample.pdf = 1.0f; sample.exited = true; } else { int component = sampler.nextDiscrete(3); float sigmaTc = _sigmaT[component]; - float xi = 1.0f - sampler.next1D(); + float tauC = _transmittance->sample(sampler, state.firstScatter)/sigmaTc; - float t = inverseOpticalDepth(h, t0, sigmaTc, xi); + float t = inverseOpticalDepth(h, t0, tauC); sample.t = min(t, maxT); - sample.weight = std::exp(-_sigmaT*densityIntegral(h, t0, sample.t)); + Vec3f tau = densityIntegral(h, t0, sample.t)*_sigmaT; + sample.weight = _transmittance->eval(tau, state.firstScatter, sample.exited); sample.exited = (t >= maxT); if (sample.exited) { - sample.pdf = sample.weight.avg(); + sample.pdf = _transmittance->surfaceProbability(tau, state.firstScatter).avg(); } else { float rho = density(h, sample.t); - sample.pdf = (rho*_sigmaT*sample.weight).avg(); - sample.weight *= rho*_sigmaS; + sample.pdf = (rho*_sigmaT*_transmittance->mediumPdf(tau, state.firstScatter)).avg(); + sample.weight *= rho*_sigmaS*_transmittance->sigmaBar(); } sample.weight /= sample.pdf; sample.t -= t0; @@ -164,17 +166,19 @@ bool AtmosphericMedium::sampleDistance(PathSampleGenerator &sampler, const Ray & return true; } -Vec3f AtmosphericMedium::transmittance(PathSampleGenerator &/*sampler*/, const Ray &ray) const +Vec3f AtmosphericMedium::transmittance(PathSampleGenerator &/*sampler*/, const Ray &ray, bool startOnSurface, + bool endOnSurface) const { Vec3f p = (ray.pos() - _center); float t0 = p.dot(ray.dir()); float t1 = ray.farT() + t0; float h = (p - t0*ray.dir()).length(); - return std::exp(-_sigmaT*densityIntegral(h, t0, t1)); + Vec3f tau = _sigmaT*densityIntegral(h, t0, t1); + return _transmittance->eval(tau, startOnSurface, endOnSurface); } -float AtmosphericMedium::pdf(PathSampleGenerator &/*sampler*/, const Ray &ray, bool onSurface) const +float AtmosphericMedium::pdf(PathSampleGenerator &/*sampler*/, const Ray &ray, bool startOnSurface, bool endOnSurface) const { if (_absorptionOnly) { return 1.0f; @@ -184,33 +188,13 @@ float AtmosphericMedium::pdf(PathSampleGenerator &/*sampler*/, const Ray &ray, b float t1 = ray.farT() + t0; float h = (p - t0*ray.dir()).length(); - Vec3f transmittance = std::exp(-_sigmaT*densityIntegral(h, t0, t1)); - if (onSurface) { - return transmittance.avg(); - } else { - return (density(h, t0)*_sigmaT*transmittance).avg(); - } - } -} -Vec3f AtmosphericMedium::transmittanceAndPdfs(PathSampleGenerator &/*sampler*/, const Ray &ray, bool startOnSurface, - bool endOnSurface, float &pdfForward, float &pdfBackward) const -{ - Vec3f p = (ray.pos() - _center); - float t0 = p.dot(ray.dir()); - float t1 = ray.farT() + t0; - float h = (p - t0*ray.dir()).length(); - - Vec3f transmittance = std::exp(-_sigmaT*densityIntegral(h, t0, t1)); - - if (_absorptionOnly) { - pdfForward = pdfBackward = 1.0f; - } else { - pdfForward = endOnSurface ? transmittance.avg() : (density(h, t1)*_sigmaT*transmittance).avg(); - pdfBackward = startOnSurface ? transmittance.avg() : (density(h, t0)*_sigmaT*transmittance).avg(); + Vec3f tau = _sigmaT*densityIntegral(h, t0, t1); + if (endOnSurface) + return _transmittance->surfaceProbability(tau, startOnSurface).avg(); + else + return (density(h, t0)*_sigmaT*_transmittance->mediumPdf(tau, startOnSurface)).avg(); } - - return transmittance; } } diff --git a/src/core/media/AtmosphericMedium.hpp b/src/core/media/AtmosphericMedium.hpp index 2dae826a..ced492d0 100644 --- a/src/core/media/AtmosphericMedium.hpp +++ b/src/core/media/AtmosphericMedium.hpp @@ -28,7 +28,7 @@ class AtmosphericMedium : public Medium inline float density(Vec3f p) const; inline float density(float h, float t0) const; inline float densityIntegral(float h, float t0, float t1) const; - inline float inverseOpticalDepth(double h, double t0, double sigmaT, double xi) const; + inline float inverseOpticalDepth(double h, double t0, double tau) const; public: AtmosphericMedium(); @@ -46,10 +46,9 @@ class AtmosphericMedium : public Medium virtual bool sampleDistance(PathSampleGenerator &sampler, const Ray &ray, MediumState &state, MediumSample &sample) const override; - virtual Vec3f transmittance(PathSampleGenerator &sampler, const Ray &ray) const override; - virtual float pdf(PathSampleGenerator &sampler, const Ray &ray, bool onSurface) const override; - virtual Vec3f transmittanceAndPdfs(PathSampleGenerator &sampler, const Ray &ray, bool startOnSurface, - bool endOnSurface, float &pdfForward, float &pdfBackward) const override; + virtual Vec3f transmittance(PathSampleGenerator &sampler, const Ray &ray, bool startOnSurface, + bool endOnSurface) const override; + virtual float pdf(PathSampleGenerator &sampler, const Ray &ray, bool startOnSurface, bool endOnSurface) const override; }; } diff --git a/src/core/media/ExponentialMedium.cpp b/src/core/media/ExponentialMedium.cpp index 25d674ad..fd462aaf 100644 --- a/src/core/media/ExponentialMedium.cpp +++ b/src/core/media/ExponentialMedium.cpp @@ -92,14 +92,13 @@ inline float ExponentialMedium::densityIntegral(float x, float dx, float tMax) c return (std::exp(-x) - std::exp(-dx*tMax - x))/dx; } -inline float ExponentialMedium::inverseOpticalDepth(float x, float dx, float sigmaT, float logXi) const +inline float ExponentialMedium::inverseOpticalDepth(float x, float dx, float tau) const { if (dx == 0.0f) { - float effectiveSigmaTc = sigmaT*std::exp(-x); - return -logXi/effectiveSigmaTc; + return tau/std::exp(-x); } else { - float denom = sigmaT + dx*std::exp(x)*logXi; - return denom <= 0.0f ? Ray::infinity() : std::log(sigmaT/denom)/dx; + float denom = 1.0f - dx*std::exp(x)*tau; + return denom <= 0.0f ? Ray::infinity() : -std::log(denom)/dx; } } @@ -117,25 +116,26 @@ bool ExponentialMedium::sampleDistance(PathSampleGenerator &sampler, const Ray & if (maxT == Ray::infinity() && dx <= 0.0f) return false; sample.t = maxT; - sample.weight = std::exp(-_sigmaT*densityIntegral(x, dx, ray.farT())); + Vec3f tau = densityIntegral(x, dx, ray.farT())*_sigmaT; + sample.weight = _transmittance->eval(tau, state.firstScatter, true); sample.pdf = 1.0f; sample.exited = true; } else { int component = sampler.nextDiscrete(3); float sigmaTc = _sigmaT[component]; - float xi = 1.0f - sampler.next1D(); - float logXi = std::log(xi); + float tauC = _transmittance->sample(sampler, state.firstScatter)/sigmaTc; - float t = inverseOpticalDepth(x, dx, sigmaTc, logXi); + float t = inverseOpticalDepth(x, dx, tauC); sample.t = min(t, maxT); - sample.weight = std::exp(-_sigmaT*densityIntegral(x, dx, sample.t)); + Vec3f tau = densityIntegral(x, dx, sample.t)*_sigmaT; + sample.weight = _transmittance->eval(tau, state.firstScatter, sample.exited); sample.exited = (t >= maxT); if (sample.exited) { - sample.pdf = sample.weight.avg(); + sample.pdf = _transmittance->surfaceProbability(tau, state.firstScatter).avg(); } else { float rho = density(x, dx, sample.t); - sample.pdf = (rho*_sigmaT*sample.weight).avg(); - sample.weight *= rho*_sigmaS; + sample.pdf = (rho*_sigmaT*_transmittance->mediumPdf(tau, state.firstScatter)).avg(); + sample.weight *= rho*_sigmaS*_transmittance->sigmaBar(); } sample.weight /= sample.pdf; @@ -146,18 +146,21 @@ bool ExponentialMedium::sampleDistance(PathSampleGenerator &sampler, const Ray & return true; } -Vec3f ExponentialMedium::transmittance(PathSampleGenerator &/*sampler*/, const Ray &ray) const +Vec3f ExponentialMedium::transmittance(PathSampleGenerator &/*sampler*/, const Ray &ray, bool startOnSurface, + bool endOnSurface) const { float x = _falloffScale*(ray.pos() - _unitPoint).dot(_unitFalloffDirection); float dx = _falloffScale*ray.dir().dot(_unitFalloffDirection); - if (ray.farT() == Ray::infinity() && dx <= 0.0f) + if (ray.farT() == Ray::infinity() && dx <= 0.0f) { return Vec3f(0.0f); - else - return std::exp(-_sigmaT*densityIntegral(x, dx, ray.farT())); + } else { + Vec3f tau = _sigmaT*densityIntegral(x, dx, ray.farT()); + return _transmittance->eval(tau, startOnSurface, endOnSurface); + } } -float ExponentialMedium::pdf(PathSampleGenerator &/*sampler*/, const Ray &ray, bool onSurface) const +float ExponentialMedium::pdf(PathSampleGenerator &/*sampler*/, const Ray &ray, bool startOnSurface, bool endOnSurface) const { if (_absorptionOnly) { return 1.0f; @@ -165,36 +168,12 @@ float ExponentialMedium::pdf(PathSampleGenerator &/*sampler*/, const Ray &ray, b float x = _falloffScale*(ray.pos() - _unitPoint).dot(_unitFalloffDirection); float dx = _falloffScale*ray.dir().dot(_unitFalloffDirection); - Vec3f transmittance = std::exp(-_sigmaT*densityIntegral(x, dx, ray.farT())); - if (onSurface) { - return transmittance.avg(); - } else { - return (density(x, dx, ray.farT())*_sigmaT*transmittance).avg(); - } + Vec3f tau = _sigmaT*densityIntegral(x, dx, ray.farT()); + if (endOnSurface) + return _transmittance->surfaceProbability(tau, startOnSurface).avg(); + else + return (density(x, dx, ray.farT())*_sigmaT*_transmittance->mediumPdf(tau, startOnSurface)).avg(); } } -Vec3f ExponentialMedium::transmittanceAndPdfs(PathSampleGenerator &/*sampler*/, const Ray &ray, bool startOnSurface, - bool endOnSurface, float &pdfForward, float &pdfBackward) const -{ - float x = _falloffScale*(ray.pos() - _unitPoint).dot(_unitFalloffDirection); - float dx = _falloffScale*ray.dir().dot(_unitFalloffDirection); - - if (ray.farT() == Ray::infinity() && dx <= 0.0f) { - pdfForward = pdfBackward = 0.0f; - return Vec3f(0.0f); - } - - Vec3f transmittance = std::exp(-_sigmaT*densityIntegral(x, dx, ray.farT())); - - if (_absorptionOnly) { - pdfForward = pdfBackward = 1.0f; - } else { - pdfForward = endOnSurface ? transmittance.avg() : (density(x, dx, ray.farT())*_sigmaT*transmittance).avg(); - pdfBackward = startOnSurface ? transmittance.avg() : (density(x, dx, 0.0f)*_sigmaT*transmittance).avg(); - } - - return transmittance; -} - } diff --git a/src/core/media/ExponentialMedium.hpp b/src/core/media/ExponentialMedium.hpp index fb0509ce..c60814b9 100644 --- a/src/core/media/ExponentialMedium.hpp +++ b/src/core/media/ExponentialMedium.hpp @@ -21,7 +21,7 @@ class ExponentialMedium : public Medium inline float density(Vec3f p) const; inline float density(float x, float dx, float t) const; inline float densityIntegral(float x, float dx, float tMax) const; - inline float inverseOpticalDepth(float x, float dx, float sigmaT, float logXi) const; + inline float inverseOpticalDepth(float x, float dx, float tau) const; public: ExponentialMedium(); @@ -39,10 +39,9 @@ class ExponentialMedium : public Medium virtual bool sampleDistance(PathSampleGenerator &sampler, const Ray &ray, MediumState &state, MediumSample &sample) const override; - virtual Vec3f transmittance(PathSampleGenerator &sampler, const Ray &ray) const override; - virtual float pdf(PathSampleGenerator &sampler, const Ray &ray, bool onSurface) const override; - virtual Vec3f transmittanceAndPdfs(PathSampleGenerator &sampler, const Ray &ray, bool startOnSurface, - bool endOnSurface, float &pdfForward, float &pdfBackward) const override; + virtual Vec3f transmittance(PathSampleGenerator &sampler, const Ray &ray, bool startOnSurface, + bool endOnSurface) const override; + virtual float pdf(PathSampleGenerator &sampler, const Ray &ray, bool startOnSurface, bool endOnSurface) const override; }; } diff --git a/src/core/media/HomogeneousMedium.cpp b/src/core/media/HomogeneousMedium.cpp index b06bb429..6e31b021 100644 --- a/src/core/media/HomogeneousMedium.cpp +++ b/src/core/media/HomogeneousMedium.cpp @@ -3,7 +3,6 @@ #include "sampling/PathSampleGenerator.hpp" #include "math/TangentFrame.hpp" -#include "math/FastMath.hpp" #include "math/Ray.hpp" #include "io/JsonObject.hpp" @@ -74,27 +73,29 @@ bool HomogeneousMedium::sampleDistance(PathSampleGenerator &sampler, const Ray & if (maxT == Ray::infinity()) return false; sample.t = maxT; - sample.weight = FastMath::exp(-_sigmaT*maxT); + sample.weight = _transmittance->eval(sample.t*_sigmaT, state.firstScatter, true); sample.pdf = 1.0f; sample.exited = true; } else { int component = sampler.nextDiscrete(3); float sigmaTc = _sigmaT[component]; - float t = -std::log(1.0f - sampler.next1D())/sigmaTc; + float t = _transmittance->sample(sampler, state.firstScatter)/sigmaTc; sample.t = min(t, maxT); sample.continuedT = t; - sample.weight = FastMath::exp(-sample.t*_sigmaT); - sample.continuedWeight = FastMath::exp(-sample.continuedT*_sigmaT); sample.exited = (t >= maxT); + Vec3f tau = sample.t*_sigmaT; + Vec3f continuedTau = sample.continuedT*_sigmaT; + sample.weight = _transmittance->eval(tau, state.firstScatter, sample.exited); + sample.continuedWeight = _transmittance->eval(continuedTau, state.firstScatter, sample.exited); if (sample.exited) { - sample.pdf = sample.weight.avg(); + sample.pdf = _transmittance->surfaceProbability(tau, state.firstScatter).avg(); } else { - sample.pdf = (_sigmaT*sample.weight).avg(); - sample.weight *= _sigmaS; + sample.pdf = (_sigmaT*_transmittance->mediumPdf(tau, state.firstScatter)).avg(); + sample.weight *= _sigmaS*_transmittance->sigmaBar(); } sample.weight /= sample.pdf; - sample.continuedWeight = _sigmaS*sample.continuedWeight/(_sigmaT*sample.continuedWeight).avg(); + sample.continuedWeight = _sigmaS*_transmittance->sigmaBar()*sample.continuedWeight/(_sigmaT*_transmittance->mediumPdf(continuedTau, state.firstScatter)).avg(); state.advance(); } @@ -104,87 +105,26 @@ bool HomogeneousMedium::sampleDistance(PathSampleGenerator &sampler, const Ray & return true; } -bool HomogeneousMedium::invertDistance(WritablePathSampleGenerator &sampler, const Ray &ray, bool onSurface) const -{ - float maxT = ray.farT(); - if (_absorptionOnly) { - if (maxT == Ray::infinity()) - return false; - return true; - } else { - Vec3f Tr = std::exp(-_sigmaT*maxT); - Vec3f pdfs = onSurface ? Tr : _sigmaT*Tr; - Vec3f cdf(pdfs.x(), pdfs.x() + pdfs.y(), pdfs.sum()); - - float target = sampler.untracked1D()*cdf.z(); - int component = target < cdf.x() ? 0 : (target < cdf.y() ? 1 : 2); - float Trc = Tr[component]; - - float xi = 1.0f - Trc; - if (onSurface) - xi -= (1.0f - Trc)*sampler.next1D(); - - sampler.putDiscrete(3, component); - sampler.put1D(xi); - } - return true; -} - -Vec3f HomogeneousMedium::transmittance(PathSampleGenerator &/*sampler*/, const Ray &ray) const +Vec3f HomogeneousMedium::transmittance(PathSampleGenerator &/*sampler*/, const Ray &ray, bool startOnSurface, + bool endOnSurface) const { if (ray.farT() == Ray::infinity()) return Vec3f(0.0f); - else { - return FastMath::exp(-_sigmaT*ray.farT()); - } + else + return _transmittance->eval(_sigmaT*ray.farT(), startOnSurface, endOnSurface); } -float HomogeneousMedium::pdf(PathSampleGenerator &/*sampler*/, const Ray &ray, bool onSurface) const +float HomogeneousMedium::pdf(PathSampleGenerator &/*sampler*/, const Ray &ray, bool startOnSurface, bool endOnSurface) const { if (_absorptionOnly) { return 1.0f; } else { - if (onSurface) - return FastMath::exp(-ray.farT()*_sigmaT).avg(); + Vec3f tau = ray.farT()*_sigmaT; + if (endOnSurface) + return _transmittance->surfaceProbability(tau, startOnSurface).avg(); else - return (_sigmaT*FastMath::exp(-ray.farT()*_sigmaT)).avg(); + return (_sigmaT*_transmittance->mediumPdf(tau, startOnSurface)).avg(); } } -Vec3f HomogeneousMedium::transmittanceAndPdfs(PathSampleGenerator &/*sampler*/, const Ray &ray, bool startOnSurface, - bool endOnSurface, float &pdfForward, float &pdfBackward) const -{ - if (ray.farT() == Ray::infinity()) { - pdfForward = pdfBackward = 0.0f; - return Vec3f(0.0f); - } else if (_absorptionOnly) { - pdfForward = pdfBackward = 1.0f; - return FastMath::exp(-_sigmaT*ray.farT()); - } else { - Vec3f weight = FastMath::exp(-_sigmaT*ray.farT()); - pdfForward = endOnSurface ? weight.avg() : (_sigmaT*weight).avg(); - pdfBackward = startOnSurface ? weight.avg() : (_sigmaT*weight).avg(); - return weight; - } -} - -bool HomogeneousMedium::invert(WritablePathSampleGenerator &sampler, const Ray &ray, bool onSurface) const -{ - if (_absorptionOnly) - return true; - - Vec3f transmittance = std::exp(-ray.farT()*_sigmaT); - Vec3f pdfs = _sigmaT*transmittance; - float target = sampler.untracked1D()*pdfs.sum(); - int component = (target < pdfs.x() ? 0 : (target < pdfs.x() + pdfs.y() ? 1 : 2)); - - float xi = 1.0f - transmittance[component]; - if (onSurface) - xi += (1.0f - xi)*sampler.untracked1D(); - sampler.putDiscrete(3, component); - sampler.put1D(xi); - - return true; -} - } diff --git a/src/core/media/HomogeneousMedium.hpp b/src/core/media/HomogeneousMedium.hpp index 72de2b1b..02023e2c 100644 --- a/src/core/media/HomogeneousMedium.hpp +++ b/src/core/media/HomogeneousMedium.hpp @@ -30,12 +30,8 @@ class HomogeneousMedium : public Medium virtual bool sampleDistance(PathSampleGenerator &sampler, const Ray &ray, MediumState &state, MediumSample &sample) const override; - virtual bool invertDistance(WritablePathSampleGenerator &sampler, const Ray &ray, bool onSurface) const; - virtual Vec3f transmittance(PathSampleGenerator &sampler, const Ray &ray) const override; - virtual float pdf(PathSampleGenerator &sampler, const Ray &ray, bool onSurface) const override; - virtual Vec3f transmittanceAndPdfs(PathSampleGenerator &sampler, const Ray &ray, bool startOnSurface, - bool endOnSurface, float &pdfForward, float &pdfBackward) const override; - virtual bool invert(WritablePathSampleGenerator &sampler, const Ray &ray, bool onSurface) const; + virtual Vec3f transmittance(PathSampleGenerator &sampler, const Ray &ray, bool startOnSurface, bool endOnSurface) const override; + virtual float pdf(PathSampleGenerator &sampler, const Ray &ray, bool startOnSurface, bool endOnSurface) const override; Vec3f sigmaA() const { return _sigmaA; } Vec3f sigmaS() const { return _sigmaS; } diff --git a/src/core/media/Medium.cpp b/src/core/media/Medium.cpp index fa6df9f9..8ee0357c 100644 --- a/src/core/media/Medium.cpp +++ b/src/core/media/Medium.cpp @@ -1,5 +1,7 @@ #include "Medium.hpp" +#include "transmittances/ExponentialTransmittance.hpp" + #include "phasefunctions/IsotropicPhaseFunction.hpp" #include "io/JsonObject.hpp" @@ -8,7 +10,8 @@ namespace Tungsten { Medium::Medium() -: _phaseFunction(std::make_shared()), +: _transmittance(std::make_shared()), + _phaseFunction(std::make_shared()), _maxBounce(1024) { } @@ -19,6 +22,8 @@ void Medium::fromJson(JsonPtr value, const Scene &scene) if (auto phase = value["phase_function"]) _phaseFunction = scene.fetchPhase(phase); + if (auto trans = value["transmittance"]) + _transmittance = scene.fetchTransmittance(trans); value.getField("max_bounces", _maxBounce); } @@ -27,6 +32,7 @@ rapidjson::Value Medium::toJson(Allocator &allocator) const { return JsonObject{JsonSerializable::toJson(allocator), allocator, "phase_function", *_phaseFunction, + "transmittance", *_transmittance, "max_bounces", _maxBounce }; } @@ -39,9 +45,9 @@ bool Medium::invertDistance(WritablePathSampleGenerator &/*sampler*/, const Ray Vec3f Medium::transmittanceAndPdfs(PathSampleGenerator &sampler, const Ray &ray, bool startOnSurface, bool endOnSurface, float &pdfForward, float &pdfBackward) const { - pdfForward = pdf(sampler, ray, endOnSurface); - pdfBackward = pdf(sampler, ray.scatter(ray.hitpoint(), -ray.dir(), 0.0f, ray.farT()), startOnSurface); - return transmittance(sampler, ray); + pdfForward = pdf(sampler, ray, startOnSurface, endOnSurface); + pdfBackward = pdf(sampler, ray.scatter(ray.hitpoint(), -ray.dir(), 0.0f, ray.farT()), endOnSurface, startOnSurface); + return transmittance(sampler, ray, startOnSurface, endOnSurface); } const PhaseFunction *Medium::phaseFunction(const Vec3f &/*p*/) const @@ -49,4 +55,9 @@ const PhaseFunction *Medium::phaseFunction(const Vec3f &/*p*/) const return _phaseFunction.get(); } +bool Medium::isDirac() const +{ + return _transmittance->isDirac(); +} + } diff --git a/src/core/media/Medium.hpp b/src/core/media/Medium.hpp index b0e65963..03554656 100644 --- a/src/core/media/Medium.hpp +++ b/src/core/media/Medium.hpp @@ -1,6 +1,8 @@ #ifndef MEDIUM_HPP_ #define MEDIUM_HPP_ +#include "transmittances/Transmittance.hpp" + #include "phasefunctions/PhaseFunction.hpp" #include "samplerecords/MediumSample.hpp" @@ -20,6 +22,7 @@ class Scene; class Medium : public JsonSerializable { protected: + std::shared_ptr _transmittance; std::shared_ptr _phaseFunction; int _maxBounce; @@ -60,11 +63,14 @@ class Medium : public JsonSerializable virtual bool sampleDistance(PathSampleGenerator &sampler, const Ray &ray, MediumState &state, MediumSample &sample) const = 0; virtual bool invertDistance(WritablePathSampleGenerator &sampler, const Ray &ray, bool onSurface) const; - virtual Vec3f transmittance(PathSampleGenerator &sampler, const Ray &ray) const = 0; - virtual float pdf(PathSampleGenerator &sampler, const Ray &ray, bool onSurface) const = 0; + virtual Vec3f transmittance(PathSampleGenerator &sampler, const Ray &ray, bool startOnSurface, + bool endOnSurface) const = 0; + virtual float pdf(PathSampleGenerator &sampler, const Ray &ray, bool startOnSurface, bool endOnSurface) const = 0; virtual Vec3f transmittanceAndPdfs(PathSampleGenerator &sampler, const Ray &ray, bool startOnSurface, bool endOnSurface, float &pdfForward, float &pdfBackward) const; virtual const PhaseFunction *phaseFunction(const Vec3f &p) const; + + bool isDirac() const; }; } diff --git a/src/core/media/VoxelMedium.cpp b/src/core/media/VoxelMedium.cpp index b7f45e23..0144ebd3 100644 --- a/src/core/media/VoxelMedium.cpp +++ b/src/core/media/VoxelMedium.cpp @@ -121,26 +121,28 @@ bool VoxelMedium::sampleDistance(PathSampleGenerator &sampler, const Ray &ray, if (_absorptionOnly) { sample.t = maxT; - sample.weight = _grid->transmittance(sampler, p, w, t0, t1, _sigmaT/wPrime); + Vec3f tau = _grid->opticalDepth(sampler, p, w, t0, t1)*(_sigmaT/wPrime); + sample.weight = _transmittance->eval(tau, state.firstScatter, true); sample.pdf = 1.0f; sample.exited = true; } else { int component = sampler.nextDiscrete(3); float sigmaTc = _sigmaT[component]; - float xi = 1.0f - sampler.next1D(); - float logXi = -std::log(xi); + float tauC = _transmittance->sample(sampler, state.firstScatter)/(sigmaTc/wPrime); - Vec2f tAndDensity = _grid->inverseOpticalDepth(sampler, p, w, t0, t1, sigmaTc/wPrime, logXi); + Vec2f tAndDensity = _grid->inverseOpticalDepth(sampler, p, w, t0, t1, tauC); sample.t = tAndDensity.x(); sample.exited = (sample.t >= t1); - float opticalDepth = sample.exited ? tAndDensity.y()/sigmaTc : logXi/sigmaTc; - sample.weight = std::exp(-_sigmaT*opticalDepth); + if (sample.exited) + tauC = tAndDensity.y(); + Vec3f tau = tauC*(_sigmaT/wPrime); + sample.weight = _transmittance->eval(tau, state.firstScatter, sample.exited); if (sample.exited) { - sample.pdf = sample.weight.avg(); + sample.pdf = _transmittance->surfaceProbability(tau, state.firstScatter).avg(); } else { float rho = tAndDensity.y(); - sample.pdf = (rho*_sigmaT*sample.weight).avg(); - sample.weight *= rho*_sigmaS; + sample.pdf = (rho*_sigmaT*_transmittance->mediumPdf(tau, state.firstScatter)).avg(); + sample.weight *= rho*_sigmaS*_transmittance->sigmaBar(); } sample.weight /= sample.pdf; sample.t /= wPrime; @@ -153,7 +155,8 @@ bool VoxelMedium::sampleDistance(PathSampleGenerator &sampler, const Ray &ray, return true; } -Vec3f VoxelMedium::transmittance(PathSampleGenerator &sampler, const Ray &ray) const +Vec3f VoxelMedium::transmittance(PathSampleGenerator &sampler, const Ray &ray, bool startOnSurface, + bool endOnSurface) const { Vec3f p = _worldToGrid*ray.pos(); Vec3f w = _worldToGrid.transformVector(ray.dir()); @@ -163,10 +166,11 @@ Vec3f VoxelMedium::transmittance(PathSampleGenerator &sampler, const Ray &ray) c if (!bboxIntersection(_gridBounds, p, w, t0, t1)) return Vec3f(1.0f); - return _grid->transmittance(sampler, p, w, t0, t1, _sigmaT/wPrime); + Vec3f tau = _grid->opticalDepth(sampler, p, w, t0, t1)*(_sigmaT/wPrime); + return _transmittance->eval(tau, startOnSurface, endOnSurface); } -float VoxelMedium::pdf(PathSampleGenerator &sampler, const Ray &ray, bool onSurface) const +float VoxelMedium::pdf(PathSampleGenerator &sampler, const Ray &ray, bool startOnSurface, bool endOnSurface) const { if (_absorptionOnly) { return 1.0f; @@ -179,38 +183,12 @@ float VoxelMedium::pdf(PathSampleGenerator &sampler, const Ray &ray, bool onSurf if (!bboxIntersection(_gridBounds, p, w, t0, t1)) return 1.0f; - Vec3f transmittance = _grid->transmittance(sampler, p, w, t0, t1, _sigmaT/wPrime); - if (onSurface) { - return transmittance.avg(); - } else { - return (_grid->density(p)*_sigmaT*transmittance).avg(); - } + Vec3f tau = _grid->opticalDepth(sampler, p, w, t0, t1)*(_sigmaT/wPrime); + if (endOnSurface) + return _transmittance->surfaceProbability(tau, startOnSurface).avg(); + else + return (_grid->density(p)*_sigmaT*_transmittance->mediumPdf(tau, startOnSurface)).avg(); } } -Vec3f VoxelMedium::transmittanceAndPdfs(PathSampleGenerator &sampler, const Ray &ray, bool startOnSurface, - bool endOnSurface, float &pdfForward, float &pdfBackward) const -{ - Vec3f p = _worldToGrid*ray.pos(); - Vec3f w = _worldToGrid.transformVector(ray.dir()); - float wPrime = w.length(); - w /= wPrime; - float t0 = 0.0f, t1 = ray.farT()*wPrime; - if (!bboxIntersection(_gridBounds, p, w, t0, t1)) { - pdfForward = pdfBackward = 1.0f; - return Vec3f(1.0f); - } - - Vec3f transmittance = _grid->transmittance(sampler, p, w, t0, t1, _sigmaT/wPrime); - - if (_absorptionOnly) { - pdfForward = pdfBackward = 1.0f; - } else { - pdfForward = endOnSurface ? transmittance.avg() : (_grid->density(p)*_sigmaT*transmittance).avg(); - pdfBackward = startOnSurface ? transmittance.avg() : (_grid->density(p)*_sigmaT*transmittance).avg(); - } - - return transmittance; -} - } diff --git a/src/core/media/VoxelMedium.hpp b/src/core/media/VoxelMedium.hpp index fbfb436f..f7657ca5 100644 --- a/src/core/media/VoxelMedium.hpp +++ b/src/core/media/VoxelMedium.hpp @@ -36,10 +36,9 @@ class VoxelMedium : public Medium virtual bool sampleDistance(PathSampleGenerator &sampler, const Ray &ray, MediumState &state, MediumSample &sample) const override; - virtual Vec3f transmittance(PathSampleGenerator &sampler, const Ray &ray) const override; - virtual float pdf(PathSampleGenerator &sampler, const Ray &ray, bool onSurface) const override; - virtual Vec3f transmittanceAndPdfs(PathSampleGenerator &sampler, const Ray &ray, bool startOnSurface, - bool endOnSurface, float &pdfForward, float &pdfBackward) const override; + virtual Vec3f transmittance(PathSampleGenerator &sampler, const Ray &ray, bool startOnSurface, + bool endOnSurface) const override; + virtual float pdf(PathSampleGenerator &sampler, const Ray &ray, bool startOnSurface, bool endOnSurface) const override; }; } diff --git a/src/core/transmittances/DavisTransmittance.cpp b/src/core/transmittances/DavisTransmittance.cpp new file mode 100644 index 00000000..a5d50e01 --- /dev/null +++ b/src/core/transmittances/DavisTransmittance.cpp @@ -0,0 +1,62 @@ +#include "DavisTransmittance.hpp" + +#include "io/JsonObject.hpp" + +namespace Tungsten { + +DavisTransmittance::DavisTransmittance() +: _alpha(1.1f) +{ +} + +void DavisTransmittance::fromJson(JsonPtr value, const Scene &scene) +{ + Transmittance::fromJson(value, scene); + value.getField("alpha", _alpha); + + if (_alpha < 1 + 1e-6f) { + std::cout << "Warning: alpha parameter of Davis transmittance has to be > 1. Clamping the current value (" << _alpha << ")" << std::endl; + _alpha = 1 + 1e-6f; + } +} + +rapidjson::Value DavisTransmittance::toJson(Allocator &allocator) const +{ + return JsonObject{Transmittance::toJson(allocator), allocator, + "type", "davis", + "alpha", _alpha + }; +} + +Vec3f DavisTransmittance::surfaceSurface(const Vec3f &tau) const +{ + return std::pow(1.0f + tau/_alpha, -_alpha); +} +Vec3f DavisTransmittance::surfaceMedium(const Vec3f &tau) const +{ + return std::pow(1.0f + tau/_alpha, -(_alpha + 1.0f)); +} +Vec3f DavisTransmittance::mediumSurface(const Vec3f &tau) const +{ + return surfaceMedium(tau); +} +Vec3f DavisTransmittance::mediumMedium(const Vec3f &tau) const +{ + return (1.0f + 1.0f/_alpha)*std::pow(1.0f + tau/_alpha, -(_alpha + 2.0f)); +} + +float DavisTransmittance::sigmaBar() const +{ + return 1.0f; +} + +float DavisTransmittance::sampleSurface(PathSampleGenerator &sampler) const +{ + return _alpha*(std::pow(1.0f - sampler.next1D(), -1.0f/_alpha) - 1.0f); +} +float DavisTransmittance::sampleMedium(PathSampleGenerator &sampler) const +{ + return _alpha*(std::pow(1.0f - sampler.next1D(), -1.0f/(1.0f + _alpha)) - 1.0f); +} + +} diff --git a/src/core/transmittances/DavisTransmittance.hpp b/src/core/transmittances/DavisTransmittance.hpp new file mode 100644 index 00000000..8b11dd38 --- /dev/null +++ b/src/core/transmittances/DavisTransmittance.hpp @@ -0,0 +1,31 @@ +#ifndef DAVISTRANSMITTANCE_HPP_ +#define DAVISTRANSMITTANCE_HPP_ + +#include "Transmittance.hpp" + +namespace Tungsten { + +class DavisTransmittance : public Transmittance +{ + float _alpha; + +public: + DavisTransmittance(); + + virtual void fromJson(JsonPtr value, const Scene &scene) override; + virtual rapidjson::Value toJson(Allocator &allocator) const override; + + virtual Vec3f surfaceSurface(const Vec3f &tau) const override final; + virtual Vec3f surfaceMedium(const Vec3f &tau) const override final; + virtual Vec3f mediumSurface(const Vec3f &tau) const override final; + virtual Vec3f mediumMedium(const Vec3f &tau) const override final; + + virtual float sigmaBar() const override final; + + virtual float sampleSurface(PathSampleGenerator &sampler) const override final; + virtual float sampleMedium(PathSampleGenerator &sampler) const override final; +}; + +} + +#endif /* DAVISTRANSMITTANCE_HPP_ */ diff --git a/src/core/transmittances/DavisWeinsteinTransmittance.cpp b/src/core/transmittances/DavisWeinsteinTransmittance.cpp new file mode 100644 index 00000000..ac8524a3 --- /dev/null +++ b/src/core/transmittances/DavisWeinsteinTransmittance.cpp @@ -0,0 +1,119 @@ +#include "DavisWeinsteinTransmittance.hpp" + +#include "math/MathUtil.hpp" + +#include "io/JsonObject.hpp" + +namespace Tungsten { + +DavisWeinsteinTransmittance::DavisWeinsteinTransmittance() +: _h(0.75f), + _c(1.0f) +{ +} + +void DavisWeinsteinTransmittance::fromJson(JsonPtr value, const Scene &scene) +{ + Transmittance::fromJson(value, scene); + value.getField("h", _h); + value.getField("c", _c); + if (_h < 0.5f || _h > 1.0f) + std::cout << "Warning: Valid range of Davis Weinstein transmittance is [0.5, 1.0]. Clamping current value (" << _h << ") to within range" << std::endl; + _h = clamp(_h, 0.5f, 1.0f); +} + +rapidjson::Value DavisWeinsteinTransmittance::toJson(Allocator &allocator) const +{ + return JsonObject{Transmittance::toJson(allocator), allocator, + "type", "davis", + "h", _h, + "c", _c + }; +} + +float DavisWeinsteinTransmittance::computeAlpha(float tau) const +{ + float beta = 2.0f*_h - 1.0f; + return std::pow(tau, 1 - beta)/(std::pow(_c, 1 + beta)); +} + +Vec3f DavisWeinsteinTransmittance::surfaceSurface(const Vec3f &tau) const +{ + float alpha = computeAlpha(tau[0]); + float Tr = std::pow(1.0f + tau[0]/alpha, -alpha); + + return Vec3f(std::isnan(Tr) ? 0 : Tr); +} +Vec3f DavisWeinsteinTransmittance::surfaceMedium(const Vec3f &tau) const +{ + float beta = 2.0f*_h - 1.0f; + float t = tau[0]; + float alpha = computeAlpha(t); + float base = 1.0f + t/alpha; + float trSurface = std::pow(base, -alpha); + + float Tr = trSurface*(beta/base - (beta - 1.0f)*alpha/t*std::log(base)); + + return Vec3f(std::isnan(Tr) ? 0 : Tr); +} +Vec3f DavisWeinsteinTransmittance::mediumSurface(const Vec3f &tau) const +{ + return surfaceMedium(tau); +} +Vec3f DavisWeinsteinTransmittance::mediumMedium(const Vec3f &tau) const +{ + float beta = 2.0f*_h - 1.0f; + float t = tau[0]; + float alpha = computeAlpha(t); + float base = 1.0f + t/alpha; + float logBase = std::log(base); + float trSurface = std::pow(base, -alpha); + float term1 = beta*(-1.0f + beta*(1.0f + t) + (-1.0f + 2.0f*beta)*t/alpha)/(t*base*base); + float term2 = ((-1.0f + beta)*beta*alpha/(t*t)*(2.0f*t + base)*logBase)/base; + float term3 = (beta - 1.0f)*alpha/t*logBase; + + float Tr = trSurface*(term1 - term2 + term3*term3); + + return Vec3f(std::isnan(Tr) ? 0 : Tr); +} + +float DavisWeinsteinTransmittance::sigmaBar() const +{ + return 1.0f; +} + +// We have no real way of analytically sampling this function, so a simple bisection will have to do +float DavisWeinsteinTransmittance::sampleSurface(PathSampleGenerator &sampler) const +{ + float xi = sampler.next1D(); + auto cdf = [this](float tau) { return 1.0f - surfaceSurface(Vec3f(tau))[0]; }; + float step = 1e6; + float result = step*2; + while (step > 1e-6) { + if (cdf(result) > xi) + result -= step; + else + result += step; + step /= 2; + } + + return result; +} +float DavisWeinsteinTransmittance::sampleMedium(PathSampleGenerator &sampler) const +{ + float xi = sampler.next1D(); + auto cdf = [this](float tau) { return 1.0f - mediumSurface(Vec3f(tau))[0]; }; + float step = 1e6; + float result = step * 2; + while (step > 1e-6) { + if (cdf(result) > xi) + result -= step; + else + result += step; + step /= 2; + } + + return result; +} + +} diff --git a/src/core/transmittances/DavisWeinsteinTransmittance.hpp b/src/core/transmittances/DavisWeinsteinTransmittance.hpp new file mode 100644 index 00000000..cddf03f0 --- /dev/null +++ b/src/core/transmittances/DavisWeinsteinTransmittance.hpp @@ -0,0 +1,34 @@ +#ifndef DAVISWEINSTEINTRANSMITTANCE_HPP_ +#define DAVISWEINSTEINTRANSMITTANCE_HPP_ + +#include "Transmittance.hpp" + +namespace Tungsten { + +class DavisWeinsteinTransmittance : public Transmittance +{ + float _h; + float _c; + + float computeAlpha(float tau) const; + +public: + DavisWeinsteinTransmittance(); + + virtual void fromJson(JsonPtr value, const Scene &scene) override; + virtual rapidjson::Value toJson(Allocator &allocator) const override; + + virtual Vec3f surfaceSurface(const Vec3f &tau) const override final; + virtual Vec3f surfaceMedium(const Vec3f &tau) const override final; + virtual Vec3f mediumSurface(const Vec3f &tau) const override final; + virtual Vec3f mediumMedium(const Vec3f &tau) const override final; + + virtual float sigmaBar() const override final; + + virtual float sampleSurface(PathSampleGenerator &sampler) const override final; + virtual float sampleMedium(PathSampleGenerator &sampler) const override final; +}; + +} + +#endif /* DAVISWEINSTEINTRANSMITTANCE_HPP_ */ diff --git a/src/core/transmittances/DoubleExponentialTransmittance.cpp b/src/core/transmittances/DoubleExponentialTransmittance.cpp new file mode 100644 index 00000000..3bdca5eb --- /dev/null +++ b/src/core/transmittances/DoubleExponentialTransmittance.cpp @@ -0,0 +1,68 @@ +#include "DoubleExponentialTransmittance.hpp" + +#include "sampling/UniformPathSampler.hpp" + +#include "math/MathUtil.hpp" + +#include "io/JsonObject.hpp" + +#include "Memory.hpp" + +namespace Tungsten { + +DoubleExponentialTransmittance::DoubleExponentialTransmittance() +: _sigmaA(0.5f), + _sigmaB(10.0f) +{ +} + +void DoubleExponentialTransmittance::fromJson(JsonPtr value, const Scene &scene) +{ + Transmittance::fromJson(value, scene); + value.getField("sigma_a", _sigmaA); + value.getField("sigma_b", _sigmaB); +} + +rapidjson::Value DoubleExponentialTransmittance::toJson(Allocator &allocator) const +{ + return JsonObject{Transmittance::toJson(allocator), allocator, + "type", "double_exponential", + "sigma_a", _sigmaA, + "sigma_b", _sigmaB, + }; +} + +Vec3f DoubleExponentialTransmittance::surfaceSurface(const Vec3f &tau) const +{ + return 0.5f*(std::exp(-_sigmaA*tau) + std::exp(-_sigmaB*tau)); +} +Vec3f DoubleExponentialTransmittance::surfaceMedium(const Vec3f &tau) const +{ + return 0.5f*(_sigmaA*std::exp(-_sigmaA*tau) + _sigmaB*std::exp(-_sigmaB*tau)); +} +Vec3f DoubleExponentialTransmittance::mediumSurface(const Vec3f &tau) const +{ + return (_sigmaA*std::exp(-_sigmaA*tau) + _sigmaB*std::exp(-_sigmaB*tau))/(_sigmaA + _sigmaB); +} +Vec3f DoubleExponentialTransmittance::mediumMedium(const Vec3f &tau) const +{ + return (sqr(_sigmaA)*std::exp(-_sigmaA*tau) + sqr(_sigmaB)*std::exp(-_sigmaB*tau))/(_sigmaA + _sigmaB); +} + +float DoubleExponentialTransmittance::sigmaBar() const +{ + return 0.5f*(_sigmaA + _sigmaB); +} + +float DoubleExponentialTransmittance::sampleSurface(PathSampleGenerator &sampler) const +{ + float t = -std::log(1.0f - sampler.next1D()); + return sampler.nextBoolean(0.5f) ? t/_sigmaA : t/_sigmaB; +} +float DoubleExponentialTransmittance::sampleMedium(PathSampleGenerator &sampler) const +{ + float t = -std::log(1.0f - sampler.next1D()); + return sampler.nextBoolean(_sigmaA/(_sigmaA + _sigmaB)) ? t/_sigmaA : t/_sigmaB; +} + +} diff --git a/src/core/transmittances/DoubleExponentialTransmittance.hpp b/src/core/transmittances/DoubleExponentialTransmittance.hpp new file mode 100644 index 00000000..9f414ffd --- /dev/null +++ b/src/core/transmittances/DoubleExponentialTransmittance.hpp @@ -0,0 +1,35 @@ +#ifndef DoubleExponentialTransmittance_HPP_ +#define DoubleExponentialTransmittance_HPP_ + +#include "Transmittance.hpp" + +namespace Tungsten { + +class DoubleExponentialTransmittance : public Transmittance +{ + float _sigmaA; + float _sigmaB; + +public: + DoubleExponentialTransmittance(); + + virtual void fromJson(JsonPtr value, const Scene &scene) override; + virtual rapidjson::Value toJson(Allocator &allocator) const override; + + virtual Vec3f surfaceSurface(const Vec3f &tau) const override final; + virtual Vec3f surfaceMedium(const Vec3f &tau) const override final; + virtual Vec3f mediumSurface(const Vec3f &tau) const override final; + virtual Vec3f mediumMedium(const Vec3f &tau) const override final; + + virtual float sigmaBar() const override final; + + virtual float sampleSurface(PathSampleGenerator &sampler) const override final; + virtual float sampleMedium(PathSampleGenerator &sampler) const override final; +}; + + +} + + + +#endif /* DoubleExponentialTransmittance_HPP_ */ diff --git a/src/core/transmittances/ErlangTransmittance.cpp b/src/core/transmittances/ErlangTransmittance.cpp new file mode 100644 index 00000000..df20da51 --- /dev/null +++ b/src/core/transmittances/ErlangTransmittance.cpp @@ -0,0 +1,69 @@ +#include "ErlangTransmittance.hpp" + +#include "sampling/UniformPathSampler.hpp" + +#include "math/MathUtil.hpp" + +#include "io/JsonObject.hpp" + +#include "Memory.hpp" + +namespace Tungsten { + +ErlangTransmittance::ErlangTransmittance() +: _lambda(5.0f) +{ +} + +void ErlangTransmittance::fromJson(JsonPtr value, const Scene &scene) +{ + Transmittance::fromJson(value, scene); + value.getField("rate", _lambda); +} + +rapidjson::Value ErlangTransmittance::toJson(Allocator &allocator) const +{ + return JsonObject{Transmittance::toJson(allocator), allocator, + "type", "erlang", + "rate", _lambda, + }; +} + +Vec3f ErlangTransmittance::surfaceSurface(const Vec3f &tau) const +{ + return 0.5f*std::exp(-_lambda*tau)*(2.0f + _lambda*tau); +} +Vec3f ErlangTransmittance::surfaceMedium(const Vec3f &tau) const +{ + return mediumSurface(tau)*_lambda*0.5f; +} +Vec3f ErlangTransmittance::mediumSurface(const Vec3f &tau) const +{ + return std::exp(-_lambda*tau)*(1.0f + _lambda*tau); +} +Vec3f ErlangTransmittance::mediumMedium(const Vec3f &tau) const +{ + return sqr(_lambda)*tau*std::exp(-_lambda*tau); +} + +float ErlangTransmittance::sigmaBar() const +{ + return _lambda*0.5f; +} + +float ErlangTransmittance::sampleSurface(PathSampleGenerator &sampler) const +{ + float xi = sampler.next1D(); + float x = 0.5f; + for (int i = 0; i < 10; ++i) { + x += (xi - (1.0f - surfaceSurface(Vec3f(x))[0]))/surfaceMedium(Vec3f(x))[0]; + x = max(x, 0.0f); + } + return x; +} +float ErlangTransmittance::sampleMedium(PathSampleGenerator &sampler) const +{ + return -1.0f/_lambda*std::log(sampler.next1D()*sampler.next1D()); +} + +} diff --git a/src/core/transmittances/ErlangTransmittance.hpp b/src/core/transmittances/ErlangTransmittance.hpp new file mode 100644 index 00000000..3d62decb --- /dev/null +++ b/src/core/transmittances/ErlangTransmittance.hpp @@ -0,0 +1,34 @@ +#ifndef ErlangTransmittance_HPP_ +#define ErlangTransmittance_HPP_ + +#include "Transmittance.hpp" + +namespace Tungsten { + +class ErlangTransmittance : public Transmittance +{ + float _lambda; + +public: + ErlangTransmittance(); + + virtual void fromJson(JsonPtr value, const Scene &scene) override; + virtual rapidjson::Value toJson(Allocator &allocator) const override; + + virtual Vec3f surfaceSurface(const Vec3f &tau) const override final; + virtual Vec3f surfaceMedium(const Vec3f &tau) const override final; + virtual Vec3f mediumSurface(const Vec3f &tau) const override final; + virtual Vec3f mediumMedium(const Vec3f &tau) const override final; + + virtual float sigmaBar() const override final; + + virtual float sampleSurface(PathSampleGenerator &sampler) const override final; + virtual float sampleMedium(PathSampleGenerator &sampler) const override final; +}; + + +} + + + +#endif /* ErlangTransmittance_HPP_ */ diff --git a/src/core/transmittances/ExponentialTransmittance.cpp b/src/core/transmittances/ExponentialTransmittance.cpp new file mode 100644 index 00000000..38d44ef5 --- /dev/null +++ b/src/core/transmittances/ExponentialTransmittance.cpp @@ -0,0 +1,57 @@ +#include "ExponentialTransmittance.hpp" + +#include "sampling/UniformPathSampler.hpp" + +#include "math/MathUtil.hpp" +#include "math/FastMath.hpp" + +#include "io/JsonObject.hpp" + +#include "Memory.hpp" + +namespace Tungsten { + +void ExponentialTransmittance::fromJson(JsonPtr value, const Scene &scene) +{ + Transmittance::fromJson(value, scene); +} + +rapidjson::Value ExponentialTransmittance::toJson(Allocator &allocator) const +{ + return JsonObject{Transmittance::toJson(allocator), allocator, + "type", "exponential", + }; +} + +Vec3f ExponentialTransmittance::surfaceSurface(const Vec3f &tau) const +{ + return FastMath::exp(-tau); +} +Vec3f ExponentialTransmittance::surfaceMedium(const Vec3f &tau) const +{ + return FastMath::exp(-tau); +} +Vec3f ExponentialTransmittance::mediumSurface(const Vec3f &tau) const +{ + return FastMath::exp(-tau); +} +Vec3f ExponentialTransmittance::mediumMedium(const Vec3f &tau) const +{ + return FastMath::exp(-tau); +} + +float ExponentialTransmittance::sigmaBar() const +{ + return 1.0f; +} + +float ExponentialTransmittance::sampleSurface(PathSampleGenerator &sampler) const +{ + return -std::log(1.0f - sampler.next1D()); +} +float ExponentialTransmittance::sampleMedium(PathSampleGenerator &sampler) const +{ + return -std::log(1.0f - sampler.next1D()); +} + +} diff --git a/src/core/transmittances/ExponentialTransmittance.hpp b/src/core/transmittances/ExponentialTransmittance.hpp new file mode 100644 index 00000000..a1731d2c --- /dev/null +++ b/src/core/transmittances/ExponentialTransmittance.hpp @@ -0,0 +1,32 @@ +#ifndef EXPONENTIALTRANSMITTANCE_HPP_ +#define EXPONENTIALTRANSMITTANCE_HPP_ + +#include "Transmittance.hpp" + +namespace Tungsten { + +class ExponentialTransmittance : public Transmittance +{ +public: + ExponentialTransmittance() = default; + + virtual void fromJson(JsonPtr value, const Scene &scene) override; + virtual rapidjson::Value toJson(Allocator &allocator) const override; + + virtual Vec3f surfaceSurface(const Vec3f &tau) const override final; + virtual Vec3f surfaceMedium(const Vec3f &tau) const override final; + virtual Vec3f mediumSurface(const Vec3f &tau) const override final; + virtual Vec3f mediumMedium(const Vec3f &tau) const override final; + + virtual float sigmaBar() const override final; + + virtual float sampleSurface(PathSampleGenerator &sampler) const override final; + virtual float sampleMedium(PathSampleGenerator &sampler) const override final; +}; + + +} + + + +#endif /* EXPONENTIALTRANSMITTANCE_HPP_ */ diff --git a/src/core/transmittances/InterpolatedTransmittance.cpp b/src/core/transmittances/InterpolatedTransmittance.cpp new file mode 100644 index 00000000..686f3748 --- /dev/null +++ b/src/core/transmittances/InterpolatedTransmittance.cpp @@ -0,0 +1,85 @@ +#include "InterpolatedTransmittance.hpp" + +#include "ExponentialTransmittance.hpp" +#include "LinearTransmittance.hpp" +#include "ErlangTransmittance.hpp" + +#include "sampling/UniformPathSampler.hpp" + +#include "math/MathUtil.hpp" + +#include "io/JsonObject.hpp" +#include "io/Scene.hpp" + +#include "Memory.hpp" + +namespace Tungsten { + +InterpolatedTransmittance::InterpolatedTransmittance() +: _trA(std::make_shared()), + _trB(std::make_shared()), + _u(0.5f) +{ +} + +void InterpolatedTransmittance::fromJson(JsonPtr value, const Scene &scene) +{ + Transmittance::fromJson(value, scene); + if (auto m = value["tr_a"]) _trA = scene.fetchTransmittance(m); + if (auto m = value["tr_b"]) _trB = scene.fetchTransmittance(m); + value.getField("ratio", _u); +} + +rapidjson::Value InterpolatedTransmittance::toJson(Allocator &allocator) const +{ + return JsonObject{Transmittance::toJson(allocator), allocator, + "type", "interpolated", + "tr_a", _trA->toJson(allocator), + "tr_b", _trB->toJson(allocator), + "ratio", _u + }; +} + +Vec3f InterpolatedTransmittance::surfaceSurface(const Vec3f &tau) const +{ + return sigmaBar()*lerp(_trA->surfaceSurface(tau)/_trA->sigmaBar(), _trB->surfaceSurface(tau)/_trB->sigmaBar(), _u); +} +Vec3f InterpolatedTransmittance::surfaceMedium(const Vec3f &tau) const +{ + return mediumSurface(tau)*sigmaBar(); +} +Vec3f InterpolatedTransmittance::mediumSurface(const Vec3f &tau) const +{ + return lerp(_trA->mediumSurface(tau), _trB->mediumSurface(tau), _u); +} +Vec3f InterpolatedTransmittance::mediumMedium(const Vec3f &tau) const +{ + Vec3f pdfA = _trA->mediumMedium(tau); + Vec3f pdfB = _trB->mediumMedium(tau); + Vec3f result; + for (int i = 0; i < 3; ++i) { + bool diracA = _trA->isDirac() && pdfA[i] > 0.0f; + bool diracB = _trB->isDirac() && pdfB[i] > 0.0f; + if (diracA ^ diracB) + result[i] = diracA ? pdfA[i] : pdfB[i]; + else + result[i] = lerp(pdfA[i], pdfB[i], _u); + } + return result; +} + +float InterpolatedTransmittance::sigmaBar() const +{ + return 1.0f/lerp(1.0f/_trA->sigmaBar(), 1.0f/_trB->sigmaBar(), _u); +} + +float InterpolatedTransmittance::sampleSurface(PathSampleGenerator &sampler) const +{ + return sampler.nextBoolean(_u) ? _trB->sampleSurface(sampler) : _trA->sampleSurface(sampler); +} +float InterpolatedTransmittance::sampleMedium(PathSampleGenerator &sampler) const +{ + return sampler.nextBoolean(_u) ? _trB->sampleMedium(sampler) : _trA->sampleMedium(sampler); +} + +} diff --git a/src/core/transmittances/InterpolatedTransmittance.hpp b/src/core/transmittances/InterpolatedTransmittance.hpp new file mode 100644 index 00000000..86f8550a --- /dev/null +++ b/src/core/transmittances/InterpolatedTransmittance.hpp @@ -0,0 +1,39 @@ +#ifndef INTERPOLATEDTRANSMITTANCE_HPP_ +#define INTERPOLATEDTRANSMITTANCE_HPP_ + +#include "Transmittance.hpp" + +#include +#include + +namespace Tungsten { + +class InterpolatedTransmittance : public Transmittance +{ + std::shared_ptr _trA; + std::shared_ptr _trB; + float _u; + +public: + InterpolatedTransmittance(); + + virtual void fromJson(JsonPtr value, const Scene &scene) override; + virtual rapidjson::Value toJson(Allocator &allocator) const override; + + virtual Vec3f surfaceSurface(const Vec3f &tau) const override final; + virtual Vec3f surfaceMedium(const Vec3f &tau) const override final; + virtual Vec3f mediumSurface(const Vec3f &tau) const override final; + virtual Vec3f mediumMedium(const Vec3f &tau) const override final; + + virtual float sigmaBar() const override final; + + virtual float sampleSurface(PathSampleGenerator &sampler) const override final; + virtual float sampleMedium(PathSampleGenerator &sampler) const override final; +}; + + +} + + + +#endif /* INTERPOLATEDTRANSMITTANCE_HPP_ */ diff --git a/src/core/transmittances/LinearTransmittance.cpp b/src/core/transmittances/LinearTransmittance.cpp new file mode 100644 index 00000000..62b8a9a0 --- /dev/null +++ b/src/core/transmittances/LinearTransmittance.cpp @@ -0,0 +1,79 @@ +#include "LinearTransmittance.hpp" + +#include "sampling/UniformPathSampler.hpp" + +#include "math/MathUtil.hpp" + +#include "io/JsonObject.hpp" + +#include "Memory.hpp" + +namespace Tungsten { + +LinearTransmittance::LinearTransmittance() +: _maxT(1.0f) +{ +} + +void LinearTransmittance::fromJson(JsonPtr value, const Scene &scene) +{ + Transmittance::fromJson(value, scene); + value.getField("max_t", _maxT); +} + +rapidjson::Value LinearTransmittance::toJson(Allocator &allocator) const +{ + return JsonObject{Transmittance::toJson(allocator), allocator, + "type", "linear", + "max_t", _maxT + }; +} + +Vec3f LinearTransmittance::surfaceSurface(const Vec3f &tau) const +{ + return 1.0f - min(tau/_maxT, Vec3f(1.0f)); +} +Vec3f LinearTransmittance::surfaceMedium(const Vec3f &tau) const +{ + Vec3f result(1.0f/_maxT); + for (int i = 0; i < 3; ++i) + if (tau[i] > _maxT) + result[i] = 0.0f; + return result; +} +Vec3f LinearTransmittance::mediumSurface(const Vec3f &tau) const +{ + Vec3f result(1.0f); + for (int i = 0; i < 3; ++i) + if (tau[i] > _maxT) + result[i] = 0.0f; + return result; +} +Vec3f LinearTransmittance::mediumMedium(const Vec3f &tau) const +{ + Vec3f result; + for (int i = 0; i < 3; ++i) + result[i] = std::abs(tau[i] - _maxT) < 1e-3f ? 1.0f : 0.0f; + return result; +} + +float LinearTransmittance::sigmaBar() const +{ + return 1.0f/_maxT; +} + +bool LinearTransmittance::isDirac() const +{ + return true; +} + +float LinearTransmittance::sampleSurface(PathSampleGenerator &sampler) const +{ + return _maxT*sampler.next1D(); +} +float LinearTransmittance::sampleMedium(PathSampleGenerator &/*sampler*/) const +{ + return _maxT; +} + +} diff --git a/src/core/transmittances/LinearTransmittance.hpp b/src/core/transmittances/LinearTransmittance.hpp new file mode 100644 index 00000000..e3458e4c --- /dev/null +++ b/src/core/transmittances/LinearTransmittance.hpp @@ -0,0 +1,33 @@ +#ifndef LINEARTRANSMITTANCE_HPP_ +#define LINEARTRANSMITTANCE_HPP_ + +#include "Transmittance.hpp" + +namespace Tungsten { + +class LinearTransmittance : public Transmittance +{ + float _maxT; + +public: + LinearTransmittance(); + + virtual void fromJson(JsonPtr value, const Scene &scene) override; + virtual rapidjson::Value toJson(Allocator &allocator) const override; + + virtual Vec3f surfaceSurface(const Vec3f &tau) const override final; + virtual Vec3f surfaceMedium(const Vec3f &tau) const override final; + virtual Vec3f mediumSurface(const Vec3f &tau) const override final; + virtual Vec3f mediumMedium(const Vec3f &tau) const override final; + + virtual bool isDirac() const override final; + + virtual float sigmaBar() const override final; + + virtual float sampleSurface(PathSampleGenerator &sampler) const override final; + virtual float sampleMedium(PathSampleGenerator &sampler) const override final; +}; + +} + +#endif /* LINEARTRANSMITTANCE_HPP_ */ diff --git a/src/core/transmittances/PulseTransmittance.cpp b/src/core/transmittances/PulseTransmittance.cpp new file mode 100644 index 00000000..50eee130 --- /dev/null +++ b/src/core/transmittances/PulseTransmittance.cpp @@ -0,0 +1,108 @@ +#include "PulseTransmittance.hpp" + +#include "sampling/UniformPathSampler.hpp" + +#include "math/MathUtil.hpp" + +#include "io/JsonObject.hpp" +#include "io/FileUtils.hpp" + +#include "Memory.hpp" + +namespace Tungsten { + +PulseTransmittance::PulseTransmittance() +: _a(0.0f), + _b(1.0f), + _numPulses(4) +{ +} + +void PulseTransmittance::fromJson(JsonPtr value, const Scene &scene) +{ + Transmittance::fromJson(value, scene); + value.getField("min", _a); + value.getField("max", _b); + value.getField("num_pulses", _numPulses); +} + +rapidjson::Value PulseTransmittance::toJson(Allocator &allocator) const +{ + return JsonObject{Transmittance::toJson(allocator), allocator, + "type", "pulse", + "min", _a, + "max", _b, + "num_pulses", _numPulses, + }; +} + +bool PulseTransmittance::isDirac() const +{ + return true; +} + +Vec3f PulseTransmittance::surfaceSurface(const Vec3f &tau) const +{ + Vec3f idxF = clamp(float(_numPulses)*(tau - _a)/(_b - _a) + 0.5f, Vec3f(0.0f), Vec3f(_numPulses)); + Vec3i idx = Vec3i(idxF); + + Vec3f height = Vec3f(_numPulses - idx)/float(_numPulses); + + Vec3f cellIntegral = height*(idxF - Vec3f(idx)); + for (int i = 0; i < 3; ++i) { + if (idx[i] > 0) + cellIntegral[i] += (idx[i] - 0.5f) - (idx[i]*(idx[i] - 1))/(2.0f*_numPulses); + else + cellIntegral[i] -= 0.5f; + } + + return 1.0f - (2.0f/_numPulses)*cellIntegral; +} +Vec3f PulseTransmittance::surfaceMedium(const Vec3f &tau) const +{ + return 2.0f/(_b - _a)*mediumSurface(tau); +} +Vec3f PulseTransmittance::mediumSurface(const Vec3f &tau) const +{ + Vec3i idx = clamp(Vec3i(float(_numPulses)*(tau - _a)/(_b - _a) + 0.5f), Vec3i(0), Vec3i(_numPulses)); + return 1.0f - Vec3f(idx)/float(_numPulses); +} +Vec3f PulseTransmittance::mediumMedium(const Vec3f &tau) const +{ + Vec3f idxF = clamp(float(_numPulses)*(tau - _a)/(_b - _a), Vec3f(0.0f), Vec3f(_numPulses)); + Vec3i idx = Vec3i(idxF); + return (1.0f/_numPulses)*Vec3f( + std::abs(idxF[0] - idx[0] - 0.5f) < 1e-3f ? 1.0f : 0.0f, + std::abs(idxF[1] - idx[1] - 0.5f) < 1e-3f ? 1.0f : 0.0f, + std::abs(idxF[2] - idx[2] - 0.5f) < 1e-3f ? 1.0f : 0.0f + ); +} + +float PulseTransmittance::sigmaBar() const +{ + return 2.0f/(_b - _a); +} + +float PulseTransmittance::sampleSurface(PathSampleGenerator &sampler) const +{ + float xi = sampler.next1D()*_numPulses*0.5f; + float delta = 1.0f/_numPulses; + + for (int i = 0; i < _numPulses; ++i) { + float h0 = 1.0f - (i + 0.0f)*delta; + float h1 = 1.0f - (i + 1.0f)*delta; + xi -= h0*0.5f; + if (xi < 0.0f) + return _a + (i + 0.0f + 0.5f*sampler.next1D())*(_b - _a)*delta; + xi -= h1*0.5f; + if (xi < 0.0f) + return _a + (i + 0.5f + 0.5f*sampler.next1D())*(_b - _a)*delta; + } + return 0.0f; +} +float PulseTransmittance::sampleMedium(PathSampleGenerator &sampler) const +{ + return _a + (0.5f + int(sampler.next1D()*_numPulses))/float(_numPulses)*(_b - _a); +} + +} diff --git a/src/core/transmittances/PulseTransmittance.hpp b/src/core/transmittances/PulseTransmittance.hpp new file mode 100644 index 00000000..76e18cee --- /dev/null +++ b/src/core/transmittances/PulseTransmittance.hpp @@ -0,0 +1,37 @@ +#ifndef PULSETRANSMITTANCE_HPP_ +#define PULSETRANSMITTANCE_HPP_ + +#include "Transmittance.hpp" + +namespace Tungsten { + +class PulseTransmittance : public Transmittance +{ + float _a, _b; + int _numPulses; + +public: + PulseTransmittance(); + + virtual void fromJson(JsonPtr value, const Scene &scene) override; + virtual rapidjson::Value toJson(Allocator &allocator) const override; + + virtual bool isDirac() const override final; + + virtual Vec3f surfaceSurface(const Vec3f &tau) const override final; + virtual Vec3f surfaceMedium(const Vec3f &tau) const override final; + virtual Vec3f mediumSurface(const Vec3f &tau) const override final; + virtual Vec3f mediumMedium(const Vec3f &tau) const override final; + + virtual float sigmaBar() const override final; + + virtual float sampleSurface(PathSampleGenerator &sampler) const override final; + virtual float sampleMedium(PathSampleGenerator &sampler) const override final; +}; + + +} + + + +#endif /* PULSETRANSMITTANCE_HPP_ */ diff --git a/src/core/transmittances/QuadraticTransmittance.cpp b/src/core/transmittances/QuadraticTransmittance.cpp new file mode 100644 index 00000000..8874f3f0 --- /dev/null +++ b/src/core/transmittances/QuadraticTransmittance.cpp @@ -0,0 +1,68 @@ +#include "QuadraticTransmittance.hpp" + +#include "sampling/UniformPathSampler.hpp" + +#include "math/MathUtil.hpp" + +#include "io/JsonObject.hpp" + +#include "Memory.hpp" + +namespace Tungsten { + +QuadraticTransmittance::QuadraticTransmittance() +: _maxT(0.75f) +{ +} + +void QuadraticTransmittance::fromJson(JsonPtr value, const Scene &scene) +{ + Transmittance::fromJson(value, scene); + value.getField("max_t", _maxT); +} + +rapidjson::Value QuadraticTransmittance::toJson(Allocator &allocator) const +{ + return JsonObject{Transmittance::toJson(allocator), allocator, + "type", "quadratic", + "max_t", _maxT + }; +} + +Vec3f QuadraticTransmittance::surfaceSurface(const Vec3f &tau) const +{ + Vec3f t = min(tau/_maxT, Vec3f(1.0f)); + return 1.0f - 2.0f*t + t*t; +} +Vec3f QuadraticTransmittance::surfaceMedium(const Vec3f &tau) const +{ + return (2.0f/_maxT)*(1.0f - min(tau/_maxT, Vec3f(1.0f))); +} +Vec3f QuadraticTransmittance::mediumSurface(const Vec3f &tau) const +{ + return 1.0f - min(tau/_maxT, Vec3f(1.0f)); +} +Vec3f QuadraticTransmittance::mediumMedium(const Vec3f &tau) const +{ + Vec3f result(1.0f/_maxT); + for (int i = 0; i < 3; ++i) + if (tau[i] > _maxT) + result[i] = 0.0f; + return result; +} + +float QuadraticTransmittance::sigmaBar() const +{ + return 2.0f/_maxT; +} + +float QuadraticTransmittance::sampleSurface(PathSampleGenerator &sampler) const +{ + return _maxT*(1.0f - std::sqrt(1.0f - sampler.next1D())); +} +float QuadraticTransmittance::sampleMedium(PathSampleGenerator &sampler) const +{ + return _maxT*sampler.next1D(); +} + +} diff --git a/src/core/transmittances/QuadraticTransmittance.hpp b/src/core/transmittances/QuadraticTransmittance.hpp new file mode 100644 index 00000000..73e6d6c0 --- /dev/null +++ b/src/core/transmittances/QuadraticTransmittance.hpp @@ -0,0 +1,31 @@ +#ifndef QuadraticTransmittance_HPP_ +#define QuadraticTransmittance_HPP_ + +#include "Transmittance.hpp" + +namespace Tungsten { + +class QuadraticTransmittance : public Transmittance +{ + float _maxT; + +public: + QuadraticTransmittance(); + + virtual void fromJson(JsonPtr value, const Scene &scene) override; + virtual rapidjson::Value toJson(Allocator &allocator) const override; + + virtual Vec3f surfaceSurface(const Vec3f &tau) const override final; + virtual Vec3f surfaceMedium(const Vec3f &tau) const override final; + virtual Vec3f mediumSurface(const Vec3f &tau) const override final; + virtual Vec3f mediumMedium(const Vec3f &tau) const override final; + + virtual float sigmaBar() const override final; + + virtual float sampleSurface(PathSampleGenerator &sampler) const override final; + virtual float sampleMedium(PathSampleGenerator &sampler) const override final; +}; + +} + +#endif /* QuadraticTransmittance_HPP_ */ diff --git a/src/core/transmittances/Transmittance.hpp b/src/core/transmittances/Transmittance.hpp new file mode 100644 index 00000000..d73be4fc --- /dev/null +++ b/src/core/transmittances/Transmittance.hpp @@ -0,0 +1,62 @@ +#ifndef TRANSMITTANCE_HPP_ +#define TRANSMITTANCE_HPP_ + +#include "sampling/PathSampleGenerator.hpp" + +#include "math/Vec.hpp" + +#include "io/JsonSerializable.hpp" + +#include "StringableEnum.hpp" + +#include +#include + +namespace Tungsten { + +class Transmittance : public JsonSerializable +{ +public: + virtual ~Transmittance() {} + + inline Vec3f eval(const Vec3f &tau, bool startOnSurface, bool endOnSurface) const + { + if (startOnSurface && endOnSurface) + return surfaceSurface(tau); + else if (!startOnSurface && !endOnSurface) + return mediumMedium(tau)/sigmaBar(); + else + return mediumSurface(tau); + } + inline float sample(PathSampleGenerator &sampler, bool startOnSurface) const + { + return startOnSurface ? sampleSurface(sampler) : sampleMedium(sampler); + } + inline Vec3f surfaceProbability(const Vec3f &tau, bool startOnSurface) const + { + return startOnSurface ? surfaceSurface(tau) : mediumSurface(tau); + } + inline Vec3f mediumPdf(const Vec3f &tau, bool startOnSurface) const + { + return startOnSurface ? surfaceMedium(tau) : mediumMedium(tau); + } + + virtual bool isDirac() const // Returns true if mediumMedium is a dirac/sum of dirac deltas + { + return false; + } + + virtual Vec3f surfaceSurface(const Vec3f &tau) const = 0; + virtual Vec3f surfaceMedium(const Vec3f &tau) const = 0; + virtual Vec3f mediumSurface(const Vec3f &tau) const = 0; + virtual Vec3f mediumMedium(const Vec3f &tau) const = 0; + + virtual float sigmaBar() const = 0; // Returns surfaceMedium(x)/mediumSurface(x) (i.e. identical to surfaceMedium(0)) + + virtual float sampleSurface(PathSampleGenerator &sampler) const = 0; + virtual float sampleMedium(PathSampleGenerator &sampler) const = 0; +}; + +} + +#endif /* TRANSMITTANCE_HPP_ */ diff --git a/src/core/transmittances/TransmittanceFactory.cpp b/src/core/transmittances/TransmittanceFactory.cpp new file mode 100644 index 00000000..08ddfb4d --- /dev/null +++ b/src/core/transmittances/TransmittanceFactory.cpp @@ -0,0 +1,26 @@ +#include "TransmittanceFactory.hpp" +#include "DoubleExponentialTransmittance.hpp" +#include "DavisWeinsteinTransmittance.hpp" +#include "InterpolatedTransmittance.hpp" +#include "ExponentialTransmittance.hpp" +#include "QuadraticTransmittance.hpp" +#include "ErlangTransmittance.hpp" +#include "LinearTransmittance.hpp" +#include "DavisTransmittance.hpp" +#include "PulseTransmittance.hpp" + +namespace Tungsten { + +DEFINE_STRINGABLE_ENUM(TransmittanceFactory, "transmittance", ({ + {"double_exponential", std::make_shared}, + {"exponential", std::make_shared}, + {"quadratic", std::make_shared}, + {"linear", std::make_shared}, + {"erlang", std::make_shared}, + {"davis", std::make_shared}, + {"davis_weinstein", std::make_shared}, + {"pulse", std::make_shared}, + {"interpolated", std::make_shared}, +})) + +} diff --git a/src/core/transmittances/TransmittanceFactory.hpp b/src/core/transmittances/TransmittanceFactory.hpp new file mode 100644 index 00000000..a9b20ba0 --- /dev/null +++ b/src/core/transmittances/TransmittanceFactory.hpp @@ -0,0 +1,17 @@ +#ifndef TRANSMITTANCEFACTORY_HPP_ +#define TRANSMITTANCEFACTORY_HPP_ + +#include "StringableEnum.hpp" + +#include +#include + +namespace Tungsten { + +class Transmittance; + +typedef StringableEnum()>> TransmittanceFactory; + +} + +#endif /* TRANSMITTANCEFACTORY_HPP_ */