diff --git a/data/tests/test_bsdf.xml b/data/tests/test_bsdf.xml index d962081b..f1c15892 100644 --- a/data/tests/test_bsdf.xml +++ b/data/tests/test_bsdf.xml @@ -1,27 +1,9 @@ - - - - - - - - - - - - - - - - - - - + + + @@ -158,6 +140,18 @@ + + + + + + + + + + @@ -200,4 +194,13 @@ + + + + + + + + + diff --git a/src/bsdfs/SConscript b/src/bsdfs/SConscript index a2665b39..681ab4e2 100644 --- a/src/bsdfs/SConscript +++ b/src/bsdfs/SConscript @@ -2,26 +2,26 @@ Import('env', 'plugins') # Basic library of smooth and rough materials plugins += env.SharedLibrary('diffuse', ['diffuse.cpp']) -#plugins += env.SharedLibrary('dielectric', ['dielectric.cpp']) -#plugins += env.SharedLibrary('conductor', ['conductor.cpp']) -#plugins += env.SharedLibrary('plastic', ['plastic.cpp']) -#plugins += env.SharedLibrary('roughdiffuse', ['roughdiffuse.cpp']) -#plugins += env.SharedLibrary('roughdielectric', ['roughdielectric.cpp']) -#plugins += env.SharedLibrary('roughconductor', ['roughconductor.cpp']) -#plugins += env.SharedLibrary('roughplastic', ['roughplastic.cpp']) +plugins += env.SharedLibrary('dielectric', ['dielectric.cpp']) +plugins += env.SharedLibrary('conductor', ['conductor.cpp']) +plugins += env.SharedLibrary('plastic', ['plastic.cpp']) +plugins += env.SharedLibrary('roughdiffuse', ['roughdiffuse.cpp']) +plugins += env.SharedLibrary('roughdielectric', ['roughdielectric.cpp']) +plugins += env.SharedLibrary('roughconductor', ['roughconductor.cpp']) +plugins += env.SharedLibrary('roughplastic', ['roughplastic.cpp']) # Materials that act as modifiers -#plugins += env.SharedLibrary('twosided', ['twosided.cpp']) -#plugins += env.SharedLibrary('mask', ['mask.cpp']) -#plugins += env.SharedLibrary('mixture', ['mixture.cpp']) -#plugins += env.SharedLibrary('coating', ['coating.cpp']) -#plugins += env.SharedLibrary('bump', ['bump.cpp']) +plugins += env.SharedLibrary('twosided', ['twosided.cpp']) +plugins += env.SharedLibrary('mask', ['mask.cpp']) +plugins += env.SharedLibrary('mixture', ['mixture.cpp']) +plugins += env.SharedLibrary('coating', ['coating.cpp']) +plugins += env.SharedLibrary('bump', ['bump.cpp']) # Other materials -#plugins += env.SharedLibrary('ward', ['ward.cpp']) -#plugins += env.SharedLibrary('phong', ['phong.cpp']) -#plugins += env.SharedLibrary('irawan', ['irawan.cpp']) -#plugins += env.SharedLibrary('difftrans', ['difftrans.cpp']) -#plugins += env.SharedLibrary('hk', ['hk.cpp']) +plugins += env.SharedLibrary('ward', ['ward.cpp']) +plugins += env.SharedLibrary('phong', ['phong.cpp']) +plugins += env.SharedLibrary('irawan', ['irawan.cpp']) +plugins += env.SharedLibrary('difftrans', ['difftrans.cpp']) +plugins += env.SharedLibrary('hk', ['hk.cpp']) Export('plugins') diff --git a/src/bsdfs/bump.cpp b/src/bsdfs/bump.cpp index 3c68c907..cdd0f635 100644 --- a/src/bsdfs/bump.cpp +++ b/src/bsdfs/bump.cpp @@ -191,28 +191,6 @@ public: return m_nested->pdf(perturbedQuery, measure); } - Spectrum sample(BSDFQueryRecord &bRec, Float &pdf, const Point2 &sample) const { - const Intersection& its = bRec.its; - Intersection perturbed; - perturbIntersection(its, perturbed); - - BSDFQueryRecord perturbedQuery(perturbed, bRec.sampler, bRec.quantity); - perturbedQuery.wi = perturbed.toLocal(its.toWorld(bRec.wi)); - perturbedQuery.typeMask = bRec.typeMask; - perturbedQuery.component = bRec.component; - Spectrum result = m_nested->sample(perturbedQuery, pdf, sample); - - if (!result.isZero()) { - bRec.sampledComponent = perturbedQuery.sampledComponent; - bRec.sampledType = perturbedQuery.sampledType; - bRec.wo = its.toLocal(perturbed.toWorld(perturbedQuery.wo)); - if (Frame::cosTheta(bRec.wo) * Frame::cosTheta(perturbedQuery.wo) <= 0) - return Spectrum(0.0f); - } - - return result; - } - Spectrum sample(BSDFQueryRecord &bRec, const Point2 &sample) const { const Intersection& its = bRec.its; Intersection perturbed; @@ -234,6 +212,28 @@ public: return result; } + Spectrum sampleXXX(BSDFQueryRecord &bRec, Float &pdf, const Point2 &sample) const { + const Intersection& its = bRec.its; + Intersection perturbed; + perturbIntersection(its, perturbed); + + BSDFQueryRecord perturbedQuery(perturbed, bRec.sampler, bRec.quantity); + perturbedQuery.wi = perturbed.toLocal(its.toWorld(bRec.wi)); + perturbedQuery.typeMask = bRec.typeMask; + perturbedQuery.component = bRec.component; + Spectrum result = m_nested->sampleXXX(perturbedQuery, pdf, sample); + + if (!result.isZero()) { + bRec.sampledComponent = perturbedQuery.sampledComponent; + bRec.sampledType = perturbedQuery.sampledType; + bRec.wo = its.toLocal(perturbed.toWorld(perturbedQuery.wo)); + if (Frame::cosTheta(bRec.wo) * Frame::cosTheta(perturbedQuery.wo) <= 0) + return Spectrum(0.0f); + } + + return result; + } + Shader *createShader(Renderer *renderer) const; std::string toString() const { diff --git a/src/bsdfs/coating.cpp b/src/bsdfs/coating.cpp index 6a2ab7e8..1061ec2b 100644 --- a/src/bsdfs/coating.cpp +++ b/src/bsdfs/coating.cpp @@ -307,7 +307,7 @@ public: } } - Spectrum sample(BSDFQueryRecord &bRec, Float &pdf, const Point2 &_sample) const { + Spectrum sampleXXX(BSDFQueryRecord &bRec, Float &pdf, const Point2 &_sample) const { bool sampleSpecular = (bRec.typeMask & EDeltaReflection) && (bRec.component == -1 || bRec.component == (int) m_components.size()-1); bool sampleNested = (bRec.typeMask & m_nested->getType() & BSDF::EAll) @@ -341,14 +341,14 @@ public: bRec.sampledType = EDeltaReflection; bRec.wo = reflect(bRec.wi); pdf = sampleNested ? probSpecular : 1.0f; - return Spectrum(R12); + return Spectrum(R12) / pdf; } else { if (R12 == 1.0f) /* Total internal reflection */ return Spectrum(0.0f); Vector wiBackup = bRec.wi; bRec.wi = wiPrime; - Spectrum result = m_nested->sample(bRec, pdf, sample); + Spectrum result = m_nested->sampleXXX(bRec, pdf, sample); bRec.wi = wiBackup; if (result.isZero()) return Spectrum(0.0f); @@ -366,21 +366,20 @@ public: if (R21 == 1.0f) /* Total internal reflection */ return Spectrum(0.0f); - if (sampleSpecular) - pdf *= 1 - probSpecular; + if (sampleSpecular) { + pdf *= 1.0f - probSpecular; + result /= 1.0f - probSpecular; + } result *= (1 - R12) * (1 - R21); if (BSDF::getMeasure(bRec.sampledType) == ESolidAngle) { /* Solid angle compression & irradiance conversion factors */ Float eta = m_extIOR / m_intIOR, etaSqr = eta*eta; - Float temp = Frame::cosTheta(bRec.wo) / Frame::cosTheta(woPrime); - - result *= etaSqr * - Frame::cosTheta(bRec.wi) / Frame::cosTheta(wiPrime) * temp; - pdf *= etaSqr * temp; - } + result *= Frame::cosTheta(bRec.wi) / Frame::cosTheta(wiPrime); + pdf *= etaSqr * Frame::cosTheta(bRec.wo) / Frame::cosTheta(woPrime); + } return result; } @@ -388,12 +387,7 @@ public: Spectrum sample(BSDFQueryRecord &bRec, const Point2 &sample) const { Float pdf; - Spectrum result = SmoothCoating::sample(bRec, pdf, sample); - - if (result.isZero()) - return Spectrum(0.0f); - else - return result / pdf; + return SmoothCoating::sampleXXX(bRec, pdf, sample); } Shader *createShader(Renderer *renderer) const; diff --git a/src/bsdfs/conductor.cpp b/src/bsdfs/conductor.cpp index 6765009e..ef92a668 100644 --- a/src/bsdfs/conductor.cpp +++ b/src/bsdfs/conductor.cpp @@ -255,7 +255,7 @@ public: fresnelConductor(Frame::cosTheta(bRec.wi), m_eta, m_k); } - Spectrum sample(BSDFQueryRecord &bRec, Float &pdf, const Point2 &sample) const { + Spectrum sampleXXX(BSDFQueryRecord &bRec, Float &pdf, const Point2 &sample) const { bool sampleReflection = (bRec.typeMask & EDeltaReflection) && (bRec.component == -1 || bRec.component == 0); diff --git a/src/bsdfs/dielectric.cpp b/src/bsdfs/dielectric.cpp index 1044f2ed..8c295cd2 100644 --- a/src/bsdfs/dielectric.cpp +++ b/src/bsdfs/dielectric.cpp @@ -387,7 +387,7 @@ public: } } - Spectrum sample(BSDFQueryRecord &bRec, Float &pdf, const Point2 &sample) const { + Spectrum sampleXXX(BSDFQueryRecord &bRec, Float &pdf, const Point2 &sample) const { bool sampleReflection = (bRec.typeMask & EDeltaReflection) && (bRec.component == -1 || bRec.component == 0); bool sampleTransmission = (bRec.typeMask & EDeltaTransmission) @@ -433,7 +433,7 @@ public: bRec.wo = reflect(bRec.wi); pdf = Fr; - return m_specularReflectance->getValue(bRec.its) * Fr; + return m_specularReflectance->getValue(bRec.its); } else { bRec.sampledComponent = 1; bRec.sampledType = EDeltaTransmission; @@ -447,7 +447,7 @@ public: /* When transporting radiance, account for the solid angle change at boundaries with different indices of refraction. */ return m_specularTransmittance->getValue(bRec.its) - * (1-Fr) * (bRec.quantity == ERadiance ? (eta*eta) : (Float) 1); + * (bRec.quantity == ERadiance ? (eta*eta) : (Float) 1); } } else if (sampleReflection) { bRec.sampledComponent = 0; diff --git a/src/bsdfs/difftrans.cpp b/src/bsdfs/difftrans.cpp index 1aacfd44..802a7ffe 100644 --- a/src/bsdfs/difftrans.cpp +++ b/src/bsdfs/difftrans.cpp @@ -99,7 +99,7 @@ public: return m_transmittance->getValue(bRec.its); } - Spectrum sample(BSDFQueryRecord &bRec, Float &pdf, const Point2 &sample) const { + Spectrum sampleXXX(BSDFQueryRecord &bRec, Float &pdf, const Point2 &sample) const { if (!(bRec.typeMask & m_combinedType)) return Spectrum(0.0f); bRec.wo = squareToHemispherePSA(sample); @@ -108,8 +108,7 @@ public: bRec.sampledComponent = 0; bRec.sampledType = EDiffuseTransmission; pdf = std::abs(Frame::cosTheta(bRec.wo)) * INV_PI; - return m_transmittance->getValue(bRec.its) - * (INV_PI * std::abs(Frame::cosTheta(bRec.wo))); + return m_transmittance->getValue(bRec.its); } void addChild(const std::string &name, ConfigurableObject *child) { diff --git a/src/bsdfs/hk.cpp b/src/bsdfs/hk.cpp index 50937cad..b33e84fb 100644 --- a/src/bsdfs/hk.cpp +++ b/src/bsdfs/hk.cpp @@ -272,7 +272,7 @@ public: return 0.0f; } - inline Spectrum sample(BSDFQueryRecord &bRec, Float &_pdf, const Point2 &_sample) const { + inline Spectrum sampleXXX(BSDFQueryRecord &bRec, Float &_pdf, const Point2 &_sample) const { AssertEx(bRec.sampler != NULL, "The BSDFQueryRecord needs to have a sampler!"); bool hasSpecularTransmission = (bRec.typeMask & EDeltaTransmission) @@ -307,7 +307,7 @@ public: bRec.wo = -bRec.wi; _pdf = hasSingleScattering ? probSpecularTransmission : 1.0f; - return eval(bRec, EDiscrete); + return eval(bRec, EDiscrete) / _pdf; } else { /* The glossy transmission/scattering component should be sampled */ bool hasGlossyReflection = (bRec.typeMask & EGlossyReflection) @@ -317,7 +317,7 @@ public: /* Sample According to the phase function lobes */ PhaseFunctionQueryRecord pRec(MediumSamplingRecord(), bRec.wi, bRec.wo); - m_phase->sample(pRec, _pdf, bRec.sampler); + m_phase->sampleXXX(pRec, _pdf, bRec.sampler); /* Store the sampled direction */ bRec.wo = pRec.wo; @@ -337,19 +337,14 @@ public: if (_pdf == 0) return Spectrum(0.0f); else - return eval(bRec, ESolidAngle); + return eval(bRec, ESolidAngle) / _pdf; } } Spectrum sample(BSDFQueryRecord &bRec, const Point2 &sample) const { - Float pdf = 0; - Spectrum result = HanrahanKrueger::sample(bRec, pdf, sample); - - if (result.isZero()) - return Spectrum(0.0f); - else - return result / pdf; + Float pdf; + return HanrahanKrueger::sampleXXX(bRec, pdf, sample); } void serialize(Stream *stream, InstanceManager *manager) const { diff --git a/src/bsdfs/irawan.cpp b/src/bsdfs/irawan.cpp index 05f2a903..2a2544c0 100644 --- a/src/bsdfs/irawan.cpp +++ b/src/bsdfs/irawan.cpp @@ -303,10 +303,10 @@ public: bRec.wo = squareToHemispherePSA(sample); bRec.sampledComponent = 0; bRec.sampledType = EGlossyReflection; - return eval(bRec, ESolidAngle) / Frame::cosTheta(bRec.wo); + return eval(bRec, ESolidAngle) * M_PI / Frame::cosTheta(bRec.wo); } - Spectrum sample(BSDFQueryRecord &bRec, Float &pdf, const Point2 &sample) const { + Spectrum sampleXXX(BSDFQueryRecord &bRec, Float &pdf, const Point2 &sample) const { bool hasSpecular = (bRec.typeMask & EGlossyReflection) && (bRec.component == -1 || bRec.component == 0) && m_ksMultiplier > 0; bool hasDiffuse = (bRec.typeMask & EDiffuseReflection) && @@ -321,7 +321,7 @@ public: bRec.sampledComponent = 0; bRec.sampledType = EGlossyReflection; pdf = Frame::cosTheta(bRec.wo) * INV_PI; - return eval(bRec, ESolidAngle); + return eval(bRec, ESolidAngle) / pdf; } void serialize(Stream *stream, InstanceManager *manager) const { diff --git a/src/bsdfs/mask.cpp b/src/bsdfs/mask.cpp index bdc8b722..1b7cef2b 100644 --- a/src/bsdfs/mask.cpp +++ b/src/bsdfs/mask.cpp @@ -140,42 +140,6 @@ public: } } - Spectrum sample(BSDFQueryRecord &bRec, Float &pdf, const Point2 &_sample) const { - Point2 sample(_sample); - Spectrum result(0.0f); - - Spectrum opacity = m_opacity->getValue(bRec.its); - Float prob = opacity.getLuminance(); - - bool sampleTransmission = bRec.typeMask & EDeltaTransmission - && (bRec.component == -1 || bRec.component == getComponentCount()-1); - bool sampleNested = bRec.component == -1 || bRec.component < getComponentCount()-1; - - if (sampleTransmission && sampleNested) { - if (sample.x <= prob) { - sample.x /= prob; - result = m_nestedBSDF->sample(bRec, pdf, sample) * opacity; - pdf *= prob; - } else { - bRec.wo = -bRec.wi; - bRec.sampledComponent = getComponentCount()-1; - bRec.sampledType = EDeltaTransmission; - pdf = 1-prob; - result = Spectrum(1.0f) - opacity; - } - } else if (sampleTransmission) { - bRec.wo = -bRec.wi; - bRec.sampledComponent = getComponentCount()-1; - bRec.sampledType = EDeltaTransmission; - pdf = 1; - result = Spectrum(1.0f) - opacity; - } else if (sampleNested) { - result = m_nestedBSDF->sample(bRec, pdf, sample) * opacity; - } - - return result; - } - Spectrum sample(BSDFQueryRecord &bRec, const Point2 &_sample) const { Point2 sample(_sample); Spectrum opacity = m_opacity->getValue(bRec.its); @@ -208,6 +172,42 @@ public: } } + Spectrum sampleXXX(BSDFQueryRecord &bRec, Float &pdf, const Point2 &_sample) const { + Point2 sample(_sample); + Spectrum result(0.0f); + + Spectrum opacity = m_opacity->getValue(bRec.its); + Float prob = opacity.getLuminance(); + + bool sampleTransmission = bRec.typeMask & EDeltaTransmission + && (bRec.component == -1 || bRec.component == getComponentCount()-1); + bool sampleNested = bRec.component == -1 || bRec.component < getComponentCount()-1; + + if (sampleTransmission && sampleNested) { + if (sample.x <= prob) { + sample.x /= prob; + result = m_nestedBSDF->sampleXXX(bRec, pdf, sample) * opacity / prob; + pdf *= prob; + } else { + bRec.wo = -bRec.wi; + bRec.sampledComponent = getComponentCount()-1; + bRec.sampledType = EDeltaTransmission; + pdf = 1-prob; + result = (Spectrum(1.0f) - opacity) / pdf; + } + } else if (sampleTransmission) { + bRec.wo = -bRec.wi; + bRec.sampledComponent = getComponentCount()-1; + bRec.sampledType = EDeltaTransmission; + pdf = 1; + result = Spectrum(1.0f) - opacity; + } else if (sampleNested) { + result = m_nestedBSDF->sampleXXX(bRec, pdf, sample) * opacity; + } + + return result; + } + void addChild(const std::string &name, ConfigurableObject *child) { if (child->getClass()->derivesFrom(MTS_CLASS(Texture)) && name == "opacity") m_opacity = static_cast(child); diff --git a/src/bsdfs/microfacet.h b/src/bsdfs/microfacet.h index 513902d5..b1126e22 100644 --- a/src/bsdfs/microfacet.h +++ b/src/bsdfs/microfacet.h @@ -240,86 +240,155 @@ public: Normal sample(const Point2 &sample, Float alphaU, Float alphaV) const { /* The azimuthal component is always selected uniformly regardless of the distribution */ - Float phiM = (2.0f * M_PI) * sample.y, - thetaM = 0.0f; - + Float cosThetaM = 0.0f, phiM = (2.0f * M_PI) * sample.y; + switch (m_type) { - case EBeckmann: - thetaM = std::atan(std::sqrt(-alphaU*alphaU * - std::log(1.0f - sample.x))); - break; - - case EGGX: - thetaM = std::atan(alphaU * std::sqrt(sample.x) / - std::sqrt(1.0f - sample.x)); + case EBeckmann: { + Float tanThetaMSqr = -alphaU*alphaU * std::log(1.0f - sample.x); + cosThetaM = 1.0f / std::sqrt(1 + tanThetaMSqr); + } break; - case EPhong: - thetaM = std::acos(std::pow(sample.x, (Float) 1 / - (alphaU + 2))); + case EGGX: { + Float tanThetaMSqr = alphaU * alphaU * sample.x / (1.0f - sample.x); + cosThetaM = 1.0f / std::sqrt(1 + tanThetaMSqr); + } break; - + + case EPhong: { + cosThetaM = std::pow(sample.x, 1 / (alphaU + 2)); + } + break; + case EAshikhminShirley: { /* Sampling method based on code from PBRT */ - Float phi, cosTheta; if (sample.x < 0.25f) { sampleFirstQuadrant(alphaU, alphaV, - 4 * sample.x, sample.y, phi, cosTheta); + 4 * sample.x, sample.y, phiM, cosThetaM); } else if (sample.x < 0.5f) { sampleFirstQuadrant(alphaU, alphaV, - 4 * (0.5f - sample.x), sample.y, phi, cosTheta); - phi = M_PI - phi; + 4 * (0.5f - sample.x), sample.y, phiM, cosThetaM); + phiM = M_PI - phiM; } else if (sample.x < 0.75f) { sampleFirstQuadrant(alphaU, alphaV, - 4 * (sample.x - 0.5f), sample.y, phi, cosTheta); - phi += M_PI; + 4 * (sample.x - 0.5f), sample.y, phiM, cosThetaM); + phiM += M_PI; } else { sampleFirstQuadrant(alphaU, alphaV, - 4 * (1 - sample.x), sample.y, phi, cosTheta); - phi = 2 * M_PI - phi; + 4 * (1 - sample.x), sample.y, phiM, cosThetaM); + phiM = 2 * M_PI - phiM; } - const Float sinTheta = std::sqrt( - std::max((Float) 0, 1 - cosTheta*cosTheta)); + } + break; + default: + SLog(EError, "Invalid distribution function!"); + } + + const Float sinThetaM = std::sqrt( + std::max((Float) 0, 1 - cosThetaM*cosThetaM)); + return Vector( + sinThetaM * std::cos(phiM), + sinThetaM * std::sin(phiM), + cosThetaM + ); + } + + + /** + * \brief Draw a sample from the microsurface normal distribution + * and return the associated probability density + * + * \param sample A uniformly distributed 2D sample + * \param alphaU The surface roughness in the tangent direction + * \param alphaV The surface roughness in the bitangent direction + * \param pdf The probability density wrt. solid angles + */ + Normal sample(const Point2 &sample, Float alphaU, Float alphaV, Float &pdf) const { + /* The azimuthal component is always selected + uniformly regardless of the distribution */ + Float cosThetaM = 0.0f; + + switch (m_type) { + case EBeckmann: { + Float tanThetaMSqr = -alphaU*alphaU * std::log(1.0f - sample.x); + cosThetaM = 1.0f / std::sqrt(1 + tanThetaMSqr); + Float cosThetaM2 = cosThetaM * cosThetaM, + cosThetaM3 = cosThetaM2 * cosThetaM; + pdf = (1.0f - sample.x) / (M_PI * alphaU*alphaU * cosThetaM3); + } + break; + + case EGGX: { + Float alphaUSqr = alphaU * alphaU; + Float tanThetaMSqr = alphaUSqr * sample.x / (1.0f - sample.x); + cosThetaM = 1.0f / std::sqrt(1 + tanThetaMSqr); + + Float cosThetaM2 = cosThetaM * cosThetaM, + cosThetaM3 = cosThetaM2 * cosThetaM, + temp = alphaUSqr + tanThetaMSqr; + + pdf = INV_PI * alphaUSqr / (cosThetaM3 * temp * temp); + } + break; + + case EPhong: { + Float exponent = 1 / (alphaU + 2); + cosThetaM = std::pow(sample.x, exponent); + pdf = (alphaU + 2) * INV_TWOPI * std::pow(sample.x, (alphaU+1) * exponent); + } + break; + + case EAshikhminShirley: { + Float phiM; + + /* Sampling method based on code from PBRT */ + if (sample.x < 0.25f) { + sampleFirstQuadrant(alphaU, alphaV, + 4 * sample.x, sample.y, phiM, cosThetaM); + } else if (sample.x < 0.5f) { + sampleFirstQuadrant(alphaU, alphaV, + 4 * (0.5f - sample.x), sample.y, phiM, cosThetaM); + phiM = M_PI - phiM; + } else if (sample.x < 0.75f) { + sampleFirstQuadrant(alphaU, alphaV, + 4 * (sample.x - 0.5f), sample.y, phiM, cosThetaM); + phiM += M_PI; + } else { + sampleFirstQuadrant(alphaU, alphaV, + 4 * (1 - sample.x), sample.y, phiM, cosThetaM); + phiM = 2 * M_PI - phiM; + } + const Float sinThetaM = std::sqrt( + std::max((Float) 0, 1 - cosThetaM*cosThetaM)), + sinPhiM = std::sin(phiM), + cosPhiM = std::cos(phiM); + + const Float exponent = alphaU * cosPhiM*cosPhiM + + alphaV * sinPhiM*sinPhiM; + pdf = std::sqrt((alphaU + 1) * (alphaV + 1)) + * INV_TWOPI * std::pow(cosThetaM, exponent); + return Vector( - sinTheta * std::cos(phi), - sinTheta * std::sin(phi), - cosTheta + sinThetaM * cosPhiM, + sinThetaM * sinPhiM, + cosThetaM ); } break; default: SLog(EError, "Invalid distribution function!"); } - - return Normal(sphericalDirection(thetaM, phiM)); + + const Float sinThetaM = std::sqrt( + std::max((Float) 0, 1 - cosThetaM*cosThetaM)); + Float phiM = (2.0f * M_PI) * sample.y; + return Vector( + sinThetaM * std::cos(phiM), + sinThetaM * std::sin(phiM), + cosThetaM + ); } - /** - * \brief Draw a sample from an isotropic microsurface normal - * distribution and return the magnitude of its 'z' component. - * - * \param sample A uniformly distributed number on [0,1] - * \param alphaU The surface roughness - */ - Float sampleIsotropic(Float sample, Float alpha) const { - switch (m_type) { - case EBeckmann: - return 1.0f / std::sqrt(1 + - std::abs(-alpha*alpha * std::log(1.0f - sample))); - - case EGGX: - return 1.0f / std::sqrt(1 + - alpha * alpha * sample / (1.0f - sample)); - - case EPhong: - return std::pow(sample, (Float) 1 / (alpha + 2)); - - default: - SLog(EError, "Invalid distribution function!"); - return 0.0f; - } - } - /** * \brief Smith's shadow-masking function G1 for each * of the supported microfacet distributions diff --git a/src/bsdfs/mixture.cpp b/src/bsdfs/mixture.cpp index 606983a7..cd47d25a 100644 --- a/src/bsdfs/mixture.cpp +++ b/src/bsdfs/mixture.cpp @@ -204,41 +204,6 @@ public: return result; } - Spectrum sample(BSDFQueryRecord &bRec, Float &pdf, const Point2 &_sample) const { - Point2 sample(_sample); - if (bRec.component == -1) { - /* Choose a component based on the normalized weights */ - size_t entry = m_pdf.sampleReuse(sample.x); - - Spectrum result = m_bsdfs[entry]->sample(bRec, pdf, sample); - if (result.isZero()) // sampling failed - return result; - - result *= m_weights[entry]; - pdf *= m_pdf[entry]; - - EMeasure measure = BSDF::getMeasure(bRec.sampledType); - for (size_t i=0; ipdf(bRec, measure) * m_pdf[i]; - result += m_bsdfs[i]->eval(bRec, measure) * m_weights[i]; - } - - bRec.sampledComponent += m_offsets[entry]; - return result; - } else { - /* Pick out an individual component */ - int requestedComponent = bRec.component; - int bsdfIndex = m_indices[requestedComponent].first; - bRec.component = m_indices[requestedComponent].second; - Spectrum result = m_bsdfs[bsdfIndex]->sample(bRec, pdf, sample) - * m_weights[bsdfIndex]; - bRec.component = bRec.sampledComponent = requestedComponent; - return result; - } - } - Spectrum sample(BSDFQueryRecord &bRec, const Point2 &_sample) const { Point2 sample(_sample); if (bRec.component == -1) { @@ -246,11 +211,11 @@ public: size_t entry = m_pdf.sampleReuse(sample.x); Float pdf; - Spectrum result = m_bsdfs[entry]->sample(bRec, pdf, sample); + Spectrum result = m_bsdfs[entry]->sampleXXX(bRec, pdf, sample); if (result.isZero()) // sampling failed return result; - result *= m_weights[entry]; + result *= m_weights[entry] * pdf; pdf *= m_pdf[entry]; EMeasure measure = BSDF::getMeasure(bRec.sampledType); @@ -275,6 +240,41 @@ public: } } + Spectrum sampleXXX(BSDFQueryRecord &bRec, Float &pdf, const Point2 &_sample) const { + Point2 sample(_sample); + if (bRec.component == -1) { + /* Choose a component based on the normalized weights */ + size_t entry = m_pdf.sampleReuse(sample.x); + + Spectrum result = m_bsdfs[entry]->sampleXXX(bRec, pdf, sample); + if (result.isZero()) // sampling failed + return result; + + result *= m_weights[entry] * pdf; + pdf *= m_pdf[entry]; + + EMeasure measure = BSDF::getMeasure(bRec.sampledType); + for (size_t i=0; ipdf(bRec, measure) * m_pdf[i]; + result += m_bsdfs[i]->eval(bRec, measure) * m_weights[i]; + } + + bRec.sampledComponent += m_offsets[entry]; + return result/pdf; + } else { + /* Pick out an individual component */ + int requestedComponent = bRec.component; + int bsdfIndex = m_indices[requestedComponent].first; + bRec.component = m_indices[requestedComponent].second; + Spectrum result = m_bsdfs[bsdfIndex]->sampleXXX(bRec, pdf, sample) + * m_weights[bsdfIndex]; + bRec.component = bRec.sampledComponent = requestedComponent; + return result; + } + } + void addChild(const std::string &name, ConfigurableObject *child) { if (child->getClass()->derivesFrom(MTS_CLASS(BSDF))) { BSDF *bsdf = static_cast(child); diff --git a/src/bsdfs/phong.cpp b/src/bsdfs/phong.cpp index 766e6e51..606d2b4a 100644 --- a/src/bsdfs/phong.cpp +++ b/src/bsdfs/phong.cpp @@ -175,7 +175,7 @@ public: return 0.0f; } - inline Spectrum sample(BSDFQueryRecord &bRec, Float &_pdf, const Point2 &_sample) const { + inline Spectrum sampleXXX(BSDFQueryRecord &bRec, Float &_pdf, const Point2 &_sample) const { Point2 sample(_sample); bool hasSpecular = (bRec.typeMask & EGlossyReflection) @@ -230,17 +230,12 @@ public: if (_pdf == 0) return Spectrum(0.0f); else - return eval(bRec, ESolidAngle); + return eval(bRec, ESolidAngle) / _pdf; } Spectrum sample(BSDFQueryRecord &bRec, const Point2 &sample) const { - Float pdf = 0; - Spectrum result = Phong::sample(bRec, pdf, sample); - - if (result.isZero()) - return Spectrum(0.0f); - else - return result / pdf; + Float pdf; + return Phong::sampleXXX(bRec, pdf, sample); } void addChild(const std::string &name, ConfigurableObject *child) { diff --git a/src/bsdfs/plastic.cpp b/src/bsdfs/plastic.cpp index d0a48ea2..a024eec5 100644 --- a/src/bsdfs/plastic.cpp +++ b/src/bsdfs/plastic.cpp @@ -182,17 +182,17 @@ public: Float probSpecular = 1.0f; if (hasSpecular && hasDiffuse) { - probSpecular = fresnel(Frame::cosTheta(bRec.wi), m_extIOR, m_intIOR); - probSpecular = (probSpecular*m_specularSamplingWeight) / - (probSpecular*m_specularSamplingWeight + - (1-probSpecular) * (1-m_specularSamplingWeight)); + Float Fr = fresnel(Frame::cosTheta(bRec.wi), m_extIOR, m_intIOR); + probSpecular = (Fr*m_specularSamplingWeight) / + (Fr*m_specularSamplingWeight + + (1-Fr) * (1-m_specularSamplingWeight)); } if (measure == EDiscrete && hasSpecular) { /* Check if the provided direction pair matches an ideal specular reflection; tolerate some roundoff errors */ if (std::abs(1 - dot(reflect(bRec.wi), bRec.wo)) < Epsilon) - return hasDiffuse ? probSpecular : 1.0f; + return probSpecular; } else if (measure == ESolidAngle && hasDiffuse) { return Frame::cosTheta(bRec.wo) * INV_PI * (hasSpecular ? (1 - probSpecular) : 1.0f); @@ -224,7 +224,7 @@ public: bRec.wo = reflect(bRec.wi); return m_specularReflectance->getValue(bRec.its) * - Fr / probSpecular; + (Fr / probSpecular); } else { bRec.sampledComponent = 1; bRec.sampledType = EDiffuseReflection; @@ -234,7 +234,7 @@ public: )); return m_diffuseReflectance->getValue(bRec.its) * - (1-Fr) / (1-probSpecular); + ((1-Fr) / (1-probSpecular)); } } else if (hasSpecular) { bRec.sampledComponent = 0; @@ -254,7 +254,7 @@ public: } } - Spectrum sample(BSDFQueryRecord &bRec, Float &pdf, const Point2 &sample) const { + Spectrum sampleXXX(BSDFQueryRecord &bRec, Float &pdf, const Point2 &sample) const { bool hasSpecular = (bRec.typeMask & EDeltaReflection) && (bRec.component == -1 || bRec.component == 0); bool hasDiffuse = (bRec.typeMask & EDiffuseReflection) @@ -276,7 +276,8 @@ public: bRec.wo = reflect(bRec.wi); pdf = probSpecular; - return m_specularReflectance->getValue(bRec.its) * Fr; + return m_specularReflectance->getValue(bRec.its) + * Fr / probSpecular; } else { bRec.sampledComponent = 1; bRec.sampledType = EDiffuseReflection; @@ -288,7 +289,7 @@ public: pdf = (1-probSpecular) * Frame::cosTheta(bRec.wo) * INV_PI; return m_diffuseReflectance->getValue(bRec.its) - * (INV_PI * Frame::cosTheta(bRec.wo) * (1-Fr)); + * (1-Fr) / (1-probSpecular); } } else if (hasSpecular) { bRec.sampledComponent = 0; @@ -304,11 +305,10 @@ public: return Spectrum(0.0f); bRec.wo = squareToHemispherePSA(sample); - + pdf = Frame::cosTheta(bRec.wo) * INV_PI; - return m_diffuseReflectance->getValue(bRec.its) - * (INV_PI * Frame::cosTheta(bRec.wo) * (1-Fr)); + return m_diffuseReflectance->getValue(bRec.its) * (1-Fr); } } diff --git a/src/bsdfs/roughconductor.cpp b/src/bsdfs/roughconductor.cpp index 21cbbb9b..ae423fff 100644 --- a/src/bsdfs/roughconductor.cpp +++ b/src/bsdfs/roughconductor.cpp @@ -298,7 +298,9 @@ public: m_alphaV->getValue(bRec.its).average()); /* Sample M, the microsurface normal */ - const Normal m = m_distribution.sample(sample, alphaU, alphaV); + Float microfacetPDF; + const Normal m = m_distribution.sample(sample, + alphaU, alphaV, microfacetPDF); /* Perfect specular reflection based on the microsurface normal */ bRec.wo = reflect(bRec.wi, m); @@ -316,14 +318,14 @@ public: * m_distribution.G(bRec.wi, bRec.wo, m, alphaU, alphaV) * dot(bRec.wi, m); - Float denominator = m_distribution.pdf(m, alphaU, alphaV) + Float denominator = microfacetPDF * Frame::cosTheta(bRec.wi); return m_specularReflectance->getValue(bRec.its) * F * (numerator / denominator); } - Spectrum sample(BSDFQueryRecord &bRec, Float &_pdf, const Point2 &sample) const { + Spectrum sampleXXX(BSDFQueryRecord &bRec, Float &pdf, const Point2 &sample) const { if (Frame::cosTheta(bRec.wi) < 0 || ((bRec.component != -1 && bRec.component != 0) || !(bRec.typeMask & EGlossyReflection))) @@ -336,7 +338,8 @@ public: m_alphaV->getValue(bRec.its).average()); /* Sample M, the microsurface normal */ - const Normal m = m_distribution.sample(sample, alphaU, alphaV); + const Normal m = m_distribution.sample(sample, + alphaU, alphaV, pdf); /* Perfect specular reflection based on the microsurface normal */ bRec.wo = reflect(bRec.wi, m); @@ -347,15 +350,23 @@ public: if (Frame::cosTheta(bRec.wo) <= 0) return Spectrum(0.0f); - /* Guard against numerical imprecisions */ - _pdf = pdf(bRec, ESolidAngle); + const Spectrum F = fresnelConductor(Frame::cosTheta(bRec.wi), + m_eta, m_k); - if (_pdf == 0) - return Spectrum(0.0f); - else - return eval(bRec, ESolidAngle); + Float numerator = m_distribution.eval(m, alphaU, alphaV) + * m_distribution.G(bRec.wi, bRec.wo, m, alphaU, alphaV) + * dot(bRec.wi, m); + + Float denominator = pdf * Frame::cosTheta(bRec.wi); + + /* Jacobian of the half-direction transform */ + pdf /= 4.0f * dot(bRec.wo, m); + + return m_specularReflectance->getValue(bRec.its) * F + * (numerator / denominator); } + void addChild(const std::string &name, ConfigurableObject *child) { if (child->getClass()->derivesFrom(MTS_CLASS(Texture))) { if (name == "alpha") diff --git a/src/bsdfs/roughdielectric.cpp b/src/bsdfs/roughdielectric.cpp index 609cf4e6..1f3d6ada 100644 --- a/src/bsdfs/roughdielectric.cpp +++ b/src/bsdfs/roughdielectric.cpp @@ -367,8 +367,8 @@ public: && (bRec.typeMask & EGlossyReflection)), hasTransmission = ((bRec.component == -1 || bRec.component == 1) && (bRec.typeMask & EGlossyTransmission)), - reflect = Frame::cosTheta(bRec.wi) - * Frame::cosTheta(bRec.wo) > 0; + reflect = Frame::cosTheta(bRec.wi) + * Frame::cosTheta(bRec.wo) > 0; /* Determine the appropriate indices of refraction */ Float etaI = m_extIOR, etaT = m_intIOR; @@ -385,10 +385,10 @@ public: return 0.0f; /* Calculate the reflection half-vector (and possibly flip it - so that it lies inside the hemisphere around the normal) */ + so that it lies inside the hemisphere around the surface normal) */ H = normalize(bRec.wo+bRec.wi) * signum(Frame::cosTheta(bRec.wo)); - + /* Jacobian of the half-direction transform */ dwh_dwo = 1.0f / (4.0f * dot(bRec.wo, H)); } else { @@ -460,8 +460,9 @@ public: #endif /* Sample M, the microsurface normal */ + Float microfacetPDF; const Normal m = m_distribution.sample(sample, - sampleAlphaU, sampleAlphaV); + sampleAlphaU, sampleAlphaV, microfacetPDF); Float F = fresnel(dot(bRec.wi, m), m_extIOR, m_intIOR), numerator = 1.0f; @@ -510,14 +511,14 @@ public: numerator *= m_distribution.eval(m, alphaU, alphaV) * m_distribution.G(bRec.wi, bRec.wo, m, alphaU, alphaV) * dot(bRec.wi, m); - - Float denominator = m_distribution.pdf(m, sampleAlphaU, sampleAlphaV) + + Float denominator = microfacetPDF * Frame::cosTheta(bRec.wi); return result * std::abs(numerator / denominator); } - Spectrum sample(BSDFQueryRecord &bRec, Float &_pdf, const Point2 &_sample) const { + Spectrum sampleXXX(BSDFQueryRecord &bRec, Float &pdf, const Point2 &_sample) const { Point2 sample(_sample); bool hasReflection = ((bRec.component == -1 || bRec.component == 0) @@ -546,15 +547,35 @@ public: #endif /* Sample M, the microsurface normal */ + Float microfacetPDF; const Normal m = m_distribution.sample(sample, - sampleAlphaU, sampleAlphaV); + sampleAlphaU, sampleAlphaV, microfacetPDF); + +#if 1 + Float ref = m_distribution.pdf(m, sampleAlphaU, sampleAlphaV); + if (std::abs(ref-microfacetPDF)/ref > Epsilon) + cout << "OOPS! ref=" << ref << ", got=" << microfacetPDF << endl; +#endif + + pdf = microfacetPDF; + + Float F = fresnel(dot(bRec.wi, m), m_extIOR, m_intIOR), + numerator = 1.0f; if (hasReflection && hasTransmission) { - Float F = fresnel(dot(bRec.wi, m), m_extIOR, m_intIOR); - if (bRec.sampler->next1D() > F) + if (bRec.sampler->next1D() > F) { choseReflection = false; + pdf *= (1-F); + } else { + pdf *= F; + } + } else { + numerator = hasReflection ? F : 1-F; } + Spectrum result; + Float dwh_dwo; + if (choseReflection) { /* Perfect specular reflection based on the microsurface normal */ bRec.wo = reflect(bRec.wi, m); @@ -564,6 +585,11 @@ public: /* Side check */ if (Frame::cosTheta(bRec.wi) * Frame::cosTheta(bRec.wo) <= 0) return Spectrum(0.0f); + + result = m_specularReflectance->getValue(bRec.its); + + /* Jacobian of the half-direction transform */ + dwh_dwo = 1.0f / (4.0f * dot(bRec.wo, m)); } else { /* Determine the appropriate indices of refraction */ Float etaI = m_extIOR, etaT = m_intIOR; @@ -580,17 +606,28 @@ public: /* Side check */ if (Frame::cosTheta(bRec.wi) * Frame::cosTheta(bRec.wo) >= 0) return Spectrum(0.0f); + + result = m_specularTransmittance->getValue(bRec.its) + * ((bRec.quantity == ERadiance) ? + ((etaI*etaI) / (etaT*etaT)) : (Float) 1); + + /* Jacobian of the half-direction transform. */ + Float sqrtDenom = etaI * dot(bRec.wi, m) + etaT * dot(bRec.wo, m); + dwh_dwo = (etaT*etaT * dot(bRec.wo, m)) / (sqrtDenom*sqrtDenom); } - /* Guard against numerical imprecisions */ - _pdf = pdf(bRec, ESolidAngle); + numerator *= m_distribution.eval(m, alphaU, alphaV) + * m_distribution.G(bRec.wi, bRec.wo, m, alphaU, alphaV) + * dot(bRec.wi, m); - if (_pdf == 0) - return Spectrum(0.0f); - else - return eval(bRec, ESolidAngle); + Float denominator = microfacetPDF * Frame::cosTheta(bRec.wi); + + pdf *= std::abs(dwh_dwo); + + return result * std::abs(numerator / denominator); } + void addChild(const std::string &name, ConfigurableObject *child) { if (child->getClass()->derivesFrom(MTS_CLASS(Texture))) { if (name == "alpha") diff --git a/src/bsdfs/roughdiffuse.cpp b/src/bsdfs/roughdiffuse.cpp index a39e9678..4ecc4df2 100644 --- a/src/bsdfs/roughdiffuse.cpp +++ b/src/bsdfs/roughdiffuse.cpp @@ -236,7 +236,7 @@ public: (Frame::cosTheta(bRec.wo) * INV_PI); } - Spectrum sample(BSDFQueryRecord &bRec, Float &pdf, const Point2 &sample) const { + Spectrum sampleXXX(BSDFQueryRecord &bRec, Float &pdf, const Point2 &sample) const { if (!(bRec.typeMask & EGlossyReflection) || Frame::cosTheta(bRec.wi) <= 0) return Spectrum(0.0f); @@ -244,7 +244,7 @@ public: bRec.sampledComponent = 0; bRec.sampledType = EGlossyReflection; pdf = Frame::cosTheta(bRec.wo) * INV_PI; - return eval(bRec, ESolidAngle); + return eval(bRec, ESolidAngle) / pdf; } void addChild(const std::string &name, ConfigurableObject *child) { diff --git a/src/bsdfs/roughplastic.cpp b/src/bsdfs/roughplastic.cpp index 1f3669bd..0a8de77f 100644 --- a/src/bsdfs/roughplastic.cpp +++ b/src/bsdfs/roughplastic.cpp @@ -305,7 +305,7 @@ public: return result; } - inline Spectrum sample(BSDFQueryRecord &bRec, Float &_pdf, const Point2 &_sample) const { + inline Spectrum sampleXXX(BSDFQueryRecord &bRec, Float &_pdf, const Point2 &_sample) const { bool hasSpecular = (bRec.typeMask & EGlossyReflection) && (bRec.component == -1 || bRec.component == 0); bool hasDiffuse = (bRec.typeMask & EDiffuseReflection) && @@ -359,17 +359,12 @@ public: if (_pdf == 0) return Spectrum(0.0f); else - return eval(bRec, ESolidAngle); + return eval(bRec, ESolidAngle) / _pdf; } Spectrum sample(BSDFQueryRecord &bRec, const Point2 &sample) const { - Float pdf = 0; - Spectrum result = RoughPlastic::sample(bRec, pdf, sample); - - if (result.isZero()) - return Spectrum(0.0f); - else - return result / pdf; + Float pdf; + return RoughPlastic::sampleXXX(bRec, pdf, sample); } void serialize(Stream *stream, InstanceManager *manager) const { diff --git a/src/bsdfs/twosided.cpp b/src/bsdfs/twosided.cpp index f6f8b022..04bc4a23 100644 --- a/src/bsdfs/twosided.cpp +++ b/src/bsdfs/twosided.cpp @@ -123,14 +123,14 @@ public: return result; } - Spectrum sample(BSDFQueryRecord &bRec, Float &pdf, const Point2 &sample) const { + Spectrum sampleXXX(BSDFQueryRecord &bRec, Float &pdf, const Point2 &sample) const { bool flipped = false; if (Frame::cosTheta(bRec.wi) < 0) { bRec.wi.z *= -1; flipped = true; } - Spectrum result = m_nestedBRDF->sample(bRec, pdf, sample); + Spectrum result = m_nestedBRDF->sampleXXX(bRec, pdf, sample); if (flipped) { bRec.wi.z *= -1; diff --git a/src/bsdfs/ward.cpp b/src/bsdfs/ward.cpp index 02c28076..70c6dcd7 100644 --- a/src/bsdfs/ward.cpp +++ b/src/bsdfs/ward.cpp @@ -262,7 +262,7 @@ public: return 0.0f; } - inline Spectrum sample(BSDFQueryRecord &bRec, Float &_pdf, const Point2 &_sample) const { + inline Spectrum sampleXXX(BSDFQueryRecord &bRec, Float &_pdf, const Point2 &_sample) const { Point2 sample(_sample); bool hasSpecular = (bRec.typeMask & EGlossyReflection) @@ -321,17 +321,12 @@ public: if (_pdf == 0) return Spectrum(0.0f); else - return eval(bRec, ESolidAngle); + return eval(bRec, ESolidAngle) / _pdf; } Spectrum sample(BSDFQueryRecord &bRec, const Point2 &sample) const { - Float pdf=0; - Spectrum result = Ward::sample(bRec, pdf, sample); - - if (result.isZero()) - return Spectrum(0.0f); - else - return result / pdf; + Float pdf; + return Ward::sampleXXX(bRec, pdf, sample); } void addChild(const std::string &name, ConfigurableObject *child) { diff --git a/src/tests/test_chisquare.cpp b/src/tests/test_chisquare.cpp index 1a0f8744..5c8c3288 100644 --- a/src/tests/test_chisquare.cpp +++ b/src/tests/test_chisquare.cpp @@ -128,6 +128,10 @@ public: EMeasure measure = ESolidAngle; if (!sampled.isZero()) measure = BSDF::getMeasure(bRec.sampledType); + + if (sampled.isZero() && sampled2.isZero()) + return boost::make_tuple(Vector(0.0f), 0.0f, measure); + Spectrum f = m_bsdf->eval(bRec, measure); pdfVal = m_bsdf->pdf(bRec, measure); Spectrum manual = f/pdfVal; @@ -178,8 +182,7 @@ public: disableFPExceptions(); #endif - return boost::make_tuple(bRec.wo, - sampled.isZero() ? 0.0f : 1.0f, measure); + return boost::make_tuple(bRec.wo, 1.0f, measure); } Float pdf(const Vector &wo, EMeasure measure) { @@ -393,7 +396,7 @@ public: Log(EInfo, "Checking the model for %i incident directions and 2D sampling", wiSamples); progress->reset(); - +#if 1 /* Test for a number of different incident directions */ for (size_t j=0; jupdate(j+1); } Log(EInfo, "The largest encountered importance weight was = %.2f", largestWeight); +#endif largestWeight = 0; if (bsdf->getComponentCount() > 1) {